*2.3. Mutation Screening*

DNA was obtained from the peripheral blood by phenol extraction or the salting-out method. Direct Sanger sequencing of twelve *LMNA* exons (canonical transcript NM\_170707.4) including flanking intronic regions was performed in 16 probands as described previously [15]. In one case, it was followed by multiplex ligation-dependent probe amplification (MRC Holland, Amsterdam, Netherlands) [18]. In 5 probands, next generation sequencing (NGS) was performed. Libraries were prepared using an Illumina: TruSeq Exome Enrichment Kit for whole exome sequencing (WES) in 1 proband, Nextera Rapid Capture Exome Library Kit with different enrichment probes for TruSight One (TSO) gene panels in 2 probands and TruSight Cardio (TSC) Sequencing Kit in 2 probands (Illumina, San Diego, CA, USA). WES and TSO libraries were sequenced on Illumina HiSeq1500 and TSC on Illumina MiSeq. Library preparation, sequencing, and data analysis were performed as described previously [19]. Variants identified in probands were followed up by Sanger sequencing of relatives' DNA using BigDye Terminator v3.1/v1.1 Cycle Sequencing Kit (Life Technologies) according to the

manufacturer's instructions and the 3500xL/3130xL Genetic Analyzer (Life Technologies, Carlsbad, CA, USA). The results were analysed with Variant Reporter 1.1 Software (Life Technologies). Non-missense mutations included insertions, deletions, nonsense mutations, or mutations affecting splicing.


*3.1. Target Problem*

The control improvement goal of this research is to resolve three problems. The first two problems are related to the limitation of genset's output power and the fluctuation of working points. The last is related to re-establishing the working points on the optimal efficiency trajectory according to the involvement of the interference torque.


3.4.2. Genetic Engineering

Biofortification through genetic engineering is an alternative approach when variation in the desired traits is not available naturally in the available germplasm, a specific micronutrient does not naturally exist in crops, and/or modifications cannot be achieved by conventional breeding [136,137]. This approach was supported by the availability of fully sequenced genomes in various crops in recent years. Along with increasing the concentration of micronutrients, this approach can also be targeted simultaneously for removal of antinutrients or inclusion of promoters that can enhance the bioavailability of micronutrients [93,106,138]. This approach had not only utilized genes associated with various metabolic pathways operated in plants, but also from bacteria and other organisms [139,140]. Development of transgenic crops requires a substantial investment during the initial stage, but this could be a sustainable approach that has the potential to target large populations, especially in developing countries [96,103,141].

Several crops have been successfully modified using a transgenic approach to overcome a micronutrient deficiency. For example, enhanced accumulation (3 to 4 times) of Fe was noted in rice via expression of the iron-storage protein, ferritin [142,143]. In rice, the national levels for Fe

and Zn biofortification nutrition targets were attained under field settings in the Philippines and Colombia [144]. The genetically engineered rice (golden rice) was developed to produce β-carotene (pro-vitamin A) to fight against vitamin A deficiency [145]. Recently, transgenic multivitamin corn was produced by the simultaneous modification of three distinct metabolic pathways to increase the levels of three vitamins, i.e., β-carotene (169-fold), ascorbate (6-fold), and folate (2-fold), in the endosperm, and this could pave the way to develop nutritionally complete cereals [146]. Using metabolic engineering, the folate concentration was increased in tomato and rice [85,147]. Storozhenko et al. [148] reported more than 100-fold increase in folate concentration in rice by overexpression of *Arabidopsis thaliana* pterin and para-aminobenzoate genes, precursors of the folate biosynthesis pathway, whereas Hossain et al. [149] reported a two- to four-fold increase in *Arabidopsis* by overexpression of the gene involved in pterin biosynthesis.

To the best of our knowledge, there are no examples of biofortified pulse crops developed through a transgenic approach for Fe, Zn, Se, I, carotenoids, or folates in the available literature. However, a genetic engineering approach has been applied in pulse crops for improvement of other nutritional profile. For example, the concentration of the essential amino acid methionine was significantly increased in transgenic common bean plants (up to 23%) by expression of a methionine-rich storage albumin from the Brazil nut [150], and in the concentration of methionine in transgenic lupins (up to 94%) by expressing a sunflower seed albumin gene [151].

In recent years, targeted gene editing technologies using artificial nucleases, zinc finger nucleases (ZFNs), transcription activator-like e ffector nucleases (TALENs), and the clustered regularly interspaced short palindromic repeat (CRISPR)/CRISPR-associated protein 9 (Cas9) system (CRISPR/Cas9) have given rise to the possibility to precisely modify genes of interest, and thus have potential application for crop improvement [152,153]. These technologies have been used in various crops including rice [154,155], wheat [156], and tomatoes [157]. Recently, CRISPR/Cas9 and TALENs technologies were used to generate mutant lines for genes involved in small RNA processing of *Glycine max* and *Medicago truncatula* [158]. Similarly, CRISPR/Cas9-mediated genome editing technology was used in cowpeas to successfully disrupt symbiotic nitrogen fixation (SNF) gene activation [159]. These findings pave the way for applicability of use of gene editing technologies for various traits of interest in legumes.


**Acknowledgements**

I am very grateful to Dr. Donald H. Kobe, Dr. Paolo Grigolini, Dr. Daniel P. Sheehan, Dr. Bruce N. Miller, and Dr. Marlan O. Scully and for many very helpful and thoughtful insights, as well as for very perceptive and valuable discussions and communications, which greatly helped my understanding of thermodynamics and statistical mechanics. Also, I am indebted to them, as well as to Dr. Bright Lowry, Dr. John Banewicz, Dr. Bruno J. Zwolinski, Dr. Roland E. Allen, Dr. Abraham Clearfield, Dr. Russell Larsen, Dr. James H. Cooke, Dr. Wolfgang Rindler, Dr. Richard McFee, Dr. Nolan Massey, and Dr. Stan Czamanski for lectures, discussions, and/or communications from which I learned very much concerning thermodynamics and statistical mechanics. I thank Dr. Stan Czamanski and Dr. S. Mort Zimmerman for the very interesting general scientific discussions over many years. I also thank Dan Zimmerman, Dr. Kurt W. Hess, and Robert H. Shelton for the very interesting general scientific discussions at times. Additionally, I thank Robert H. Shelton for very helpful advice concerning diction.


2.2.2. Random Forest

The random forest algorithm is based on decision trees, where the data are split into different branches that are created on the basis of specific data subsets. Random forest consists of creating multiple decision trees and randomizing the set of features these trees are fed into [48]. This approach is a trademark for what is known as ensemble learning. The response of each tree is then compared, and, in the case of classification, the mode of the outputs is considered as the categorical prediction. The diverse nature of the random forest algorithm, i.e., the use of multiple classifiers to find a robust prediction, allows for lowvariance responses, which is a desired characteristic in any machine learning method [49]. The predictions are also expected to be unbiased.

A particularity of this method is the identification and ranking of the most relevant features in the datasets with respect to the categorical responses. This can be useful as a complement for the study of the effect of single features on the output response.

The main disadvantage of this algorithm is the large computational time required for implementation, which increases with the number of trees to build (defined by the user). More details regarding the algorithm can be found in [50].


**Appendices and nomenclature**

which were all lower than the demands for lipid of aquatic animals [8]. For example, the recommended total lipid level in the diet for shrimp is in general higher than 6.5% [4]. Although external equipment could be used to settle the excessive part of biofloc, how to treat this deposit containing high content organic matter and

The efficiency for producing biofloc also needs to be elevated, if biofloc is used as a gradient of formulated feedstuff for replacement of fishmeal or soybean meal or used to feed to aquatic animals directly as a food with whole nutritional gradients usually contained in artificial feedstuff. The productivity of biofloc recent is not

Moreover, the improvement of biofloc palatability should be researched, which is important to the utilization of biofloc either eaten in situ or used as a food source [76]. Attractants or feeding promoting agents, like garlicin, betaine (trimethylglycine), trimethylamine oxide (TMAO), and s,s-dimethyl-β-propionic acid thetine (DMPT), could be taken into consideration as additives during development of biofloc in situ or preparing process for biofloc pellets. Thus, the effects of those agents on biofloc

Biofloc technology (BFT) supplies a novel solution for this issue without huge water exchange, even zero water exchange. In general, ammonia would be removed quickly within several hours in a BFT system. Moreover, because of the very high nutritional value for fish and shrimp, bioflocs, the by-product of BFT, could also be reused as a complemented food in situ or a gradient for feedstuff to replace expensive fishmeal, and biofloc also could be processed to formulate diet to feed fish and shrimp directly. However, some aspects with regard to the effective use of biofloc as a food source for fish and shrimp, such as high lipid content, productivity, and

This work is supported by the development funds of the Chinese central government to guide local science and technology (2017CT5013); the Sci-Tech program of Hunan province, China (2016NK2132); and the science and research program of

the Education Department of Hunan province, China (16C1085, 18B394).

The author also thanks the support from Collaborative Innovation Center (Hunan) for Efficient and Health Production of Fisheries, Hunan Engineering Research Center of Aquatic Organism Resources and Environmental Ecology, Hunan Engineering Research Center of Aquatic Organism Resources and Environmental Ecology, and Academician Workstation (Fisheries) of Hunan Province.

bacteria, part of which may be pathogens, was also a problem [39].

*Emerging Technologies, Environment and Research for Sustainable Aquaculture*

attraction to fish and shrimp should be studied in detail, respectively.

adequate for those uses in practical operations.

palatability, need to be further researched in detail.

The author declares no conflict of interest.

**Notes/thanks/other declarations**

**9. Conclusions**

**Acknowledgements**

**Conflict of interest**

**16**


**Examination and Modification of Multi-Factor Model in Explaining Stock Excess Return with Hybrid Approach in Empirical Study of Chinese Stock Market**

### **Jian Huang \* and Huazhang Liu**

Division of Business Management, Beijing Normal University-HongKong Baptist University United International College, Zhuhai 519087, China<EMAIL_ADDRESS>or<EMAIL_ADDRESS>**\*** Correspondence<EMAIL_ADDRESS>or<EMAIL_ADDRESS>Tel.: +86-186-6454-9728

Received: 28 April 2019; Accepted: 22 May 2019; Published: 25 May 2019

**Abstract:** To search significant variables which can illustrate the abnormal return of stock price, this research is generally based on the Fama-French five-factor model to develop a multi-factor model. We evaluated the existing factors in the empirical study of Chinese stock market and examined for new factors to extend the model by OLS and ridge regression model. With data from 2007 to 2018, the regression analysis was conducted on 1097 stocks separately in the market with computer simulation based on Python. Moreover, we conducted research on factor cyclical pattern via chi-square test and developed a corresponding trading strategy with trend analysis. For the results, we found that except market risk premium, each industry corresponds di fferently to the rest of six risk factors. The factor cyclical pattern can be used to predict the direction of seven risk factors and a simple moving average approach based on the relationships between risk factors and each industry was conducted in back-test which suggested that SMB (size premium), CMA (investment growth premium), CRMHL (momentum premium), and AMLH (asset turnover premium) can gain positive return.

**Keywords:** multi-factor model; risk factors; OLS and ridge regression model; python; chi-square test


**Hadar Kimelman 1,2 and Moshe Shemesh 1,\***

Received: 29 August 2019; Accepted: 23 September 2019; Published: 29 September 2019

**Abstract:** Live probiotic bacteria obtained with food are thought to have beneficial effects on a mammalian host, including their ability to reduce intestinal colonization by pathogens. To ensure the beneficial effects, the probiotic cells must survive processing and storage of food, its passage through the upper gastrointestinal tract (GIT), and subsequent chemical ingestion processes until they reach their target organ. However, there is considerable loss of viability of the probiotic bacteria during the drying process, in the acidic conditions of the stomach, and in the high bile concentration in the small intestine. *Bacillus subtilis*, a spore-forming probiotic bacterium, can effectively maintain a favorable balance of microflora in the GIT. *B. subtilis* produces a protective extracellular matrix (ECM), which is shared with other probiotic bacteria; thus, it was suggested that this ECM could potentially protect an entire community of probiotic cells against unfavorable environmental conditions. Consequently, a biofilm-based bio-coating system was developed that would enable a mutual growth of *B. subtilis* with different lactic acid bacteria (LAB) through increasing the ECM production. Results of the study demonstrate a significant increase in the survivability of the bio-coated LAB cells during the desiccation process and passage through the acidic environment. Thus, it provides evidence about the ability of *B. subtilis* in rescuing the desiccation-sensitive LAB, for instance, *Lactobacillus rhamnosus*, from complete eradication. Furthermore, this study demonstrates the antagonistic potential of the mutual probiotic system against pathogenic bacteria such as *Staphylococcus aureus*. The data show that the cells of *B. subtilis* possess robust anti-biofilm activity against *S. aureus* through activating the antimicrobial lipopeptide production pathway.

**Keywords:** beneficial biofilm; bio-coating; *B. subtilis*; probiotics; extracellular matrix; pathogen elimination


**1. Introduction**

We use recent advances in macroeconometrics and financial econometrics to investigate the macroeconomic effects of oil price shocks and oil price uncertainty. In doing so, they appeal to the real options theory (also known as investment under uncertainty), which predicts that firms are likely to delay making irreversible investment decisions in the face of uncertainty about the price of oil, particularly when the cash flow from investments is contingent on the oil price.

In particular, we investigate the effects of oil price uncertainty on total factor productivity growth, and in doing so, we control for two different channels for TFP growth—benefits from the quality of the national innovation system and from adopting new technologies. We use an unbalanced panel for 28 European Union countries over the period from 1990 to 2018. Consistent with the earlier literature, we find that oil price uncertainty has a negative and statistically significant effect on TFP growth even after we control for technological advancements and the effects of globalization.

According to the Solow growth model [1], the aggregate value added (or GDP) growth can be decomposed into contributions from the aggregate capital input, *K*, aggregate labour input, *L*, and aggregate total factor productivity, TFP, as follows:

$$
\Delta l n \text{GDP} = \nu\_K \Delta l nK + \nu\_L \Delta l nL + \Delta l n \text{TFP} \tag{1}
$$

where *νK* is the output elasticity of capital and *νL* is the output elasticity of labour. Under the assumption of constant return to scale *νK* + *νL* = 1, TFP driven by factors, such as technological progress, that are not tied to explicit input usage. In theory, TFP growth captures technical change and overall efficiency.

We consider two separate channels of technological change: benefits from the National Innovation System and benefits from adopting new technologies. Even within the European Union (EU), there are countries with full-cycle national innovation systems and countries who mainly adopt rather than invent modern technologies. The former countries are referred to as "innovation leaders" and the latter as "innovation followers", using (here) the Innovation Union Scoreboard terminology in a slightly different manner. Thus, we consider different sets of key indicators to assess technological progress in different groups.

However, not only technological progress affects TFP growth. Globalization and openness of the modern economies suggest that there is an impact of the global business cycle on productivity. The linkage between the stock market volatility and business cycle was investigated in the literature (see, for example [2] and literature review in [3]). In addition, Elder and Serletis [4] and Serletis and Xu [5], among others, argue that volatility in oil prices has had a negative and statistically significant effect on several investment measures, durables consumption, and aggregate output.

Thus, we investigate the effects of oil price volatility on TFP growth, controlling for variables measuring technological progress and openness of the economy. Our study examines two central research questions. Firstly, we propose a regression model for TFP growth, combining different indicators depending on different innovation strategy types. Secondly, we assess uncertainty in commodity markets.

The rest of the paper is organized as follows. Sections 2 and 3 present the methodology used in the empirical investigation. Section 4 presents the data and the empirical results. Section 5 provides a discussion of the results and addresses the policy implications. The final section concludes the paper.

#### **2. Research Model**

We use the Conference Board Total Economy Database (TED) as a source of data for TFP growth (The Conference Board, 2019. The Conference Board Total Economy Database ™, April 2019, http://www.conference-board.org/data/economydatabase/ accessed on: 10 March 2020). The database contains time series data for more than 120 countries, covering the period since 1990.

According to Jorgenson and co-authors [6,7], capital services and labour input are measured as translog aggregates of heterogeneous capital and labour types. Under the assumptions of competitive markets, full input utilization, and constant returns to scale, the contribution of each input to output equals the share of input cost to total cost—see Measuring Productivity [8].

TED capital is decomposed into Information and Communication Technology (ICT) capital and non-ICT-capital. Labour is decomposed into pure employment quantity and labour quality. One can find a detailed description of the sources and methods used to construct all TED variables in materials provided on the official website.

According to endogenous growth theory, TFP can be modelled as a function of a country's innovation capacity (see, for example, [9,10]). Innovation capacity depends upon the size and quality of the national innovation system, openness to international trade, the degree of technological specialization, and the ability to adopt and commercialize new-to-the-world technology.

We specify the TFP growth equation to reflect the role of the national innovation system as a source of productivity growth, the impact of technology spillovers, the role of international trade and FDI, and "the health" of the global economy as follows:

$$
\Delta \ln \text{TFP}\_{it} = \sum\_{j} \beta^{j} X^{j}\_{it} + \alpha\_{0} + \mathfrak{a}\_{i} + \delta Z\_{t} + \varepsilon\_{it} \tag{2}
$$

where *εit* is an error term, *i* denotes a country, *t* denotes time, *αi* captures country fixed effects, and *Zt* accounts for the "health" of the global economy at time *t*.

As TFP growth is measured from the supply side, period effects in TFP growth models can capture global demand changes. We use oil price volatility to capture the 'health' of the economy and changes in demand instead of period fixed effect.

There are several channels through which oil prices may affect productivity. Firstly, at the end of 2019, before the pandemic, the global economy consumed around 100 million barrels of crude oil per day, compared to about 40 million barrels at the end of the 1970s (https://www.iea.org/reports/world-energy-outlook-2019. accessed on: 1 December 2020). In the case of stable oil prices, firms can more accurately plan their expenses and investments. Rising oil prices reduce the availability of inputs and lead to output decreases (see Hamilton [11]). The impact of oil price movements on GDP and several macroeconomic and financial variables for the United States economy has been widely investigated in the literature (see, for example, Barsky and Kilian [12], Kilian and Vigfusson [13], and Azad and Serletis [14], among others). In fact, Azad and Serletis [14] find the linkage between oil price uncertainty and macroeconomic indicators of emerging economies to be significant.

Secondly, financial and commodity markets are good indicators of the state of the global economy. The correlation between commodity (notably, crude oil) prices and equity prices after the global crisis has been established in the literature (see Lombardi and Ravazzolo [15]).

High volatility of financial and commodity markets causes fear among portfolio investors, and may increase risk aversion. In turn, it may reduce investments in risky innovation projects and slow down technological progress.

For EU countries, uncertainty in commodity markets is even more important, as these countries are pure importers of fossil fuels and other raw materials. According to Eurostat, (https://ec.europa.eu/eurostat/cache/infographs/energy/bloc-2c.html accessed on: 3 March 2021) the dependency rate (the share of net imports in gross inland energy consumption) was equal to 58% in 2018. The dependency rate on energy imports has increased since 2000, when it was just 56%, despite the increasing efforts to achieve the renewable energy directive target by 2020. Several barriers that were found for the United States and Russia (see [16]) are also relevant for EU countries (see, for example, [17]).

The vector *Xj* in Equation (2) captures a country's technological advancement and the effect of globalization.

R&D is one of the leading indicators of a national innovation system and potential readiness for technological advances. As R&D expenditure comprises all expenditure on research and development in business enterprises, the government sector, higher education and non-profit firms, it captures all national efforts on product and process innovations quite well. It is therefore suited as a proxy for technological innovation. In the EU, it is one of the main indicators of achieving the strategic goal to grow through innovations. We examine the impact of the growth rate of R&D expenditure on TFP growth to reduce heteroscedasticity in our model.

Globalization is one of the main trends of the last three decades, considered in this study. Public policy, which was changing in the 1990s and 2000s towards a free-market economic system, and communications technology innovations have been identified as the two primary driving factors of globalization. Nowadays, there is growing uncertainty in industrialized countries regarding whether globalization means more opportunity or more risk. Trump, Brexit and increasing populism are named direct consequences of this development (see [18]). However, we can't deny the impact of globalization on economic growth in general, and TFP growth in particular. One of our objectives in this study is to investigate the linkage between globalization and TFP growth.

There are several ways to measure globalization. For example, the KOF Globalization Index measures the economic, social, and political dimensions of globalization. (https:// www.kof.ethz.ch/en/forecasts-and-indicators/indicators/kof-globalisation-index.html. accessed on: 03 March 2021) We use factors included in the KOF Globalisation Index from an economic perspective, such as trade openness and free capital movement. Note that European countries are among the leaders in the globalization ranking.

Participation in international trade may be regarded as a driver of productivity increases (see [19]). Firstly, only efficient and highly productive firms can be successful in world markets. So, the increase in exports can be regarded as a sign of increasing productivity. On the other hand, EU firms mainly import raw materials and intermediate products

to produce higher value-added products. So, the value of imports is also associated with productivity, at least in highly industrialized countries.

Secondly, for trading on the world market, the developed infrastructure is essential. To compete on global markets, firms need to be efficient in logistics, inventive in business processes, and successful in adopting advanced technologies. Thus, we consider the growth rate of both exports and imports to be a good candidate in explaining productivity dynamics for EU countries.

As the impact of FDI on TFP growth has been supported by many researchers, we examine FDI inflows and outflows. FDI inflows can bring technological, marketing and organizational innovations, at least in theory, affecting TFP growth. However, for the EU countries, the role of FDI outflows as a source of technology spillovers seems to also be important, as a sizable proportion of the EU outward flows typically are destined for the United States and other European countries (for example, Switzerland and Iceland) ( https://ec.europa.eu/eurostat/web/products-datasets/-/bop\_fdi6\_geo accessed on: 3 March 2021). Moreover, transnational corporations (TNCs) benefit from investing in emerging countries [20], which can indirectly affect the TFP of the home country of the transnational corporations. Globalization is also characterized by an ongoing fragmentation of production [21]. More efficient production chains, achieved by lowering transaction costs while investing in countries with lower factor costs, are assumed to increase productivity.

The technological spillovers can also be captured through flows of payments for the use of intellectual property. We use charges for the benefit of intellectual property rights paid and received by a country to reflect the scale of new technology commercialization on the one hand and the scale of new technology adoption on the other hand.

#### **3. Oil Price Volatility**

To measure oil price volatility, we follow the procedure used in [22]. Our first step is to detrend the Brent oil price series and obtain the cyclical component, denoted as *Brentcycl*. The volatility of oil price is then measured by the conditional variance of the forecast of the cyclical component of the oil price as follows:

$$Brent^{cyl}{}\_{t} = \gamma\_0 + \gamma\_0 \sigma\_l + \varepsilon\_l \sigma^2{}\_{t} = \omega + \kappa \varepsilon^2{}\_{t-1} + \beta \sigma^2{}\_{t-1} \tag{3}$$

where *σ*2 *t* is the one-period ahead forecast variance based on past information (conditional variance), *εt* is a conditionally normal innovation. Thus, we use the GARCH-M(1,1) model, as in for example, [4,23,24], to measure oil price volatility. The GARCH-based measure of volatility is common in the empirical literature (see for example [25,26]).

#### **4. Empirical Evidence**

#### *4.1. Data Description and Analysis*

We retrieve data from TED, Eurostat, the OECD database, and the World Bank database. We use an unbalanced panel of 28 EU countries (including all EU members as of the end of 2019) for the period 1990–2018. Data on R&D expenditures, international trade and royalties are converted to constant PPP dollars. The growth rate is measured in log changes multiplied by 100.

To visualize the dynamics of TFP at a country level, we compute a country-specific TFP index. We create the index setting the year 1990 TFP level equal to 100 for each country. We show the TFP dynamics in panel (a) of Figure 1 for the EU15 countries and in panel (b) of Figure 1 for the EU13 countries, using data from the TED database.

**Figure 1.** (**a**) TFP dynamics by country, EU15 countries, 1990–2018. (**b**) TFP dynamics by country, EU13 countries, 1990–2018.

The dynamics of productivity vary significantly among the 'old' members of the EU. Ireland has the highest productivity gains, followed by Finland, Sweden and Germany (see panel (a) of Figure 1). The stagnating productivity trends in Spain and Italy are evident, and Greece shows a sharp decline in productivity in the 2010s.

The majority of the former communist states have shown an upward productivity trend since the early 2000s. For example, Romania is one of the most fast-growing EU countries in terms of TFP (see panel (b) of Figure 1). However, the TFP indexes for Bulgaria and Croatia are far below 100 in 2018, meaning that TFP in these countries was lower in 2018 compared to 1990. (It is to be noted that the TFP series for Croatia is based on TED estimates using data for Yugoslavia. http://www.conference-board.org/data/ economydatabase/ accessed on: 10 March 2020) Notably, all countries display a noticeable decrease in productivity in 2009; however, the recovery rates are different.

To measure the uncertainty in commodity markets, we estimate the unknown parameters of the GARCH-M(1,1) model (3), using monthly observations of the Brent oil prices (see Figure 2).

**Figure 2.** Trend and cyclical components of the Brent spot oil price. (FOB, dollars per barrel).

We obtain a one-period ahead forecast of monthly volatility from the GARCH-M model (3). The volatility *VolT* of the year *T* is measured as the monthly average of the one-year window:

$$\operatorname{Vol}\_{T} = \frac{1}{12} \sum\_{j=1}^{12} \partial\_{\overline{T},j}^{2}$$

where *σ*ˆ 2*T*,*j* is an estimate of the GARH component of Equation (3) for year *T* at month *j*. We plot the obtained data in Figure 3. The sharp peak corresponds to the most significant uncertainty in oil and financial markets during the global financial crisis. This peak is closely associated with the consecutive steep decline in TFP growth of EU countries (see Figure 1).

**Figure 3.** Oil price uncertainty from the GARCH-in-Mean model.

The model key variables' description and summary statistics are reported in Table 1.

**Table 1.** Variable's description and descriptive statistics.

Note: Oil price volatility is defined for each of the 28 years and is a common factor for any country. The estimated values from the GARCH-M model are divided by 100.

> Before running the regression model, we proceed with unit root tests for panel data. We use the the Levin, Lin, and Chu [27] test, assuming a common unit root proccess, and the Im, Pesaran, and Shin [28] ADF test, assuming individual unit root process in the considered time-series. Results are reported in Table 2.

**Table 2.** Unit root test.

#### *4.2. Estimation Results*

The empirical analysis is conducted for 28 EU countries. Firstly, we estimate our model (2) with a single variable in vector X, namely the growth rate of real R&D expenditures, and oil price volatility, assuming the cross-section fixed effect. The Haussman test has shown that the random effects do not over-perform the fixed effects specification. The fixed effects specification suits better our analysis as we are interested in the country-specific part of the growth rate of TFP (parameter *αi* in Equation (2)).

Then, we check whether the growth rate of international trade and changes in FDI inflows and outflows could improve the model specification. To reflect the impact of knowledge spillovers on total factor productivity growth, we include the royalty variable. It is calculated as the sum of receipts and payments for the use of intellectual property rights.

> trade

Table 3 summarizes the estimated impact of R&D expenditures, internationaland royalties, changes in FDI flows and oil price volatility on TFP growth.

**Table 3.** The estimated effect of global and domestic factors on TFP growth for EU28 (refers to EU members as of the end of 2019).

Notes: White cross-section standard errors in parentheses, \*\*\* *p*-value < 0.01, \*\* *p*-value < 0.05, \* *p*-value < 0.10.

> In order to analyze whether the effects of technological innovations, openness and knowledge spillovers differ between countries with different historical and economic development, the countries are divided into two groups: "old" and "new" EU members, which is relatively common in empirical analyses.

> To test for equality of means of TFP growth rates between the EU15 and EU13 countries for the 1990–2018 period, we run the Satterthwaite–Welch *t*-test, which allows for unequal variances. We cannot reject the null hypothesis that the mean growth rate is equal for these two groups at the 5% significance level (the test statistic is equal to 1.9 with a *p*-value of 0.056). However, we assume that the impact of different factors on total factor productivity differs between countries with mature and developing markets.

## *4.3. Discussion*

All "old" EU countries (except Finland) have negative country-specific "fixed" growth rates of TFP, although, in the case of Austria and Germany, the coefficients are very close to zero. All "new" EU countries have positive country-specific 'fixed' growth rates of TFP, with Latvia and Lithuania having the highest values. The effect is partly due to a "low base" TFP in the former Soviet countries.

R&D expenditure growth has an impact on productivity in all EU countries. The openness of the economy also helps to boost productivity growth in all EU countries. This is in line with the findings of the "Globalization Report 2018: Who Benefits Most from Globalization?" It was shown that "for the third time in a row, as in 2014 and 2016, when measured in terms of real gross domestic product (GDP) per capita, industrialized countries continue to be the biggest winners of increasing globalization, while developing and emerging economies lag behind".

Our findings of the positive impact of FDI outflows on productivity growth are consistent with the results obtained by Altomonte and Ottaviano [29] in their study of the role of international production sharing in EU productivity at the micro-level. Royalties are more important for the EU13 countries than for the EU15 countries. Among the EU13 countries, almost all countries except Hungary, pay a lot more for the intellectual property rights than they receive from other countries.

## **5. Conclusions**

Total factor productivity measures the overall efficiency of labour and capital in the production process. During the examined period, from the early 1990s to the late 2010s, the EU countries showed a growth trend in the total factor productivity on average.

We assume that the total factor productivity is influenced by several groups of factors. Firstly, it is research and development carried out in a given country. For developed countries (most EU countries belong to this group), innovation is the most important driver of economic growth. We use internal R&D expenditure as a proxy variable for assessing the scale of research and innovation. We find that growth in R&D spending is associated with an increase in the TFP growth rate, all other things being equal.

The second group of factors includes indicators of the openness of the economy. These include the volume of foreign trade and inbound and outbound foreign direct investment. High competition in international markets increases the demands on companies to be more efficient and productive. Thus, it contributes to an increase in the overall productivity of the economy. We find a significant positive effect of international trade on TFP growth.

The impact of FDI inflow on productivity in developing countries is well documented in the literature. The positive effect of FDI-trade linkage lies in its contribution to integrating the host country into the world economy. The same is relevant for the former communist European countries which experienced a transition to a capitalist economic system during the 1990s. The impact of FDI outflow on the productivity of the home country is less examined in the empirical literature. Our estimates show that, controlling other variables, the impact of FDI inflow on TFP growth is not significant; however, the impact of FDI outflow is positive and significant.

However, it is not only innovation, technology and the overall efficiency of companies that drive productivity growth. Uncertainty on world commodity and financial markets plays a significant role and can have negative effects on economic growth. The volatility of the oil market, which is closely related to the volatility of other markets, inhibits the growth of total factor productivity, as shown in the study.

In this paper, we have focused on the effects of oil price shocks and oil price uncertainty on total factor productivity. Assessing the importance of oil price shocks, by simultaneously evaluating the effects monetary policy shocks, fiscal shocks, and other measures of risk and uncertainty, is an area for potentially productive research.

**Author Contributions:** Conceptualization, S.B. and A.S.; methodology, S.B.; software, S.B.; validation, S.B. and A.S.; formal analysis, S.B.; investigation, S.B. and A.S.; writing—original draft preparation, S.B. and A.S.; writing—review and editing, A.S.; visualization, A.S.; supervision, A.S.; project administration, A.S.; funding acquisition, S.B. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the RUDN University Strategic Academic Leadership Program.

**Institutional Review Board Statement:** Not applicable.

**Acknowledgments:** We would like to thank the Editor and two referees for comments that greatly improved the paper.

**Conflicts of Interest:** The authors declare no conflict of interest.


2. Unsteady flow due to sudden movement of the plate

#### 2.1. Governing equations

An infinite length plate is set in a stationary fluid as an initial condition. Let us consider the situation that the infinite length plate suddenly moves along its parallel direction at a constant speed uw. This problem was first solved by Stokes [2] in his famous treatment of the pendulum. Since Lord Rayleigh [3] also treated this flow, it is often called the Rayleigh problem in the literature. One takes that x is the plate movement direction and y is distance from the plate. Since the velocity component perpendicular to the plate v is zero, the momentum equation is simplified and is shown as a diffusion equation

$$\frac{\partial u}{\partial t} = \nu \frac{\partial^2 u}{\partial y^2} \tag{1}$$

Here, u is the velocity component parallel to the plate direction, t is the time, and ν is the kinematic viscosity. The boundary conditions for this partial differential equation are as follows:

$$\begin{cases} y = 0: & \quad \mu = \mu\_w \\ y \to \ast : & \quad \mu \to 0 \end{cases} \tag{2}$$

In order to reduce the partial differential equation to an ordinary equation, the following dimensionless velocity U and the similar variable η are introduced

$$
\mu = \mu\_w \mathcal{U}(\eta), \quad \eta = \frac{y}{2\sqrt{\nu t}} \tag{3}
$$

Then, the following ordinary differential equation can be obtained

$$\frac{d^2\mathcal{U}}{d\eta^2} + 2\eta \frac{d\mathcal{U}}{d\eta} = 0\tag{4}$$

The boundary condition for the ordinary differential equation is as follows using the similar variable η instead of y:

$$\begin{cases} \eta = 0 : \quad \mathcal{U} = 1\\ \eta \to \stackrel{\sim}{\sim} : \quad \mathcal{U} \to 0 \end{cases} \tag{5}$$

As a consequence, one needs to solve this boundary value problem. The theoretical solution can be easily obtained and expressed using the error function

$$\mathcal{U} = \frac{\mu}{\mu\_w} = 1 - \text{erf}(\eta) = 1 - \frac{2}{\sqrt{\pi}} \int\_0^{\eta} \exp\left(-\xi^2\right) d\xi \tag{6}$$

The velocity profile is shown in Figure 1.

be handled and it is valid when the flow velocity is much slower than the sound velocity and/ or the temperature difference in the fluid is small enough to consider the thermal expansion coefficient is independent to the temperature. The former situation is known the low Mach

Another simplification on the incompressible flows is the reduction of dimension due to the characteristic of similarity and periodicity. For the boundary layer flows such as the Blasius flow, the stagnation-point flow, and the von Kármán rotating disk flow have the similar solution where the flow transition from laminar to turbulence does not occur. In those cases, a combined dimensionless variable (similar variable) η is introduced and the velocity distribution can be only a function of η. While for the onset of instability such as the Rayleigh-Bénard convection, the Bénard-Marangoni convection, and the Taylor-Couette flow, the periodic characteristic of flow structure is observed. At the stage of onset of instability, the non-linear term is negligible and therefore the function of flow field is separated into the amplitude part and periodic part, respectively. This makes the effort on numerical analysis to reduce significantly

This chapter consists of three main bodies. First, a numerical technique for solving the boundary value problem called the first Stokes problem or the Rayleigh problem [1] is introduced. The differential equation is transferred into an ordinary equation and it is solved by a finite difference method using the Jacobi method. Second, similar solution of natural convection heat transfer heated from a vertical plate with uniform heat flux is introduced together with the method how to obtain the system of ordinary differential equations. The obtained Nusselt numbers are compared with some previous studies. Third, for example, of the linear stability analysis, one shows that the HSMAC method can be applied to obtain the critical values for the onset of secondary flow such

as the Taylor-Couette flow. The Eigen functions of flow and pressure fields are visualized.

An infinite length plate is set in a stationary fluid as an initial condition. Let us consider the situation that the infinite length plate suddenly moves along its parallel direction at a constant speed uw. This problem was first solved by Stokes [2] in his famous treatment of the pendulum. Since Lord Rayleigh [3] also treated this flow, it is often called the Rayleigh problem in the literature. One takes that x is the plate movement direction and y is distance from the plate. Since the velocity component perpendicular to the plate v is zero, the momentum equation is

> ∂u ∂t ¼ ν

∂2 u

Here, u is the velocity component parallel to the plate direction, t is the time, and ν is the kinematic viscosity. The boundary conditions for this partial differential equation are as

∂y2 (1)

number approximation, while the latter one the Boussinesq approximation.

244 Finite Element Method - Simulation, Numerical Analysis and Solution Techniques

and also to contribute the augmentation of accuracy of the results.

2. Unsteady flow due to sudden movement of the plate

2.1. Governing equations

follows:

simplified and is shown as a diffusion equation

Figure 1. Velocity profile.

### 2.2. Numerical method for solving the ordinary differential equation using finite difference method

For numerical solution, it is necessary to define the range of η, as recognized from Figure 1, η = 4 is enough. Hence, the boundary condition shown below is used instead of Eq. (5)

$$\begin{cases} \eta = 0 : & \mathcal{U} = 1 \\ \eta = 4 : & \mathcal{U} = 0 \end{cases} \tag{7}$$

η ¼ 0 : U1 ¼ 1 η ¼ 4 : UN ¼ 0

In the following, the case of N = 7 is considered, for example. By substituting i = 2–6 into Eq. (8),

1

0

U2 U3 U4 U5 U6 1

�α2U1

Numerical Analysis of the Incompressible Fluid Flow and Heat Transfer

http://dx.doi.org/10.5772/intechopen.72263

1

CCCCCCCCCA

1

0

U2 U3 U4 U5 U6 1

CCCCCCCCCA

(11)

CCCCCCCCCA

BBBBBBBBB@

�β6U7

0 0 0

0

BBBBBBBBB@

1

CCCCCCCCCA

CCCCCCCCCA ¼

CCCCCCCCCA

This kind of tridiagonal matrix is often seen and can be solved by a direct numerical method, such as Tomas method. However, the rank of the matrix is usually extremely large and one

In general, the rank of the matrix appearing in computational fluid dynamics (CFD) is large and iterative methods such as Jacobi, Gauss-Seidel, or successive over relaxation (SOR) methodare employed. In this subsection, the Jacobi method is explained. The matrix can be divided into

0

BBBBBBBBB@

�α2U1

�βN�1UN

0 0 0

0

BBBBBBBBB@

In the Jacobi method, only the diagonal part is put in the left-hand side (n + 1 step), while the

BBBBBBBBB@

(9)

247

(10)

(

�2 β2 0 00 α3 � 2 β3 0 0 0 α4 � 2 β4 0 0 0 α5 � 2 β5 000 α6 � 2

introduces an iterative method for solving the king-size matrix.

1

0

U2 U3 U4 U5 U6 1

CCCCCCCCCA þ

CCCCCCCCCA

1

0

U2 U3 U4 U5 U6 1

CCCCCCCCCA ¼

CCCCCCCCCA

lower and upper parts are moved to the right-hand side (n step)

BBBBBBBBB@

BBBBBBBBB@

the following simultaneous equation is obtained:

0

BBBBBBBBB@

2.3. Iterative method for matrix solver

0

BBBBBBBBB@

0

þ

BBBBBBBBB@

three parts of lower, diagonal, and upper as follows:

As illustrated in Figure 2, in which vertical and horizontal axes are exchanged from Figure 1, one needs to obtain each value of dimensionless velocity numerically. The approximated velocity profile is expressed by connecting these values smoothly. For simplicity, the intervals between neighboring two points are the same and it is noted as Δη. When the second-order central difference method is used, Eq. (4) is as follows:

$$\frac{\mathcal{U}\_{i+1} - 2\mathcal{U}\_i + \mathcal{U}\_{i-1}}{\left(\Delta \eta\right)^2} + 2\eta\_i \frac{\mathcal{U}\_{i+1} - \mathcal{U}\_{i-1}}{2\left(\Delta \eta\right)} = 0, \quad (i = 2, 3, \dots \text{ $N - 1$ })$$

Here, N is total number of grids and in this chapter, the first grid point starts from 1 as its definition. The above equation becomes

$$\underbrace{\left\{1-\eta\_{i}(\Delta\eta)\right\}}\_{a\_{i}}\mathcal{U}\_{i-1}-2\mathcal{U}\_{i}+\underbrace{\left\{1+\eta\_{i}(\Delta\eta)\right\}}\_{f\_{i}}\mathcal{U}\_{i+1}=0,\quad(i=2,3,\cdots N-1).\tag{8}$$

Here, ηi ¼ ð Þ i � 1 Δη � � and αi and βi are coefficients determined by the number of grids. The boundary condition (7) is modified

Figure 2. Equidistant grids discretized.

Numerical Analysis of the Incompressible Fluid Flow and Heat Transfer http://dx.doi.org/10.5772/intechopen.72263 247

$$\begin{cases} \eta = 0: & \mathcal{U}\_1 = 1 \\ \eta = 4: & \mathcal{U}\_N = 0 \end{cases} \tag{9}$$

In the following, the case of N = 7 is considered, for example. By substituting i = 2–6 into Eq. (8), the following simultaneous equation is obtained:

$$
\begin{pmatrix} -2 & \beta\_2 & 0 & 0 & 0 \\ \alpha\_3 & -2 & \beta\_3 & 0 & 0 \\ 0 & \alpha\_4 & -2 & \beta\_4 & 0 \\ 0 & 0 & \alpha\_5 & -2 & \beta\_5 \\ 0 & 0 & 0 & \alpha\_6 & -2 \end{pmatrix} \begin{pmatrix} U\_2 \\ U\_3 \\ U\_4 \\ U\_5 \\ U\_6 \end{pmatrix} = \begin{pmatrix} -\alpha\_2 \boldsymbol{\mu} \boldsymbol{I}\_1 \\ 0 \\ 0 \\ 0 \\ 0 \\ -\beta\_6 \boldsymbol{\mu} \boldsymbol{I}\_7 \end{pmatrix} \tag{10}
$$

This kind of tridiagonal matrix is often seen and can be solved by a direct numerical method, such as Tomas method. However, the rank of the matrix is usually extremely large and one introduces an iterative method for solving the king-size matrix.

#### 2.3. Iterative method for matrix solver

2.2. Numerical method for solving the ordinary differential equation using finite

η = 4 is enough. Hence, the boundary condition shown below is used instead of Eq. (5)

�

central difference method is used, Eq. (4) is as follows:

246 Finite Element Method - Simulation, Numerical Analysis and Solution Techniques

2 þ 2ηi

Ui�1 � 2Ui þ 1 þ ηi

Uiþ1 � 2Ui þ Ui�1 ð Þ Δη

definition. The above equation becomes

1 � ηi ð Þ Δη � � |fflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflffl} αi

boundary condition (7) is modified

Figure 2. Equidistant grids discretized.

Here, ηi ¼ ð Þ i � 1 Δη

For numerical solution, it is necessary to define the range of η, as recognized from Figure 1,

η ¼ 0 : U ¼ 1 η ¼ 4 : U ¼ 0

As illustrated in Figure 2, in which vertical and horizontal axes are exchanged from Figure 1, one needs to obtain each value of dimensionless velocity numerically. The approximated velocity profile is expressed by connecting these values smoothly. For simplicity, the intervals between neighboring two points are the same and it is noted as Δη. When the second-order

Uiþ1 � Ui�1

Here, N is total number of grids and in this chapter, the first grid point starts from 1 as its

ð Þ Δη � � |fflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflffl} βi

2ð Þ Δη ¼ 0, ið Þ ¼ 2; 3; ⋯N � 1

� � and αi and βi are coefficients determined by the number of grids. The

Uiþ1 ¼ 0, ið Þ ¼ 2; 3; ⋯N � 1 : (8)

(7)

difference method

In general, the rank of the matrix appearing in computational fluid dynamics (CFD) is large and iterative methods such as Jacobi, Gauss-Seidel, or successive over relaxation (SOR) methodare employed. In this subsection, the Jacobi method is explained. The matrix can be divided into three parts of lower, diagonal, and upper as follows:

$$
\begin{pmatrix} 0 & 0 & 0 & 0 & 0 \\ \alpha\_3 & 0 & 0 & 0 & 0 \\ 0 & \alpha\_4 & 0 & 0 & 0 \\ 0 & 0 & \alpha\_5 & 0 & 0 \\ 0 & 0 & 0 & \alpha\_6 & 0 \end{pmatrix} \begin{pmatrix} U\_2 \\ U\_3 \\ U\_4 \\ U\_5 \\ U\_6 \\ U\_6 \end{pmatrix} + \begin{pmatrix} -2 & 0 & 0 & 0 & 0 \\ 0 & -2 & 0 & 0 & 0 \\ 0 & 0 & -2 & 0 & 0 \\ 0 & 0 & 0 & -2 & 0 \\ 0 & 0 & 0 & 0 & -2 \end{pmatrix} \begin{pmatrix} U\_2 \\ U\_3 \\ U\_4 \\ U\_5 \\ U\_6 \end{pmatrix} \\
$$

$$
+ \begin{pmatrix} 0 & \beta\_2 & 0 & 0 & 0 \\ 0 & 0 & \beta\_3 & 0 & 0 \\ 0 & 0 & 0 & \beta\_4 & 0 \\ 0 & 0 & 0 & 0 & \beta\_5 \\ 0 & 0 & 0 & 0 & 0 \end{pmatrix} \begin{pmatrix} U\_2 \\ U\_3 \\ U\_4 \\ U\_5 \\ U\_6 \\ U\_7 \end{pmatrix} = \begin{pmatrix} -a\_2 \, \boldsymbol{I}\_1 \\ 0 \\ 0 \\ 0 \\ 0 \\ -\beta\_{N-1} \, \boldsymbol{I}\_N \end{pmatrix} \tag{11}
$$

In the Jacobi method, only the diagonal part is put in the left-hand side (n + 1 step), while the lower and upper parts are moved to the right-hand side (n step)

�20 0 0 0 0 � 20 00 0 0 � 200 000 � 2 0 0000 � 2 0 BBBBBB@ 1 CCCCCCA |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} D U2 U3 U4 U5 U6 0 BBBBBB@ 1 CCCCCCA nþ1 |fflfflfflfflfflfflffl{zfflfflfflfflfflfflffl} U ! nþ1 ¼ �α2U1 0 0 0 �βN�1UN 0 BBBBBB@ 1 CCCCCCA |fflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflffl} b ! � 00000 α3 00 0 0 0 α4 00 0 0 0 α5 0 0 000 α6 0 0 BBBBBB@ 1 CCCCCCA |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} L U2 U3 U4 U5 U6 0 BBBBBB@ 1 CCCCCCA n |fflfflfflfflffl{zfflfflfflfflffl} U ! n � 0 β2 000 0 0 β3 0 0 000 β4 0 0000 β5 00000 0 BBBBBB@ 1 CCCCCCA |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} U U2 U3 U4 U5 U6 0 BBBBBB@ 1 CCCCCCA n |fflfflfflfflffl{zfflfflfflfflffl} U ! n (12)

Here, n is the old iteration step and n + 1 is the new iteration step. Hence, the following equation is repeatedly used:

$$
\overrightarrow{\boldsymbol{U}}^{n+1} = \boldsymbol{D}^{-1} \left[ \overrightarrow{\boldsymbol{b}} - (\boldsymbol{L} + \boldsymbol{U}) \overrightarrow{\boldsymbol{U}}^{n} \right] \tag{13}
$$

shown in Eqs. (15)–(17) together with the boundary condition (18). Here, one defines that x axis is in the vertical direction and its velocity component is u, and y axis is in the direction

> ∂2 u

y ¼ 0 : u ¼ v ¼ 0, q ¼ �kð Þ ∂T=∂y

Here, β is the thermal expansion coefficient, g is the acceleration due to gravity, α is the thermal

First, dimensionless variables, such as velocity and temperature, are set as follows using the

vaxa yaua |ffl{zffl} ½ � 1

> ∂2 U

∂V ∂Y ¼ 0

, V ¼ v va

∂Y2 þ gβð Þ Tb � T∞ xa ua 2 |fflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflffl} ½ � 3

, θ ¼ T � Tb Ta

> þ gβTaxa ua 2 |fflfflffl{zfflfflffl} ½ � 4

θ

, U ¼ u ua

∂U ∂X þ

y ! ∞ : u ! 0, T ! T∞

∂2 T

∂y ¼ 0 (15)

Numerical Analysis of the Incompressible Fluid Flow and Heat Transfer

http://dx.doi.org/10.5772/intechopen.72263

∂y2 þ gβð Þ T � T∞ (16)

∂y2 (17)

(18)

249

(19)

∂u ∂x þ ∂v

perpendicular to the vertical plate and its velocity component is v.

u ∂u ∂x þ v ∂u ∂y ¼ ν

�

unknown reference value denoted with subscripts a and b:

Equation (19) is substituted into Eqs. (15)–(18), and one gets

vaxa yaua |ffl{zffl} ½ � 1

V ∂U ∂Y ¼ νxa ya 2ua |ffl{zffl} ½ � 2

, Y ¼ y ya

X ¼ x xa

U ∂U ∂X þ

diffusivity, k is the thermal conductivity, and T is the temperature.

u ∂T ∂x þ v ∂T ∂y ¼ α

Continuity of mass

Momentum equation

Energy equation

Boundary equation

3.3. Non-dimensionalization

This is equivalent to the following equation:

$$\mathcal{U}\mathcal{U}\_{i}^{n+1} = \frac{1}{2} (\alpha\_i \mathcal{U}\_{i-1}^n + \beta\_i \mathcal{U}\_{i+1}^n), \quad (i = 2, 3, 4, 5, 6) \tag{14}$$

By using Eq. (9), Eq. (14) is computed repeatedly and then the value of each grid gradually converges to a certain solution. The Gauss-Seidel and SOR methods are known as the faster convergence method.

### 3. Similarity solution for natural convection heated from a vertical plate

#### 3.1. Introduction

In this section, let us consider the natural convection heat transfer for a vertical plate heated with uniform heat flux in the wide range of Prandtl number from zero to infinity. In order to explain the numerical method as how to solve the governing equations, one assumes that the flow and temperature fields formed in the vicinity of the heated plate have a similarity and then one introduces the finite difference method to obtain numerical results.

#### 3.2. Governing equations

One assumes that the flow is incompressible laminar and boundary layer equations are used in this analysis. The governing equations with presuming the Boussinesq approximation are shown in Eqs. (15)–(17) together with the boundary condition (18). Here, one defines that x axis is in the vertical direction and its velocity component is u, and y axis is in the direction perpendicular to the vertical plate and its velocity component is v.

Continuity of mass

�20 0 0 0 0 � 20 00 0 0 � 200 000 � 2 0 0000 � 2

0 β2 000 0 0 β3 0 0 000 β4 0 0000 β5 00000

equation is repeatedly used:

convergence method.

3.1. Introduction

3.2. Governing equations

1

0

248 Finite Element Method - Simulation, Numerical Analysis and Solution Techniques

BBBBBB@

U2 U3 U4 U5 U6

U2 U3 U4 U5 U6 1

nþ1

¼

0

BBBBBB@

1

CCCCCCA

�

0

BBBBBB@

00000 α3 00 0 0 0 α4 00 0 0 0 α5 0 0 000 α6 0

1

0

BBBBBB@

U2 U3 U4 U5 U6

1

n

CCCCCCA

(13)

CCCCCCA

�βN�1UN

Here, n is the old iteration step and n + 1 is the new iteration step. Hence, the following

Un iþ1

By using Eq. (9), Eq. (14) is computed repeatedly and then the value of each grid gradually converges to a certain solution. The Gauss-Seidel and SOR methods are known as the faster

3. Similarity solution for natural convection heated from a vertical plate

In this section, let us consider the natural convection heat transfer for a vertical plate heated with uniform heat flux in the wide range of Prandtl number from zero to infinity. In order to explain the numerical method as how to solve the governing equations, one assumes that the flow and temperature fields formed in the vicinity of the heated plate have a similarity and

One assumes that the flow is incompressible laminar and boundary layer equations are used in this analysis. The governing equations with presuming the Boussinesq approximation are

�ð Þ L þ U U !n � �

� �, ið Þ ¼ 2; 3; 4; 5; 6 (14)

¼ D�1 b !

i�1 þ βi

CCCCCCA

n

1

CCCCCCA

U !nþ1

2 αiUn

then one introduces the finite difference method to obtain numerical results.

(12)

CCCCCCA

1

0

BBBBBB@

CCCCCCA

This is equivalent to the following equation:

Unþ1 i ¼ 1

0

BBBBBB@

0

BBBBBB@

�

$$\frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} = 0 \tag{15}$$

Momentum equation

$$
\mu \frac{\partial u}{\partial x} + v \frac{\partial u}{\partial y} = \nu \frac{\partial^2 u}{\partial y^2} + \text{g}\beta (T - T\_{\text{\textquotedblleft}}) \tag{16}
$$

Energy equation

$$
\mu \frac{\partial T}{\partial x} + v \frac{\partial T}{\partial y} = \alpha \frac{\partial^2 T}{\partial y^2} \tag{17}
$$

Boundary equation

$$\begin{cases} y = 0: \quad u = v = 0, \quad q = -k(\partial T/\partial y) \\ y \to \rightsquigarrow \quad u \to 0, \quad T \to T\_{\rightsquigarrow} \end{cases} \tag{18}$$

Here, β is the thermal expansion coefficient, g is the acceleration due to gravity, α is the thermal diffusivity, k is the thermal conductivity, and T is the temperature.

#### 3.3. Non-dimensionalization

First, dimensionless variables, such as velocity and temperature, are set as follows using the unknown reference value denoted with subscripts a and b:

$$X = \frac{x}{x\_a}, \quad Y = \frac{y}{y\_a}, \quad \mathcal{U} = \frac{u}{u\_a}, \quad V = \frac{v}{v\_a}, \quad \Theta = \frac{T - T\_b}{T\_a} \tag{19}$$

Equation (19) is substituted into Eqs. (15)–(18), and one gets

$$\frac{\partial \mathcal{U}}{\partial \mathcal{X}} + \underbrace{\frac{\upsilon\_a \mathbf{x}\_a}{\underbrace{y\_a u\_a}} \frac{\partial \mathcal{V}}{\partial \mathcal{Y}}}\_{[1]} = \mathbf{0}$$

$$\mathcal{U}\frac{\partial \mathcal{U}}{\partial \mathcal{X}} + \underbrace{\frac{\upsilon\_a \mathbf{x}\_a}{y\_a u\_a} \mathcal{V}\frac{\partial \mathcal{U}}{\partial \mathcal{Y}}}\_{[1]} = \underbrace{\frac{\nu \mathbf{x}\_a}{y\_a^2 u\_a} \frac{\partial^2 \mathcal{U}}{\partial \mathcal{Y}^2}}\_{[2]} + \underbrace{\frac{g\beta (T\_b - T\_\Leftrightarrow) \mathbf{x}\_a}{u\_a^2}}\_{[3]} + \underbrace{\frac{g\beta T\_a \mathbf{x}\_a}{u\_a^2} \Theta}\_{[4]}$$

$$\begin{aligned} \underbrace{\mathcal{U}\frac{\partial\Theta}{\partial\mathcal{X}} + \underbrace{\frac{v\_{a}x\_{a}}{y\_{a}u\_{a}}}\_{[1]}}\_{[1]}\underbrace{\delta\Theta}\_{\partial Y} &= \underbrace{\frac{\alpha x\_{a}}{y\_{a}^{2}u\_{a}}}\_{[5]}\frac{\partial^{2}\Theta}{\partial Y^{2}} \\\\ \left\{ \begin{aligned} Y &= 0 &: \quad \mathcal{U} = V = 0, \quad & \underbrace{(qy\_{a})/(kT\_{a})}\_{[6]} = -\partial\theta/\partial Y \\\\ Y &\to \infty &: \quad \mathcal{U} = 0, \quad & \theta = \underbrace{(T\_{\infty} - T\_{b})/T\_{a}}\_{[7]} \end{aligned} \right. \end{aligned}$$

Y ¼ 0 : U ¼ V ¼ 0, ∂θ=∂Y ¼ �1

xa Ra∗ ð Þ Pr �1=5 , U ¼ u

k Ra∗ ð Þ Pr �1=5

qxa

Furthermore, one assumes that the velocity and temperature fields has a similarity along the direction of vertical plate, so one puts X = 1. These equations are useful for analyzing low

> þ Pr d2 U

� V dθ dη þ d2 θ

η ¼ 0 : U ¼ V ¼ 0, dθ=dη ¼ �1

∗ ¼ gβqx4

ανk , Pr ¼ ν

, V ¼ v α x Rax ∗ ð Þ Pr 1 5 ,

α

α xa Ra∗ ð Þ Pr 2=5 ,

Numerical Analysis of the Incompressible Fluid Flow and Heat Transfer

http://dx.doi.org/10.5772/intechopen.72263

(29)

dη2 þ θ ¼ 0 (30)

dη2 ¼ 0 (31)

(27)

251

(28)

(32)

(33)

Y ! ∞ : U ¼ 0, θ ¼ 0

Ra∗ ð Þ Pr 1=5 , θ ¼ T � T∞

dV dη ¼ 1 5 η dU dη � 3U

dη � V dU dη

η ! ∞ : U ¼ θ ¼ 0

The dimensionless variables and non-dimensional numbers are defined as follows:

, U ¼ u α x Rax ∗ ð Þ Pr 2 5

U dV

U 5 η dθ dη � θ

5

5 , Rax

The local Nusselt number can be obtained by the following derivation:

η ¼ y x Rax ∗ ð Þ Pr �1

θ ¼ T � T∞ qx k Rax ∗ ð Þ Pr �1

The dimensionless variables are summarized as follows:

, Y ¼ y

X ¼ x xa

V ¼ v α xa

Prandtl number cases and summarized as follows:

Low Prandtl number

Continuity of mass

Momentum equation

Energy equation

Boundary conditions

At the moment stage, xa is recognized as the height of the vertical plate. Putting [3] = 0, and one obtains Tb ¼ T∞.Hence [7] becomes θ ¼ 0.

Putting [6] = 1, and one gets qya kTa ¼ 1 ) Ta ¼ qya k Putting [5] = 1, and one gets αxa ya 2ua ¼ 1 ) ya ¼ αxa ua � �1=2 Putting [1] = 1, vaxa yaua ¼ 1 ) va ¼ yaua xa

$$\text{Putting } [4] = 1, \frac{\mathfrak{g}\mathcal{T}\_{\mathfrak{s}}\mathfrak{x}\_{\mathfrak{s}}}{u\_{\mathfrak{s}}^{2}} = 1 \quad \Rightarrow \quad \mathfrak{u}\_{\mathfrak{a}} = \left(\mathfrak{g}\mathcal{T}\_{\mathfrak{a}}\mathbf{x}\_{\mathfrak{a}}\right)^{1/2} = \left(\mathfrak{g}\mathcal{T}\_{\mathfrak{k}}^{\mathfrak{g}\mathfrak{z}}\mathbf{x}\_{\mathfrak{a}}\right)^{1/2} = \left(\mathfrak{g}\mathcal{T}\_{\mathfrak{k}}^{\mathfrak{g}}\left(\frac{\alpha\mathbf{x}\_{\mathfrak{a}}}{u\_{\mathfrak{a}}}\right)^{1/2}\mathbf{x}\_{\mathfrak{a}}\right)^{1/2}$$

$$\mathbf{u} \cdot \mathbf{u}\_d = \left\{ \left( \mathbf{g} \beta \frac{\eta}{k} \right)^2 \alpha \mathbf{x}\_d^{\ 3} \right\}^{1/5} = \left( \frac{\mathbf{g} \beta \mathbf{q} \mathbf{x}\_d^4}{k a^2} \right)^{2/5} \frac{\alpha}{\mathbf{x}\_d} = (\mathbf{R} \mathbf{a}^\* \mathbf{P} \mathbf{r})^{2/5} \frac{\alpha}{\mathbf{x}\_d}, \quad \because \mathbf{R} \mathbf{a}^\* = \frac{\mathbf{g} \beta \mathbf{q} \mathbf{x}\_d^{\ 4}}{k \alpha \mathbf{v}} \tag{20}$$

$$\mathbf{y}\_d = \left(\frac{a\mathbf{x}\_d}{\mathbf{u}\_d}\right)^{1/2} = \left(\frac{a\mathbf{x}\_d}{\left(Ra^\*Pr\right)^{2/5}\frac{a}{\mathbf{x}\_d}}\right)^{1/2} = \left(\frac{\mathbf{x}\_d^{\;2}}{\left(Ra^\*Pr\right)^{2/5}}\right)^{1/2} = \mathbf{x}\_d \left(Ra^\*Pr\right)^{-1/5} \tag{21}$$

$$T\_a = \frac{qy\_a}{k} = \frac{q}{k} \frac{x\_a}{(Ra^\*Pr)^{1/5}} = \frac{q\chi\_a}{k} \left(Ra^\*Pr\right)^{-1/5} \tag{22}$$

$$w\_a = \frac{y\_a u\_a}{\chi\_a} = \frac{\frac{\chi\_a}{(Ra^\*Pr)^{1/5}} (Ra^\*Pr)^{2/5} \frac{\alpha}{\chi\_a}}{\chi\_a} = (Ra^\*Pr)^{1/5} \frac{\alpha}{\chi\_a} \tag{23}$$

In the above process, finally one obtains the dimensionless equations as follows:

$$
\frac{\partial \mathcal{U}}{\partial X} + \frac{\partial V}{\partial Y} = 0 \tag{24}
$$

$$
\Delta U \frac{\partial \mathcal{U}}{\partial X} + V \frac{\partial \mathcal{U}}{\partial Y} = Pr \frac{\partial^2 \mathcal{U}}{\partial Y^2} + \theta \tag{25}
$$

$$U\frac{\partial\theta}{\partial X} + V\frac{\partial\theta}{\partial Y} = \frac{\partial^2\theta}{\partial Y^2} \tag{26}$$

Numerical Analysis of the Incompressible Fluid Flow and Heat Transfer http://dx.doi.org/10.5772/intechopen.72263 251

$$\begin{cases} Y = 0: \quad \mathcal{U} = V = 0, \quad \partial\theta/\partial Y = -1\\ Y \to \rightsquigarrow \quad \mathcal{U} = 0, \quad \theta = 0 \end{cases} \tag{27}$$

The dimensionless variables are summarized as follows:

$$\begin{aligned} X &= \frac{\mathbf{x}}{\mathbf{x}\_d}, \quad Y = \frac{y}{\mathbf{x}\_d (Ra^\*Pr)^{-1/5}}, \quad \mathcal{U} = \frac{u}{\frac{\alpha}{\mathbf{x}\_d} (Ra^\*Pr)^{2/5}},\\ V &= \frac{v}{\frac{\alpha}{\mathbf{x}\_d} (Ra^\*Pr)^{1/5}}, \quad \Theta = \frac{T - T\_{\infty}}{\frac{q\mathbf{x}\_d}{k} (Ra^\*Pr)^{-1/5}} \end{aligned} \tag{28}$$

Furthermore, one assumes that the velocity and temperature fields has a similarity along the direction of vertical plate, so one puts X = 1. These equations are useful for analyzing low Prandtl number cases and summarized as follows:

#### Low Prandtl number

Continuity of mass

U ∂θ ∂X þ

250 Finite Element Method - Simulation, Numerical Analysis and Solution Techniques

8 >>>><

>>>>:

Putting [6] = 1, and one gets qya

Putting [5] = 1, and one gets αxa

ua

∴ ua ¼ gβ

ya ¼ αxa ua � �1=2

Putting [1] = 1, vaxa

Putting [4] = 1, gβTaxa

vaxa yaua |ffl{zffl} ½ � 1

Y ! ∞ : U ¼ 0, θ ¼ ð Þ T∞ � Tb =Ta |fflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflffl}

Y ¼ 0 : U ¼ V ¼ 0, qya

At the moment stage, xa is recognized as the height of the vertical plate.

kTa ¼ 1 ) Ta ¼ qya

xa

2ua ¼ 1 ) ya ¼ αxa

¼ gβqxa

xa

xa Ra∗ ð Þ Pr 1=5 ¼ qxa

Ra∗ ð Þ Pr 1=5 Ra∗ ð Þ Pr 2=5 α

xa

∂V

∂Y ¼ Pr ∂2

∂Y ¼ ∂2

Putting [3] = 0, and one obtains Tb ¼ T∞.Hence [7] becomes θ ¼ 0.

ya

2 ¼ 1 ) ua ¼ gβTaxa

¼ αxa

Ta ¼ qya

va ¼ yaua xa ¼

Ra∗ ð Þ Pr 2=5 α

k ¼ q k

xa

In the above process, finally one obtains the dimensionless equations as follows:

U ∂U

∂U ∂X þ

∂X þ V ∂U

U ∂θ

∂X þ V ∂θ

!1=2

yaua ¼ 1 ) va ¼ yaua

αxa 3 � �1=5

q k � �2 V ∂θ ∂Y ¼ αxa ya 2ua |ffl{zffl} ½ � 5

∂2 θ ∂Y2

¼ �∂θ=∂Y

� �=ð Þ kTa |fflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflffl} ½ � 6

½ � 7

k

� �1=2 ¼ gβ qya

4 kα2 � �2=5

ua � �1=2

> α xa

¼ xa

xa

U

θ

k xa � �1=2 ¼ gβ q

2 Ra∗ ð Þ Pr 2=5

!1=2

¼ Ra∗ ð Þ Pr 2=5 α

¼ Ra∗ ð Þ Pr 1=5 α

k αxa ua � �1=2

xa

xa

¼ xa Ra∗ ð Þ Pr �1=5 (21)

4 kαν (20)

(23)

� �1=2

, ∵Ra∗ ¼ gβqxa

k Ra∗ ð Þ Pr �1=5 (22)

xa

∂Y ¼ 0 (24)

∂Y2 þ θ (25)

∂Y2 (26)

$$\frac{dV}{d\eta} = \frac{1}{5} \left( \eta \frac{d\mathcal{U}}{d\eta} - \mathfrak{G}\mathcal{U} \right) \tag{29}$$

Momentum equation

$$
\mathcal{U}\frac{dV}{d\eta} - V\frac{d\mathcal{U}}{d\eta} + \mathcal{P}r\frac{d^2\mathcal{U}}{d\eta^2} + \theta = 0\tag{30}
$$

Energy equation

$$\frac{dI}{5}\left(\eta \frac{d\theta}{d\eta} - \theta\right) - V\frac{d\theta}{d\eta} + \frac{d^2\theta}{d\eta^2} = 0\tag{31}$$

Boundary conditions

$$\begin{cases} \eta = 0: \quad \mathcal{U} = V = 0, \quad d\theta/d\eta = -1\\ \eta \to \stackrel{\sim}{\sim}: \quad \mathcal{U} = \theta = 0 \end{cases} \tag{32}$$

The dimensionless variables and non-dimensional numbers are defined as follows:

$$\begin{aligned} \eta &= \frac{y}{\text{x}(Ra\_{\text{x}}^{\*}Pr)^{-\frac{1}{5}}}, \quad U = \frac{u}{\frac{\alpha}{\text{x}}(Ra\_{\text{x}}^{\*}Pr)^{\frac{2}{5}}}, \quad V = \frac{v}{\frac{\alpha}{\text{x}}(Ra\_{\text{x}}^{\*}Pr)^{\frac{1}{5}}},\\ \Theta &= \frac{T - T\_{\text{w}}}{\frac{qX}{k}(Ra\_{\text{x}}^{\*}Pr)^{-\frac{1}{5}}}, \quad Ra\_{\text{x}}^{\*} = \frac{g\beta qx^{4}}{\alpha vk}, \quad Pr = \frac{\nu}{\alpha} \end{aligned} \tag{33}$$

The local Nusselt number can be obtained by the following derivation:

$$Nu\_x = \frac{h\_x x}{k} = \frac{qx}{(T\_w - T\_\circ)k} = \frac{qx}{T\_a \theta\_w k} = \underbrace{\frac{\left(\frac{g\beta q}{\alpha^2 k}\right)^{1/5}}{qx^{1/5}}}\_{T\_x - 1} \frac{k}{\theta\_w k} \tag{34}$$

$$= \left(\frac{g\beta q}{\alpha^2 k}\right)^{1/5} \frac{x^{4/5}}{\theta\_w} = \left(\frac{g\beta qx^4}{\alpha^2 k}\right)^{1/5} \frac{1}{\theta\_w} = (Ra\_x^\* Pr)^{1/5} \frac{1}{\theta\_w}$$

$$Nu\_x = (Ra\_x^\* Pr)^{\frac{1}{5}} (\theta|\_{\eta=0})^{-1} = (Ra\_x Nu\_x Pr)^{\frac{1}{5}} (\theta|\_{\eta=0})^{-1} \tag{35}$$

Therefore, the local Nusselt number can be obtained just from the dimensionless temperature at the wall using Eq. (36)

$$\frac{N\mu\_{\chi}}{(R a\_{\mathfrak{x}} Pr)^{\frac{1}{\mathfrak{x}}}} = \left(\left.\Theta\right|\_{\mathfrak{y}=0}\right)^{-\frac{5}{4}}\tag{36}$$

Nux ð Þ Rax 1 4 ¼ θj η¼0 �5 4

the cases of Pr ≥ 1, while the right-hand side ones the cases of Pr ≤ 1

of discrepancy exists. In this study, the boundary condition for Pr ! 0

η ! ∞ : U ¼ θ ¼ 0

Figure 3 shows the numerical result for the various Prandtl number cases. The upper figures indicate the vertical velocity and lower ones the temperature. The left-hand side figures show

Numerical Analysis of the Incompressible Fluid Flow and Heat Transfer

http://dx.doi.org/10.5772/intechopen.72263

Table 1 shows the summary of the local Nusselt number for various Prandtl number cases together with the reference of Churchill and Ozoe for comparison [4]. The agreement is quite good except for the extreme cases such as Pr ! 0 and ∞. In such extreme cases, a small amount

η ¼ 0 : dU=dη ¼ V ¼ 0, dθ=dη ¼ �1

Figure 3. Vertical velocity and temperature distributions for various Prandtl numbers. The left-hand side indicates high

Prandtl number cases while the right-hand side low Prandtl number cases.

3.4. Numerical results

and that for Pr ! ∞

(42)

253

#### High Prandtl number

If the Prandtl number is higher than unity, the following equations are useful:

Continuity of mass

$$\frac{dV}{d\eta} = \frac{1}{5} \left( \eta \frac{d\mathcal{U}}{d\eta} - 3\mathcal{U} \right) \tag{37}$$

Momentum equation

$$
\rho \mathcal{U} \frac{dV}{d\eta} - V \frac{d\mathcal{U}}{d\eta} + Pr \frac{d^2 \mathcal{U}}{d\eta^2} + Pr \theta = 0 \tag{38}
$$

Energy equation

$$\frac{dI}{5}\left(\eta\frac{d\theta}{d\eta}-\theta\right)-V\frac{d\theta}{d\eta}+\frac{d^2\theta}{d\eta^2}=0\tag{39}$$

Boundary conditions

$$\begin{cases} \eta = 0: \quad \mathcal{U} = V = 0, \quad d\theta/d\eta = -1\\ \eta \to \stackrel{\sim}{\sim}: \quad \mathcal{U} = \theta = 0 \end{cases} \tag{40}$$

The dimensionless variables and non-dimensional numbers are defined as follows:

$$\begin{aligned} \eta &= \frac{y}{\mathfrak{x}(Ra\_{\mathfrak{x}}^{\*})^{-\frac{1}{5}}}, \quad \mathcal{U} = \frac{u}{\frac{\alpha}{\mathfrak{x}}(Ra\_{\mathfrak{x}}^{\*})^{\frac{2}{5}}}, \quad \mathcal{V} = \frac{v}{\frac{\alpha}{\mathfrak{x}}(Ra\_{\mathfrak{x}}^{\*})^{\frac{1}{5}}},\\ \Theta &= \frac{T - T\_{\infty}}{\frac{q\chi}{k}(Ra\_{\mathfrak{x}}^{\*})^{-\frac{1}{5}}}, \quad Ra\_{\mathfrak{x}}^{\*} = \frac{g\beta q \mathfrak{x}^{4}}{\alpha vk}, \quad Pr = \frac{\nu}{\alpha} \end{aligned} \tag{41}$$

Numerical Analysis of the Incompressible Fluid Flow and Heat Transfer http://dx.doi.org/10.5772/intechopen.72263 253

$$\frac{N\mu\_{\rm x}}{(Ra\_{x})^{\frac{1}{4}}} = \left(\theta|\_{\eta=0}\right)^{-\frac{1}{4}}\tag{42}$$

#### 3.4. Numerical results

Nux ¼ hxx

252 Finite Element Method - Simulation, Numerical Analysis and Solution Techniques

¼ gβq α2k � �1=5 x4=5

Nux ¼ Rax

at the wall using Eq. (36)

High Prandtl number

Continuity of mass

Momentum equation

Energy equation

Boundary conditions

k ¼ qx

θw

5 θj η¼0 � ��1

Nux ð Þ RaxPr 1 4 ¼ θj η¼0 � ��5 4

If the Prandtl number is higher than unity, the following equations are useful:

dV dη ¼ 1 5 η dU dη � 3U � �

dη � V dU dη

þ Pr d2 U

> � V dθ dη þ d2 θ

η ¼ 0 : U ¼ V ¼ 0, dθ=dη ¼ �1

∗ ¼ gβqx4

, V ¼ v α x Rax ∗ ð Þ1 5 ,

α

ανk , Pr ¼ ν

η ! ∞ : U ¼ θ ¼ 0

The dimensionless variables and non-dimensional numbers are defined as follows:

, U ¼ u α x Rax ∗ ð Þ2 5

U dV

U 5 η dθ dη � θ � �

�

η ¼ y x Rax ∗ ð Þ�1 5

θ ¼ T � T∞ qx k Rax ∗ ð Þ�1 5 , Rax

∗ ð Þ Pr 1

ð Þ Tw � T∞ k ¼ qx

¼ gβqx4 α2k � �1=5 1

Therefore, the local Nusselt number can be obtained just from the dimensionless temperature

Taθwk ¼

θw

¼ ð Þ RaxNuxPr 1

gβq α2k � �1=5 k

¼ Ra∗

qx1=5 |fflfflfflfflffl{zfflfflfflfflffl} Ta �1

xPr � �1=5 1

5 θj η¼0 � ��1

qx θwk

θw

dη2 þ Prθ ¼ 0 (38)

dη2 ¼ 0 (39)

(34)

(35)

(36)

(37)

(40)

(41)

Figure 3 shows the numerical result for the various Prandtl number cases. The upper figures indicate the vertical velocity and lower ones the temperature. The left-hand side figures show the cases of Pr ≥ 1, while the right-hand side ones the cases of Pr ≤ 1

Table 1 shows the summary of the local Nusselt number for various Prandtl number cases together with the reference of Churchill and Ozoe for comparison [4]. The agreement is quite good except for the extreme cases such as Pr ! 0 and ∞. In such extreme cases, a small amount of discrepancy exists. In this study, the boundary condition for Pr ! 0

$$\begin{cases} \eta = 0: \quad d\mathcal{U}/d\eta = V = 0, \quad d\theta/d\eta = -1\\ \eta \to \stackrel{\sim}{\sim}: \quad \mathcal{U} = \theta = 0 \end{cases}$$

and that for Pr ! ∞

Figure 3. Vertical velocity and temperature distributions for various Prandtl numbers. The left-hand side indicates high Prandtl number cases while the right-hand side low Prandtl number cases.

Table 1. Local Nusselt number for various values of Prandtl number (the upper: present results, the lower: Churchill and Ozoe [4]).

$$\begin{cases} \eta = 0: \quad \mathcal{U} = V = 0, \quad d\theta/d\eta = -1\\ \eta \to \stackrel{\sim}{\sim}: \quad d\mathcal{U}/d\eta = \theta = 0. \end{cases}$$

are used. Owing to this kind of special treatments for the boundary condition of such extreme cases, one can obtain accurate numerical results for the system of ordinary equations. The results between the solution of the present method and that of Le Fevre [5] for the case of constant temperature of heated wall are identical to each other. The value for Pr ! ∞ is 0.5027 and that for Pr = 0 is 0.6004.

### 4. Linear stability of Taylor-Couette flow

#### 4.1. Governing equations

In the text book of Chandrasekar [6], various examples of the linear stability analysis such as the Rayleigh-Bénard convection, the Taylor-Couette flow, and the Rayleigh-Taylor instability were studied extensively. More recently, Koschmieder [7] described the research focusing on the Bénard cells and the Taylor vortices. In this section, only the Taylor-Coette flow is considered. Figure 4 shows the schematic model considered for the Taylor-Couette flow. In this section, the fluid flow inside of the co-axial double cylindrical enclosure is assumed to be incompressible Newtonian, isothermal and axisymmetric. The gray part represents the computational domain. It is known that the stationary secondary flow is generated at a certain condition under the influence of centrifugal force due to the rotation of primary basic flow which is in azimuthal direction. The continuity of mass and momentum equations are shown in the cylindrical coordinate system as follows:

$$\frac{\partial u\_r}{\partial r} + \frac{u\_r}{r} + \frac{\partial u\_z}{\partial z} = 0 \tag{43}$$

∂uz ∂t þ ur

Figure 4. Schematic model for the Taylor-Couette flow.

4.2. Basic state and linearization

Azimuthal velocity

Pressure

follows:

∂uz ∂r þ uz ∂uz ∂z ¼ � 1 r ∂p ∂z þ ν

uθð Þ¼� r

p rð Þ¼ ; z

r2 2Ω2 � r2

r2 1 � r2 2 ∂2 uz ∂r2 þ

Here, it is indicated that r is the radial, θ is the azimuthal, and z is the axial components.

the basic states for the azimuthal component of velocity and pressure are as follows:

The cylindrical enclosure is long enough to neglect the top and bottom ends. In that situation,

1Ω1

ð rf g uθð Þr

r

Here, Ω1 is the angular velocity at the inner cylinder, Ω2 is the angular velocity at the outer cylinder, p is the pressure, r is the density, and g is the acceleration due to gravity. In order to derive disturbance equations for the linear stability, the three components of velocity and pressure are represented as a summation of basic state and infinitesimal disturbance as

r þ r2 1r2

2

1 r ∂uz ∂r þ ∂2 uz ∂z2

� �

Numerical Analysis of the Incompressible Fluid Flow and Heat Transfer

http://dx.doi.org/10.5772/intechopen.72263

255

2ð Þ Ω2 � Ω1 r2 1 � r2 2

1

dr � rgz þ p0 (48)

� g (46)

r (47)

$$\frac{\partial u\_r}{\partial t} + u\_r \frac{\partial u\_r}{\partial r} + u\_z \frac{\partial u\_r}{\partial z} - \frac{u\_\theta}{r} = -\frac{1}{\rho} \frac{\partial p}{\partial r} + \nu \left( \frac{\partial^2 u\_r}{\partial r^2} + \frac{1}{r} \frac{\partial u\_r}{\partial r} - \frac{u\_r}{r^2} + \frac{\partial^2 u\_r}{\partial z^2} \right) \tag{44}$$

$$\frac{\partial u\_{\theta}}{\partial t} + u\_r \frac{\partial u\_{\theta}}{\partial r} + u\_z \frac{\partial u\_{\theta}}{\partial z} + \frac{u\_r u\_{\theta}}{r} = \nu \left( \frac{\partial^2 u\_{\theta}}{\partial r^2} + \frac{1}{r} \frac{\partial u\_{\theta}}{\partial r} - \frac{u\_{\theta}}{r^2} + \frac{\partial^2 u\_{\theta}}{\partial z^2} \right) \tag{45}$$

Numerical Analysis of the Incompressible Fluid Flow and Heat Transfer http://dx.doi.org/10.5772/intechopen.72263 255

Figure 4. Schematic model for the Taylor-Couette flow.

$$\frac{\partial u\_z}{\partial t} + u\_r \frac{\partial u\_z}{\partial r} + u\_z \frac{\partial u\_z}{\partial z} = -\frac{1}{\rho} \frac{\partial p}{\partial z} + \nu \left( \frac{\partial^2 u\_z}{\partial r^2} + \frac{1}{r} \frac{\partial u\_z}{\partial r} + \frac{\partial^2 u\_z}{\partial z^2} \right) - g \tag{46}$$

Here, it is indicated that r is the radial, θ is the azimuthal, and z is the axial components.

#### 4.2. Basic state and linearization

The cylindrical enclosure is long enough to neglect the top and bottom ends. In that situation, the basic states for the azimuthal component of velocity and pressure are as follows:

Azimuthal velocity

$$\overline{u}\_{\theta}(r) = -\frac{r\_2^2 \Omega\_2 - r\_1^2 \Omega\_1}{r\_1^2 - r\_2^2} r + \frac{r\_1^2 r\_2^2 (\Omega\_2 - \Omega\_1)}{r\_1^2 - r\_2^2} \frac{1}{r} \tag{47}$$

Pressure

η ¼ 0 : U ¼ V ¼ 0, dθ=dη ¼ �1

Table 1. Local Nusselt number for various values of Prandtl number (the upper: present results, the lower: Churchill and

0.456

0.4564 0.456

0.5234 0.524

N/A

0.5495 0.550

0.5631 0.5627

are used. Owing to this kind of special treatments for the boundary condition of such extreme cases, one can obtain accurate numerical results for the system of ordinary equations. The results between the solution of the present method and that of Le Fevre [5] for the case of constant temperature of heated wall are identical to each other. The value for Pr ! ∞ is 0.5027

In the text book of Chandrasekar [6], various examples of the linear stability analysis such as the Rayleigh-Bénard convection, the Taylor-Couette flow, and the Rayleigh-Taylor instability were studied extensively. More recently, Koschmieder [7] described the research focusing on the Bénard cells and the Taylor vortices. In this section, only the Taylor-Coette flow is considered. Figure 4 shows the schematic model considered for the Taylor-Couette flow. In this section, the fluid flow inside of the co-axial double cylindrical enclosure is assumed to be incompressible Newtonian, isothermal and axisymmetric. The gray part represents the computational domain. It is known that the stationary secondary flow is generated at a certain condition under the influence of centrifugal force due to the rotation of primary basic flow which is in azimuthal direction. The continuity of mass and momentum equations are shown in the

∂z ¼ 0 (43)

(44)

(45)

∂ur ∂r þ ur r þ ∂uz

2 r ¼ � 1 r ∂p ∂r þ ν ∂2 ur ∂r2 þ

> 1 r ∂uθ ∂r � uθ r2 þ ∂2 uθ ∂z2

∂2 uθ ∂r2 þ 1 r ∂ur ∂r � ur r2 þ ∂2 ur ∂z2

η ! ∞ : dU=dη ¼ θ ¼ 0:

Pr 0 0.01 0.1 1 10 100 ∞

0.5970 0.597

0.6694 0.670

N/A 0.4564

254 Finite Element Method - Simulation, Numerical Analysis and Solution Techniques

4. Linear stability of Taylor-Couette flow

cylindrical coordinate system as follows:

∂ur ∂r þ uz ∂ur ∂z � uθ

> ∂uθ ∂r þ uz ∂uθ ∂z þ uruθ r ¼ ν

∂ur ∂t þ ur

> ∂uθ ∂t þ ur

and that for Pr = 0 is 0.6004.

0.7107 0.6922

Nux ð Þ Rax 1 = 4

Nux ð Þ RaxPr 1= 4

Ozoe [4]).

4.1. Governing equations

$$\overline{p}(r,z) = \int \frac{\rho \{\overline{u}\_{\theta}(r)\}^2}{r} \, dr - \rho \mathbf{g}z + p\_0 \tag{48}$$

Here, Ω1 is the angular velocity at the inner cylinder, Ω2 is the angular velocity at the outer cylinder, p is the pressure, r is the density, and g is the acceleration due to gravity. In order to derive disturbance equations for the linear stability, the three components of velocity and pressure are represented as a summation of basic state and infinitesimal disturbance as follows:

$$\begin{aligned} u\_{\theta}(r,z,t) &= \overline{u}\_{\theta}(r) + \overline{v}'(r,z,t), & u\_{r}(r,z,t) &= \overline{u}'(r,z,t), & u\_{z}(r,z,t) &= \overline{w}'(r,z,t), \\ p(r,z,t) &= \overline{p}(r,z) + \overline{p}'(r,z,t) \end{aligned} \tag{49}$$

After neglecting the second-order disturbance, the following linearized equations are obtained:

$$\frac{\partial \mu'}{\partial r} + \frac{\mu'}{r} + \frac{\partial w'}{\partial z} = 0 \tag{50}$$

R ¼ r r2

ReΩ ¼ Ω1r2

, Uθ ¼ uθ

2 ν

the computational results:

Ta ¼ 4Ω1

2 r4 1 ν2

Ω1r2

The boundary conditions are as follows:

, k ¼ r2a, η ¼ r1

(

, U~ ; V~ ; W~ � � ¼ ð Þ u~; ~v; w~

, μ ¼ Ω2 Ω1

r2

1 � μ � � 1 � 4μ � �

found in the recent papers published by the present author [9, 10]

Figure 5. The staggered grids in the radius direction together with the points of each variable definition.

1 � η2 ð Þ2 ¼ 4Re 2

Ω1r2

, P~ ¼ ~p

, S ¼ s ν=r2

R ¼ η : U~ ¼ V~ ¼ W~ ¼ 0 Inner wall ð Þ R ¼ 1 : U~ ¼ V~ ¼ W~ ¼ 0 Outer wall ð Þ:

After Chandrasekar [6], the following two non-dimensional numbers are introduced to verify

Ω

In this section, it is assumed that S = 0. This indicates that the secondary flow caused by the centrifugal instability is stationary and it contains toroidal vortices. To deal with the simultaneous ordinary differential equations for the boundary value problem, a one-dimensional staggered grid system is employed as shown in Figure 5. All the equations are discretized by the fourth order central difference method with a given wavenumber k using the HSMAC method [8] during which ReΩ is obtained by the Newton method. The following equations are used for correction of the pressure and velocity simultaneously. Here, the subscript i indicates grid location, while the superscripts m and n indicate the iteration of the corrections for the convergence of Eq. (56) and the time step, respectively. The more detailed explanation can be

rνΩ1 ,

η4 1 � μ � � 1 � 4μ � �

1 � η2 ð Þ2 , κ ¼ 1 � μ=η2

2 , D ¼ d

Numerical Analysis of the Incompressible Fluid Flow and Heat Transfer

dR , D∗ ¼ d

http://dx.doi.org/10.5772/intechopen.72263

dR þ 1 R (60)

257

(61)

1 � μ (62)

$$\frac{\partial \overline{u'} }{\partial t} = -\frac{1}{\rho} \frac{\partial p'}{\partial r} + \nu \left( \frac{\partial^2 \underline{u'}}{\partial r^2} + \frac{1}{r} \frac{\partial \underline{u'}}{\partial r} - \frac{\underline{u'}}{r^2} + \frac{\partial^2 \underline{u'}}{\partial z^2} \right) + \frac{2 \overline{u}\_\theta \underline{v'}}{r} \tag{51}$$

$$\frac{\partial \overline{v}'}{\partial t} = \nu \left( \frac{\partial^2 \overline{v}'}{\partial r^2} + \frac{1}{r} \frac{\partial \overline{v}'}{\partial r} - \frac{\overline{v}'}{r^2} + \frac{\partial^2 \overline{v}'}{\partial z^2} \right) - \left( \frac{d \overline{u}\_\theta}{dr} + \frac{\overline{u}\_\theta}{r} \right) \underline{u}' \tag{52}$$

$$\frac{\partial w^{\prime}}{\partial t} = -\frac{1}{\rho} \frac{\partial p^{\prime}}{\partial z} + \nu \left( \frac{\partial^2 w^{\prime}}{\partial r^2} + \frac{1}{r} \frac{\partial w^{\prime}}{\partial r} + \frac{\partial^2 w^{\prime}}{\partial z^2} \right) \tag{53}$$

By considering the periodicity of the secondary flow which could be happened, each component of infinitesimal disturbance is assumed to be given in the following form. Here, a is the axial wavenumber (real number) and s is angular frequency (complex number)

$$\frac{\stackrel{\text{u}'}{\tilde{u}(r)}}{\tilde{u}(r)} = \frac{\stackrel{\text{v}'}{\tilde{v}(r)}}{\tilde{v}(r)} = \frac{\stackrel{\text{w}'}{\tilde{w}(r)}}{\tilde{p}(r)} = \frac{\stackrel{\text{p}'}{\tilde{p}(r)}}{\tilde{p}(r)} = \exp\left(iaz + st\right) \tag{54}$$

#### 4.3. Linear stability analysis

The dimensionless simultaneous ordinary equations are summarized as follows:

Basic velocity

$$\overline{\mathcal{U}}\mathcal{U}\_{\theta}(\mathcal{R}) = \frac{\mu - \eta^2}{1 - \eta^2} \mathcal{R} + \frac{\eta^2 \left(1 - \mu\right)}{1 - \eta^2} \frac{1}{\mathcal{R}} \tag{55}$$

Disturbance equations for amplitude functions

$$D\_\*\tilde{\mathcal{U}} + i\mathcal{k}\tilde{\mathcal{W}} = 0\tag{56}$$

$$
\Delta \tilde{s} \tilde{\boldsymbol{L}} = -D \tilde{\boldsymbol{P}} + \left( \boldsymbol{D} \boldsymbol{D}\_{\ast} - \boldsymbol{k}^{2} \right) \tilde{\boldsymbol{\mathcal{U}}} + \mathrm{Re}\_{\Omega} \frac{2 \mathcal{U}\_{\theta}}{R} \tilde{\boldsymbol{V}} \tag{57}
$$

$$
\tilde{SV} = \left(\mathbf{D}\mathbf{D}\_\* - k^2\right)\tilde{V} - \mathrm{Re}\_\Omega \left(\mathbf{D}\_\* \overline{\mathbf{U}}\_\theta\right)\tilde{\mathbf{U}}\tag{58}
$$

$$\mathcal{S}\tilde{\mathcal{W}} = -\mathrm{i}k\tilde{P} + \left(D\_{\ast}D - k^{2}\right)\tilde{\mathcal{W}}\tag{59}$$

Here, the dimensionless variables and non-dimensional numbers are as follows.  Linear stability analysis

Basic velocity

ð Þ r; z; t , urð Þ¼ r; z; t u

∂u 0 ∂r þ u 0 r þ ∂w0

> ∂2 u 0 ∂r2 þ 1 r ∂u 0 ∂r � u 0 r2 þ ∂2 u 0 ∂z2

axial wavenumber (real number) and s is angular frequency (complex number)

The dimensionless simultaneous ordinary equations are summarized as follows:

Uθð Þ¼ R μ � η2

SU~ ¼ �DP~ þ DD∗ � k

SV~ ¼ DD∗ � k

0 ~v rð Þ ¼ w0

After neglecting the second-order disturbance, the following linearized equations are obtained:

∂2 w0 ∂r2 þ

By considering the periodicity of the secondary flow which could be happened, each component of infinitesimal disturbance is assumed to be given in the following form. Here, a is the

0

1 � η2 R þ η2 1 � μ

2 U~ þ ReΩ

2 V~ � ReΩ D∗Uθ

Here, the dimensionless variables and non-dimensional numbers are as follows. The outer

1 � η2

1

D∗U~ þ ikW~ ¼ 0 (56)

2Uθ

SW~ ¼ �ikP~ þ D∗D � k2 W~ (59)

w r ~ ð Þ ¼ p

1 r ∂w0 ∂r þ ∂2 w0 ∂z2

0

ð Þ r; z; t , uzð Þ¼ r; z; t w0

þ 2uθv 0

uθ r

u 0

~p rð Þ ¼ exp ð Þ iaz þ st (54)

R (55)

R V~ (57)

U~ (58)

� duθ dr þ

∂z ¼ 0 (50)

ð Þ r; z; t (49)

ð Þ r; z; t ,

r (51)

(52)

(53)

$$\begin{cases} R = \eta: & \tilde{\mathcal{U}} = \tilde{V} = \tilde{W} = 0 \quad \text{(Inner wall)}\\ R = 1: & \tilde{\mathcal{U}} = \tilde{V} = \tilde{W} = 0 \quad \text{(Outer wall)}. \end{cases} \tag{61}$$

After Chandrasekar [6], the following two non-dimensional numbers are introduced to verify the computational results:

$$Ta = \frac{4\Omega\_1^2 r\_1^4}{\nu^2} \frac{(1-\mu)\left(1-4\mu\right)}{\left(1-\eta^2\right)^2} = 4\text{Re}\_\Omega^2 \frac{\eta^4 \left(1-\mu\right)\left(1-4\mu\right)}{\left(1-\eta^2\right)^2}, \quad \kappa = \frac{1-\mu/\eta^2}{1-\mu} \tag{62}$$

In this section, it is assumed that S = 0. This indicates that the secondary flow caused by the centrifugal instability is stationary and it contains toroidal vortices. To deal with the simultaneous ordinary differential equations for the boundary value problem, a one-dimensional staggered grid system is employed as shown in Figure 5. All the equations are discretized by the fourth order central difference method with a given wavenumber k using the HSMAC method [8] during which ReΩ is obtained by the Newton method. The following equations are used for correction of the pressure and velocity simultaneously. Here, the subscript i indicates grid location, while the superscripts m and n indicate the iteration of the corrections for the convergence of Eq. (56) and the time step, respectively. The more detailed explanation can be found in the recent papers published by the present author [9, 10]

Figure 5. The staggered grids in the radius direction together with the points of each variable definition.

$$\begin{split} \, ^{m+1}\hat{P}\_{i}^{n+1} &= ^{m}\hat{P}\_{i}^{n+1} + ^{m}(\delta\hat{P})\_{i}^{n+1} \\ &= ^{m}\hat{\boldsymbol{U}}\_{i}^{n+1} - \frac{\, ^{m}\hat{\boldsymbol{U}}\_{i+1}^{n+1} + 2\mathcal{T}^{m}\hat{\boldsymbol{U}}\_{i}^{n+1} - 2\mathcal{T}^{m}\hat{\boldsymbol{U}}\_{i-1}^{n+1} + \, ^{m}\hat{\boldsymbol{U}}\_{i-2}^{n+1}}{24(\Delta\boldsymbol{R})} + \frac{\, ^{m}\hat{\boldsymbol{U}}\_{i+1}^{n+1} + \theta^{m}\hat{\boldsymbol{U}}\_{i}^{n+1} + \theta^{m}\hat{\boldsymbol{U}}\_{i-1}^{n+1} - \, ^{m}\hat{\boldsymbol{U}}\_{i-2}^{n+1}}{16R\_{\rm R}} + \, \text{ik} \cdot ^{m}\hat{\boldsymbol{W}}\_{i}^{n+1} \end{split} \tag{63}$$

$$\mathbf{U}^{m+1} \tilde{\mathbf{U}}\_i^{n+1} = \prescript{m}{}{\mathbf{U}}\_i^{n+1} + \frac{\Delta \mathbf{r}}{\Delta \mathbf{R}} \cdot \prescript{m}{}{\left(\delta \tilde{\mathbf{P}}\right)^{n+1}\_i} \quad \overset{m+1}{\ }{\mathbf{U}}\_{i-1}^{n+1} = \prescript{m}{}{\tilde{\mathbf{U}}}\_{i-1}^{n+1} - \frac{\Delta \mathbf{r}}{\Delta \mathbf{R}} \cdot \prescript{m}{}{\left(\delta \tilde{\mathbf{P}}\right)^{n+1}\_i} \tag{64}$$

$$\boldsymbol{^{m+1}i\tilde{W}^{n+1}\_{i}} = {}^{m}\tilde{W}^{n+1}\_{i} - i\mathbf{k}(\Delta\boldsymbol{\pi}) \cdot {}^{m} \left(\delta\tilde{P}\right)^{n+1}\_{i} \tag{65}$$

Table 2 shows the computational results for various rotation speeds at η = 0.5. When μ > 0.25, the basic flow is always stable due to the Rayleigh's criterion. The present results exhibit slightly smaller values of Taylor number than those of Chandrasekar. Figures 6 and 7 show the amplitude functions and Eigen functions, respectively, for the case of μ = 0 (the outer cylinder is stationary), and Figures 8 and 9 show the case of μ = �0.5 (the outer cylinder rotates with half angular velocity in opposite direction to the inside rotation).

The simultaneous ordinary equations from (56) to (59) were divided into the real and imaginary parts. However, only four equations among the eight equations are necessary to solve in this problem because of the symmetricity and anti-symmetricity of the complex variables. In Figures 6 and 8, the real part of U, ~ V , ~ P~ and the imaginary part of W~ are shown. For the visualization shown in Figures 7 and 9, the Stokes stream function Ψ is defined as follows:

U~ ℜ � cos ð Þ¼ kZ 1

the similar manner using the trigonometric functions.

function, azimuthal velocity, and pressure.

Figure 6. Amplitude functions (η = 0.5, μ = 0, k = 6.325).

R ∂Ψ

Here, the subscripts ℜ and ℑ represent the real part and the imaginary part, respectively. The visualization of other variables, such as the azimuthal velocity and the pressure, are treated in

Figure 7. Visualization of Eigen functions for two wavelengths (η = 0.5, μ = 0, k = 6.325). From left to right, Stokes stream

∂Z , W~ ℑ � sin ð Þ¼ kZ 1

R ∂Ψ

Numerical Analysis of the Incompressible Fluid Flow and Heat Transfer

http://dx.doi.org/10.5772/intechopen.72263

259

∂R (66)

Table 2. Computational results and comparison with Chandrasekar (η = 0.5).

Figure 6. Amplitude functions (η = 0.5, μ = 0, k = 6.325).

mþ1 P~ nþ1 i ¼ mP~

> ¼ mP~ nþ1 i �

nþ1 i þ mðδP~Þ

mþ1 U~ nþ1

defined as follows:

�mU~ nþ1

nþ1 i

258 Finite Element Method - Simulation, Numerical Analysis and Solution Techniques

iþ1 þ 27mU~ nþ1

i ¼ mU~ nþ1

i þ

mþ1 W~ nþ1

variables. In Figures 6 and 8, the real part of U,

Table 2. Computational results and comparison with Chandrasekar (η = 0.5).

i � 27mU~ nþ1

Δτ

with half angular velocity in opposite direction to the inside rotation).

ΔR � m δP~ nþ1

i ¼ mW~ nþ1

i�1 þ mU~ nþ1 i�2 24ðΔRÞ þ �mU~ nþ1

Δτf2=ðΔRÞ

i , mþ1

Table 2 shows the computational results for various rotation speeds at η = 0.5. When μ > 0.25, the basic flow is always stable due to the Rayleigh's criterion. The present results exhibit slightly smaller values of Taylor number than those of Chandrasekar. Figures 6 and 7 show the amplitude functions and Eigen functions, respectively, for the case of μ = 0 (the outer cylinder is stationary), and Figures 8 and 9 show the case of μ = �0.5 (the outer cylinder rotates

The simultaneous ordinary equations from (56) to (59) were divided into the real and imaginary parts. However, only four equations among the eight equations are necessary to solve in this problem because of the symmetricity and anti-symmetricity of the complex

shown. For the visualization shown in Figures 7 and 9, the Stokes stream function Ψ is

κ μ Critical wave number Critical Ta number Wavenumber Ta number 0 1/4 6.286 15316 6.4 15332 0.4 1/6 6.293 19518 6.4 19542 0.6 2/17 6.299 22617 6.4 22644 1.0 0 6.325 33062 6.4 33100 4/3 �1/8 6.403 53210 6.4 53280 1.6 �1/4 6.715 98520 6.4 99072 1.8 �4/11 7.819 197715 7.8 199540 1.9 �9/21 8.733 288761 8.6 293630 2.0 �1/2 9.602 417734 9.6 428650

~ V ,

Present (201 grids) Chandrasekar [6]

iþ1 þ 9mU~ nþ1

i�1 ¼ mU~ nþ1

2 þ k 2 g

U~ nþ1

i � ikð Þ� Δτ m δP~ nþ1

i þ 9mU~ nþ1

16RPi

i�1 � Δτ

i�1 � mU~ nþ1 i�2

ΔR � m δP~ nþ1

~ P~ and the imaginary part of W~ are

i (65)

þ ik � mW~ nþ1 i

i (64)

(63)

Figure 7. Visualization of Eigen functions for two wavelengths (η = 0.5, μ = 0, k = 6.325). From left to right, Stokes stream function, azimuthal velocity, and pressure.

$$
\tilde{M}\_{\mathfrak{R}} \cdot \cos(kZ) = \frac{1}{R} \frac{\partial \Psi}{\partial Z}, \quad \tilde{W}\_{\mathfrak{F}} \cdot \sin(kZ) = \frac{1}{R} \frac{\partial \Psi}{\partial R} \tag{66}
$$

Here, the subscripts ℜ and ℑ represent the real part and the imaginary part, respectively. The visualization of other variables, such as the azimuthal velocity and the pressure, are treated in the similar manner using the trigonometric functions.

References

Springer; 1999

zine. 1911;21:697-711

78(794):1680-1695

Japan; August 2014

Journal of Heat Transfer. 1973;95(4):540-541

Laboratory of University of California. 1975. p. LA-5852

Papers III

1961

[1] Schlichting H, Gersten K. Boundary Layer Theory. 8th Revised and Enlarged ed. Berlin:

Numerical Analysis of the Incompressible Fluid Flow and Heat Transfer

http://dx.doi.org/10.5772/intechopen.72263

261

[2] Stokes GG. On the effect of the internal friction of fluids on the motion of pendulums. Transactions of the Cambridge Philosophical Society. 1856;9(Part II):8-106 or Collected

[3] Rayleigh L. On the motion of solid bodies through viscous liquids. Philosophical Maga-

[4] Churchill SW, Ozoe H. A correlation for laminar free convection from a vertical plate.

[5] Le Fevre EJ. Laminar free convection from a vertical plane surface. In: Proceedings of the 9th International Congress of Applied Mechanics; Brussels; Vol. 4; 1956. pp. 168-174 [6] Chandrasekar S. Hydrodynamic and Hydromagnetic Stability. Oxford: Clarendon Press;

[7] Koschmieder EL. Bénard cells and Taylor vortices. In: Cambridge Monographs on Mechanics and Applied Mathematics. United State of America: Cambridge University Press; 1993

[8] Hirt CW, Nichols BD, Romero NC. SOLA: A Numerical Solution Algorithm for Transient Fluid Flows. Los Alamos, New Mexico, United State of America: Los Alamos Scientific

[9] Tagawa T, Egashira R. Fluid flow of a liquid metal in a cylinder driven by a rotating magnetic field. Transactions of the Japan Society of Mechanical Engineers, Part B. 2012;

[10] Tagawa T. Numerical investigation of Bénard-Marangoni convection of paramagnetic liquid in annular layers. In: The 15th International Heat Transfer Conference; Kyoto,

Figure 8. Amplitude functions (η = 0.5, μ = 0.5, k = 9.602).

Figure 9. Visualization of Eigen functions for two wavelengths (η = 0.5, μ = 0.5, k = 9.602). From left to right, Stokes stream function, azimuthal velocity, and pressure.

### Author details

Toshio Tagawa

Address all correspondence to<EMAIL_ADDRESS>
Department of Aerospace Engineering, Tokyo Metropolitan University, Japan


*Article* **Study on the Design Stage from a Dimensional and Energetic Point of View for a Marine Technical Water Generator Suitable for a Medium Size Container Ship**

**Catalin Faitar 1,2,\* and Eugen Rusu 1**

**Abstract:** The purpose of this study is to provide an overview of the development of the modern low-speed marine two-stroke diesel engine from the point of view of the technical water cooling plant, taking into account and starting from the market requirements for power and speed, with information and design options relevant to the entire shipping industry. Thus, through the ideas of this project, we analyze notions and relevant aspects of systems related to marine slow turning engines, including the basic thermodynamic structure of the technical water system designed for a marine engine. This study also presents the design criteria that define the size and design concept of the engine structure components, with a focus on the technical water cooling installation. The concepts for the main engine hot parts served by the technical water cooling installation play a vital role in the marine technical water generator. The choices the engine designer must make regarding basic auxiliary systems, such as fuel injection and exhaust valve actuation, are important factors to keep in mind when installing a technical water generator onboard a ship. The automation and control systems that govern the modern electronic engine, which drive the supply pump of the technical water cooling system can provide a simplified view of the engine development process. In order to point out the contributions of this study, it is important to focus on the calculations used to determine the main parameters of a technical water generator especially designed for a midsized container ship.

**Keywords:** technical; desalination; marine; seawater; operation; generator


**Author details**

Retno Indrati

Department of Food and Agricultural Product Technology, Universitas Gadjah Mada, Yogyakarta, Indonesia

\*Address all correspondence to<EMAIL_ADDRESS>
© 2021 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/ by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.


**Cavitation-Enhanced Drug Delivery and Immunotherapy**

Editors

**Brandon Helfield Shashank Sirsi James Kwan Michael Gray**

Basel • Beijing • Wuhan • Barcelona • Belgrade • Novi Sad • Cluj • Manchester

*Editors* Brandon Helfield Department of Physics, Department of Biology Concordia University Montreal, QC, Canada

Shashank Sirsi Erik Jonsson School of Engineering and Computer Science The University of Texas at Dallas Richardson, TX, USA

James Kwan Department of Engineering Science University of Oxford Oxford, UK

Michael Gray Institute of Biomedical Engineering University of Oxford Oxford, UK

*Editorial Office* MDPI St. Alban-Anlage 66 4052 Basel, Switzerland

This is a reprint of articles from the Special Issue published online in the open access journal *Pharmaceutics* (ISSN 1999-4923) (available at: https://www.mdpi.com/journal/pharmaceutics/ special issues/Microbubbles Delivery).

For citation purposes, cite each article independently as indicated on the article page online and as indicated below:

Lastname, A.A.; Lastname, B.B. Article Title. *Journal Name* **Year**, *Volume Number*, Page Range.

**ISBN 978-3-0365-9434-7 (Hbk) ISBN 978-3-0365-9435-4 (PDF) doi.org/10.3390/books978-3-0365-9435-4**

© 2023 by the authors. Articles in this book are Open Access and distributed under the Creative Commons Attribution (CC BY) license. The book as a whole is distributed by MDPI under the terms and conditions of the Creative Commons Attribution-NonCommercial-NoDerivs (CC BY-NC-ND) license.