# On the Onset of Hydrodynamic Instability with Convective Heat Transfer Through a Rotating Curved Rectangular Duct

On the Onset of Hydrodynamic Instability with Convective Heat Transfer Through a Rotating Curved Rectangular Duct

Samir Chandra Ray Mohammad Sanjeed Hasan Rabindra Nath Mondal*

Department of Mathematics, Bangabandhu Sheikh Mujibur Rahman Science and Technology University, Gopalganj 8100, Bangladesh

Department of Mathematics, Jagannath University, Dhaka 1100, Bangladesh

Corresponding Author Email:
rnmondal71@yahoo.com
Page:
31-44
|
DOI:
https://doi.org/10.18280/mmep.070105
5 November 2019
|
Revised:
7 January 2020
|
Accepted:
12 January 2020
|
Available online:
31 March 2020
| Citation

OPEN ACCESS

Abstract:

The present paper addresses numerical prediction of hydrodynamic instability with convective heat transfer through a rotating curved rectangular duct of curvature $0.1$. The bottom wall of the duct is heated while cooling from the ceiling. Numerical calculations are carried out by using a spectral method and covering a wide range of the Taylor number $0\le Tr\le 2000$ for the constant pressure gradient force, the Dean number, $Dn=1000$. First, solution structure of the steady solutions is investigated. As a result, three branches of asymmetric steady solutions with two- to ten-vortex solutions are obtained by using Newton-Raphson iteration method. Then unsteady solutions are obtained by time evolution calculations and flow transitions are well justified by obtaining the phase space and power spectrum of the solutions. It is found that chaotic flow turns into steady-state flow through periodic oscillating flow, if Tr is increased. Streamlines and isotherms are also obtained at several values of Tr, and it is found that the unsteady flow consists of two- to ten-vortex solutions. The present study shows that combined action of the centrifugal-Coriolis-buoyancy forces contribute to generate the vorticity. The present study exposes the role of secondary vortices on convective heat transfer, which shows that convective heat transfer is significantly enhanced by the secondary flow; and the chaotic flow, which occurs at small Tr but at large Dn, enhances heat transfer more effectively than the steady-state or periodic solutions.

Keywords:

rotating curved duct, steady solutions, time evolution calculation, Taylor number, secondary flow

1. Introduction

Fluid flow through curved ducts and channels has been extensively studied over a wide range of applications because of their enormous applications in fluids engineering such as in turbo-machinery, refrigeration, air conditioning systems, heat exchangers, gas turbines, rocket engine, internal combustion engines and also in aviation engineering. In the past decades, there have been a lot of works on fluid flow and heat transfer of rotating arrangement. The scholars are mentioned to Dierich and Nikrityuk [1], Chatterjee et al. [2], Khaled [3] and Wei-Cheng et al. [4] for some distinguished reviews on heat transfer and fluid flows in different sectors of engineering. In a curved duct, it is strongly anticipated that centrifugal forces are originated in the flow on account of curvature provoking an opposite directional revolving vortex rate acted on the axial direction through the duct that generates the properties of spiraling motion in the bending pathway acquainted as secondary flow. In a convinced precise flow agreement and beyond, an extra couple of counter-rotating vortices develop at the exterior wall of the duct which is widely entitled as Dean Vortices [5]. Many theoretical and experimental studies have been conducted by considering this flow; for example, the articles by Berger et al. [6], Mondal et al. [7, 8] and Chen et al. [9] may be referenced.

Hydrodynamic instability in the curved duct is a prevalent phenomenon which is resulted from imbalance amongst the radial pressure gradient and centrifugal forces. This instability is subjected to Dean number (Dn) [10] and physical parameters of the curved duct’s aspect ratio and curvature ratio. Several studies are conducted to investigate hydrodynamic stability of fluid flow in the curved ducts. Sugiyama et al. [11] investigated laminar flow in the curved ducts under different aspect ratios and curvature ratios, experimentally. Winters [12] studied the symmetric and asymmetric steady solutions in the curved duct via finite element method. Different branches of the steady solutions on linear hydrodynamic stability of laminar flow in curved duct are numerically studied by Yanase et al. [13]. Hasan et al. [14] have interpreted the effect of steady solution structures as well as the linear stability for rotating curved duct. Chandratilleke and Nursubyakto [15] indicated the remarkable effect of the aspect ratio on Dean hydrodynamic instability by numerical simulation. In a numerical study, Fellouah et al. [16] investigated the flow behavior in $180^{\circ}$  curved duct at various curvature ratios and aspect ratios.

The fluid flow in a rotating curved duct generates centrifugal and Coriolis force. Such rotating passages are used in many engineering applications e.g. in cooling system for conductors of electrical generators. For isothermal flows of a constant property fluid, the Coriolis force tends to generate vortices while centrifugal force is purely hydrostatic by Zhang et al. [17]. Mondal et al. [18] performed comprehensive numerical study on fully developed bifurcation structure and stability of two-dimensional (2D) flow through a curved duct with square cross section and found a close relationship between the unsteady solutions and the bifurcation diagram of steady solutions. When a temperature induced variation of fluid density occurs for non-isothermal flows, centrifugal type buoyancy forces can contribute to the generation of vorticity [19]. These effects counteract each other in a non-linear manner depending on the direction of wall, resistance coefficient and the flow domain. Therefore, the effect of counteract of the system is more subtle and complicated and yields new; richer features of flow and heat transfer in general, bifurcation and stability in particular, for non-isothermal flows.

Unsteady solutions of fully developed curved duct flows were initiated by Yanase and Nishiyama [20] for a rectangular cross section. In that study, they investigated unsteady solutions for the case where dual solutions exist. The time-dependent behavior of the flow in a curved rectangular duct of large aspect ratio was investigated, in detail, by Yanase et al. [13] numerically. They performed time-evolution calculations of the unsteady solutions with and without symmetry condition and showed that periodic oscillations appear with symmetry condition while aperiodic time variation without symmetry condition. Hasan et al. [21] have analyzed the flow velocity and temperature profiles for periodic, multi-periodic and chaotic flow in a curved square duct. Wang and Liu [22] performed numerical as well as experimental investigations of periodic oscillations for the fully developed flow in a curved square duct. Flow visualization in the range of Dean numbers from 50 to 500 was conducted in their experiment. Mondal et al. [7] applied spectral method to study non-isothermal flow through a stationary curved rectangular duct of aspect ratios 1 to 3, and showed that the steady-state flow turns into chaotic flow through various flow instabilities if the aspect ratio is increased. Very recently, Islam et al. [23] conducted spectral numerical study of Dean-Taylor flow through a rotating coiled rectangular duct, where they performed unsteady solutions of the flow for both rotating and co-rotating cases and discussed the role of secondary vortices on convective heat transfer. However, solution structure as well as transient behavior of the unsteady solution is not yet resolved for the flow through a rotating curved rectangular duct with bottom wall heating and cooling from the ceiling, which motivated the present study to fill up this gap.

The significant inflictions of flow are to increase the thermal passage between two walls because it is probable that the Dean flow may contribute to transport energy and then soar up heat flow between these two walls. Norouzi et al. [24] investigated inertial and creeping flow of a second-order fluid with convective heat transfer in a curved square duct by using finite difference method. Chandratilleke et al. [25] presented a numerical investigation to examine the secondary vortex motion and heat transfer process in fluid flow through curved rectangular ducts of aspect ratios 1 to 6. The study formulated an improved simulation model based on 3-dimensional vortex structures for describing secondary flow and its thermal characteristics. Zhang et al. [26] adopted finite volume method and SIMPLE algorithm to perform laminar and turbulent flow through a curved square duct at low Reynolds number. There was a good agreement between their numerical result and the experimental data. Heat transfer between the cooled and heated side walls as well as temperature gradients for a wide range of Taylor number have been enumerated by Hasan et al. [27]. Wu et al. [28] performed numerical study of the secondary flow characteristics in a curved square duct by using spectral method, where the walls of the duct except the outer wall rotate around the centre of curvature and an azimuthal pressure gradient was imposed. Recently, Mondal et al. [29] investigated combined effects of centrifugal and Coriolis instability of the isothermal/non-isothermal flows through a rotating curved rectangular duct numerically. The secondary flow characteristics in a curved square duct were investigated experimentally by using visualization method by Yamamoto et al. [30]. Recently, Razavi et al. [31] investigated flow characteristics, heat transfer and entropy generation in a rotating curved duct by using control volume method. The effects of Dean number, wall heat flux and force ratio on the entropy generation were presented in that paper. Very recently, Li et al. [32] conducted a combined experimental and numerical study on 3D flow development in a curved rectangular duct with varying curvature. Effects of curvature, Reynolds number and aspect ratio on hydrodynamic instability were discussed in that paper to accurately predict the core of secondary base vortices. To the best of the authors' knowledge, however, there has not yet been done any work studying the solution structure as well as effects of secondary vortices on convective heat transfer for the non-isothermal flow through a rotating curved rectangular duct whose bottom wall is heated and cooling from the ceiling. But from the scientific as well as engineering point of view it is quite interesting because this type of flow is often encountered in engineering applications such as in gas turbines, metallic industry etc.

Examining the unique features of secondary flow and heat transfer, the main objective of the present study is to discuss the solution structure of the steady solutions and to investigate time-dependent behavior of the unsteady solutions through a rotating curved rectangular duct whose bottom wall is heated and cooling from the ceiling.

2. Mathematical Formulations

Consider a hydro-dynamically and thermally fully developed two-dimensional (2D) flow of viscous incompressible fluid through a rotating coiled rectangular duct, whose height and width are $2h$ and $2d$, respectively. The coordinate system with the relevant notation is shown in Figure 1. The system rotates at a constant angular velocity ${{\Omega }_{T}}$ around the $y'$ axis. It is assumed that the bottom wall of the duct is heated while cooling from the ceiling, the inner and outer walls being thermally insulated. It is also assumed that the flow is uniform in the axial direction, which is driven by a constant pressure gradient along the center-line of the duct as shown in Figure 1. The dimensional variables are made non-dimensional by using the representative length $d$, the representative velocity ${{U}_{0}}=\frac{\upsilon }{d}$, where $\upsilon$ is the kinematic viscosity of the fluid. $\delta$ is the curvature of the duct defined as $\delta =\frac{d}{L}$.

## 1.png

Figure 1. Coordinate system of the curved rectangular duct

By assumption of rotationally symmetric flow, stream functions for cross-sectional velocities have the following form

\left. \begin{align} & u=\frac{1}{r}\frac{\partial \psi }{\partial y}=\frac{1}{1+\delta x}\frac{\partial \psi }{\partial y} \\ & v=\frac{1}{r}\frac{\partial \psi }{\partial x}=-\frac{1}{1+\delta x}\frac{\partial \psi }{\partial x} \\\end{align} \right\}     (1)

Eq. (1) satisfies the continuity equation. Now, stream-wise velocity $\left( w \right)$, cross-sectional stream-function $\left( \psi \right)$ and $T$ are defined based on the Navier-Stokes equation as follows:

\begin{align} & (1+\delta x)\frac{\partial w}{\partial t}=Dn-\frac{1}{2}\frac{\partial (w,\psi )}{\partial (x,y)}-\frac{{{\delta }^{2}}w}{1+\delta x}+(1+\delta x){{\Delta }_{2}}w \\ & -\frac{1}{2}\frac{\delta }{(1+\delta x)}\frac{\partial \psi }{\partial y}w+\delta \frac{\partial w}{\partial x}-\delta Tr\frac{\partial \psi }{\partial y} \\\end{align}    (2)

\begin{align} & \left( {{\Delta }_{2}}-\frac{\delta }{1+\delta x}\frac{\partial }{\partial x} \right)\frac{\partial \psi }{\partial t}=-\frac{1}{2}\frac{1}{(1+\delta x)}\frac{\partial ({{\Delta }_{2}}\psi ,\psi )}{\partial (x,y)}+\Delta _{2}^{2}\psi \\ & +\frac{1}{2}\frac{\delta }{{{\left( 1+\delta x \right)}^{2}}}\times \left[ \frac{\partial \psi }{\partial y}\left( 2{{\Delta }_{2}}\psi -\frac{3\delta }{1+\delta x}\frac{\partial \psi }{\partial x}+\frac{{{\partial }^{2}}\psi }{\partial {{x}^{2}}} \right) \right. \\ & \left. -\frac{\partial \psi }{\partial x}\frac{{{\partial }^{2}}\psi }{\partial x\partial y} \right]+\frac{\delta }{{{\left( 1+\delta x \right)}^{2}}}\times \left[ 3\delta \frac{{{\delta }^{2}}\psi }{\partial {{x}^{2}}}-\frac{3{{\delta }^{2}}}{1+\delta x}\frac{\partial \psi }{\partial x} \right] \\ & -\frac{2\delta }{1+\delta x}\frac{\partial }{\partial x}{{\Delta }_{2}}\psi +\frac{1}{2}w\frac{\partial w}{\partial y}-Gr(1+\delta x)\frac{\partial T}{\partial x}-\frac{1}{2}Tr\frac{\partial \psi }{\partial y}, \\\end{align}    (3)

$\frac{\partial T}{\partial t}=\frac{1}{\Pr }\left( {{\Delta }_{2}}T+\frac{\delta }{1+\delta x}\frac{\partial T}{\partial x} \right)-\frac{1}{(1+\delta x)}\frac{\partial (T,\psi )}{\partial (x,y)}$     (4)

where, ${{\Delta }_{2}}\equiv \frac{{{\partial }^{2}}}{\partial {{x}^{2}}}+4\frac{{{\partial }^{2}}}{\partial {{y}^{2}}},\,\,\,\,\,\,\,\frac{\partial \left( f,g \right)}{\partial \left( x,y \right)}\equiv \frac{\partial f}{\partial x}\frac{\partial g}{\partial y}-\frac{\partial f}{\partial y}\frac{\partial g}{\partial x}.$

The equations for $w$, $\psi$ and $T$ are actually benefited for numerical computation. The non-dimensional parameters $Dn$, the Dean number; $Gr$, the Grashof number; $Tr$, the Taylor number, and $Pr$, the Prandtl number, which appear in Eqns. (2) - (4) are defined as:

$D n=\frac{G d^{3}}{\mu \nu} \sqrt{\frac{2 d}{L}}$, $G r=\frac{\beta g \Delta T d^{3}}{v^{2}}$, $T r=\frac{2 \sqrt{2 \delta} \Omega_{T} d^{3}}{v \delta}$, $\mathrm{Pr}=\frac{v}{\kappa}$

The rigid boundary conditions for $w$ and $\psi$ are used as

\begin{align}& w(\pm 1,y)=w(x,\pm 1)=\psi (\pm 1,y)=\psi (x,\pm 1)=\frac{\partial \psi }{\partial x}(\pm 1,y) \\ & =\frac{\partial \psi }{\partial y}(x,\pm 1)=0 \\\end{align}     (5)

and the temperature $T$ is assumed to be constant on the walls as

$T\left( x,1 \right)=1,\text{ }T\left( x,-1 \right)=-1,\text{ }T\left( 1,\pm y \right)=y$     (6)

There is a class of solutions which satisfy the following symmetry condition with respect to the horizontal plane $y=0$.

\left. \begin{align}& w(x,y,t)\Rightarrow w(-x,y,t), \\& \psi (x,y,t)\Rightarrow -\psi (-x,y,t), \\& T(x,y,t)=-T(-x,y,t) \\\end{align} \right\}     (7)

The solution which satisfies the condition (6) is called a symmetric solution, and that which does not an asymmetric solution. In the present study, only $Tr$ vary $\left( 0\le Tr\le 2000 \right)$ while $Dn$, $Gr$, $Pr$ and $\delta$ are fixed as $Dn=1000$, $Gr=100$, $Pr=7.0$ (water) and $\delta =0.1$.

3. Numerical Calculations

3.1 Method of numerical calculation

In order to solve the Eqns. (2)-(4) numerically, the spectral method is used. This is the method which is thought to be the best numerical method for solving the Navier-Stokes as well as energy equations [33]. By this method the variables are expanded in a series of functions consisting of the Chebyshev polynomials. That is, the expansion functions ${{\phi }_{n}}\left( x \right)$ and ${{\psi }_{n}}\left( x \right)$ are expressed as

\left. \begin{align}(8) & {{\varphi }_{n}}(x)=(1-{{x}^{2}}){{C}_{n}}(x), \\ & {{\psi }_{n}}(x)={{(1-{{x}^{2}})}^{2}}{{C}_{n}}(x) \\\end{align} \right\}     (8)

where, ${{C}_{n}}(x)=\cos \left( n{{\cos }^{-1}}(x) \right)$ is the ${{n}^{th}}$ order Chebyshev polynomial. $w\left( x,y,t \right)$, $\psi \left( x,y,t \right)$ and $T\left( x,y,t \right)$ are expanded in terms of the expansion functions ${{\phi }_{n}}\left( x \right)$ and ${{\psi }_{n}}\left( x \right)$ as:

$\left. \begin{matrix} w(x,y,t)=\sum\limits_{m=0}^{M}{\sum\limits_{n=0}^{N}{{{w}_{m\,n}}\,(t)\,{{\varphi }_{m}}\,(x){{\varphi }_{n}}\,(y)\,}} \\ \psi (x,y,t)=\sum\limits_{m=0}^{M}{\sum\limits_{n=0}^{N}{{{\psi }_{m\,n}}\,(t)\,{{\psi }_{m}}\,(x)\,{{\psi }_{n}}\,(y)}} \\ T(x,y,t)=\sum\limits_{m=0}^{M}{\sum\limits_{n=0}^{N}{{{T}_{m\,n}}\,(t)\,{{\varphi }_{m}}\,(x)\,{{\varphi }_{n}}(y)}-y} \\\end{matrix} \right\}$     (9)

where, $M$ and $N$ are the truncation numbers in the $x$- and $y$-directions respectively, and ${{w}_{mn}}$, ${{\psi }_{mn}}$ and ${{T}_{mn}}$ are the coefficients of expansion. To obtain the steady solution, the time derivative terms $\frac{\partial w}{\partial t},$ $\frac{\partial \psi }{\partial t},$ and $\frac{\partial T}{\partial t}$ are taken to be zero and the expansion series (9) with coefficients ${{w}_{mn}}$, ${{\psi }_{mn}}$ and ${{T}_{mn}}$ being time independent, are substituted into the basic Eqns. (2), (3) and (4), and the collocation method [33] is applied. As a result, a set of nonlinear algebraic equations for ${{w}_{mn}}$, ${{\psi }_{mn}}$ and ${{T}_{mn}}$ are obtained. The collocation points $({{x}_{i}},{{y}_{j}})$ are taken to be

\left. \begin{align} & {{x}_{i}}=\cos \left[ \pi \left( 1-\frac{i}{M+2} \right) \right],i=1,\ldots ,M+1 \\ & {{y}_{j}}=\cos \left[ \pi \left( 1-\frac{j}{N+2} \right) \right],j=1,\ldots ,N+1 \\\end{align} \right\}     (10)

Steady solutions are obtained by the Newton-Rapshon iteration method assuming that all the coefficients to be time independent. To avoid difficulty near the point of inflection for the steady solutions, we use the arc-length method. In the arc-length method, the arc-lengths play a central role in the formulation. The arc-length equation is

$\sum\limits_{m=0}^{M}{\sum\limits_{n=0}^{N}{\left\{ {{\left( \frac{d{{w}_{mn}}}{ds} \right)}^{2}}+{{\left( \frac{d{{\psi }_{mn}}}{ds} \right)}^{2}}+{{\left( \frac{d{{T}_{mn}}}{ds} \right)}^{2}} \right\}}}=1$     (11)

The convergence is assured by taking sufficiently small ${{\varepsilon }_{p}}$$\left( {{\varepsilon }_{p}}<{{10}^{-10}} \right)$ defined as

{{\varepsilon }_{p}}=\sum\limits_{m=0}^{M}{\sum\limits_{n=0}^{N}{\left[ \begin{align} & {{(w_{mn}^{\left( p+1 \right)}-w_{mn}^{p})}^{2}}+{{(\psi _{mn}^{\left( p+1 \right)}-\psi _{mn}^{p})}^{2}} \\ & +{{(T_{mn}^{\left( p+1 \right)}-T_{mn}^{p})}^{2}} \\\end{align} \right]}}     (12)

Finally, in order to calculate the unsteady solutions, the Crank-Nicolson and Adams-Bashforth methods together with the function expansion (9) and the collocation methods are applied to Eqns. (2)-(4).

3.2 Resistance coefficient

The resistant coefficient $\lambda$ is used as the representative quantity of the flow state. It is also called the hydraulic resistance coefficient, and is generally used in fluids engineering, defined as

$\frac{{{P}_{1}}^{*}-{{P}_{2}}^{*}}{\Delta {{z}^{*}}}=\frac{\lambda }{{{d}_{h}}^{*}}\frac{1}{2}\rho {{\langle {{w}^{*}}\rangle }^{2}}$     (13)

where, quantities with an asterisk denote dimensional ones, $\left\langle {} \right\rangle$ stands for the mean over the cross section of the duct and ${{d}_{h}}^{*}$ is the hydraulic diameter. The mean axial velocity $\langle {{w}^{*}}\rangle$ is calculated by

$\langle {{w}^{*}}\rangle =\frac{\nu }{4\sqrt{2\delta }d}\int_{-1}^{1}{dx}\int_{-1}^{1}{w(x,y,t)}dy$     (14)

Since $\frac{{{P}_{1}}-{{P}_{2}}}{\Delta z*}=G$, $\lambda$ is related to the mean non-dimensional axial velocity $\langle w\rangle$ as

$\lambda =\frac{16\sqrt{2\delta }Dn}{3{{\left\langle w \right\rangle }^{2}}}$     (15)

where, $\langle w\rangle =\sqrt{2\delta }d{\langle {{w}^{*}}\rangle }/{\nu }\;$.

4. Results and Discussions

In this study, we first investigate solution structure of the steady solutions by Newton-Raphson iteration method for the curvature $\delta =0.1$ and discuss pattern variation of secondary flows on various branches of steady solutions. After an extensive survey, three branches of steady solutions are obtained over the Taylor number $0\le Tr\le 2000$. A bifurcation diagram of steady solution is shown in Figure 2 for $Dn=1000$ and $Gr=100$ using $\lambda$, the representative quantity of the flow state. The three steady solution branches are named the first steady solution branch (Branch 1, black solid line), the second steady solution branch (Branch 2, purple solid line) and the third steady solution branch (Branch 3, blue solid line) respectively. The solution branches are obtained by the path continuation technique with various initial guesses and are distinguished by the nature and number of secondary vortices appearing in the cross section of the duct. It is observed that there is no bifurcating relationship among the three branches of steady solutions. In the following, the three branches of steady solutions, obtained for $\delta =0.1$, are discussed in brief.

## 2.png

Figure 2. Solution structure of steady solutions for $Dn=1000$, $Gr=100$ and $0\le Tr\le 2000$

## 3.png

Figure 3. (a) First steady solution branch for $Dn=1000$ and $Gr=100$(the route of the branch: $a\to b\to c$). (b) Enlargement of Figure 3 (a)

## 4.png

Figure 4. Streamlines of secondary flow (top) and axial flow (middle) and isotherm (bottom) on the first steady solution branch for $Dn=1000$, $Gr=100$ at various values of $Tr$

The first steady solution branch for $\delta =0.1$ is plotted exclusively in Figure 3 designated by black solid line for $249\le Tr\le 2000$. Figure 3 (a) shows that the branch starts from point ‘a’ $\left( Tr=2000 \right)$ and extends to the direction of decreasing $Tr$ up to $Tr=249$ where it experiences a smooth turning and then goes to the direction of increasing $Tr$ up to $Tr=2000$. An enlargement of Figure 3 (a) is shown in Figure 3 (b), where we see that the branch has two sub-branches which closely overlaps each other after turning at point b. Streamlines of secondary flow (top) and axial flow distribution (middle) and isotherms (bottom) on the first steady solution branch are shown in Figure 4 at several values of $Tr$. As seen in Figure 4, the branch consists of asymmetric two-vortex solutions.

## 5.png

Figure 5. (a) Second steady solution branch for $Dn=1000$, $Gr=100$and $0\le Tr\le 850$ (the route of the branch: $\left( a\to b\to c\to d\to e\to f\to g\to h\to i \right)$. (b) Enlargement of Figure 5(a), (c) Enlargement of Figure 5(a) around $Tr=200$

The second steady solution branch for $\delta =0.1$ is exclusively shown in Figure 5(a) by depicting purple solid line for $0\le Tr\le 850$. As seen in Figure 5(a), the branch starts at Tr =0 and extends to the direction of increasing Tr up to Tr = 850 experiencing many turnings on its way. To have a clear view about the turnings of the branch, enlargements of the second branch are shown in Figure 5(b) and in Figure 5(c). Typical contours of secondary flow patterns (top), axial flow distribution (middle) and temperature profiles (bottom) on the second steady solution branch are shown in Figure 6 at several values of Tr, where we see that the branch consists of two- to eight-vortex solutions. Temperature distribution shows that the stream lines are consistent with the secondary vortices and axial flow distribution and heat is transferred from the bottom heated wall to the fluid. Here the contours of $w$, $\psi$ and $T$ are drawn with the increment $\Delta w=4.0$, $\Delta \psi =0.8$, and $\Delta T=0.3$. The same increment of $w$, $\psi$, $T$ are used for all figures in this study, if not specified. The right hand side of each duct box for $w$ and $\psi$ are in the outside direction of the curvature. The solid lines $\left( \psi \ge 0,\text{ T}\ge \text{0} \right)$ show that the secondary flow is in the clockwise direction while the dotted ones $\left( \psi <0 \right)$ in the counter clockwise direction.

## 6.png

Figure 6. Streamlines of secondary flow (top) and axial flow (middle) and isotherms (bottom) on the second steady solution branch for $Dn=1000$, $Gr=100$ at various values of $Tr$

## 7.png

Figure 7. (a) Third steady solution branch for $Dn=1000$, $Gr=100$ and $0\le Tr\le 1200$ (the route of the branch: $a\to b\to c\to d\to e\to f\to g\to h\to i\to j\to k$). (b) Enlargement of Figure 7(a) at large $Tr$, (c) Enlargement of Figure 7(a) at small $Tr$

## 8.png

Figure 8. Streamlines of secondary flow (top) and axial flow (middle) and isotherms (bottom) on the third steady solution branch for $Dn=1000$, $Gr=100$ at various values of $Tr$

The third steady solution branch for $\delta =0.1$ is plotted exclusively in Figure 7(a) for $0\le Tr\le 1200$ designated by blue solid line. As seen in Figure 7(a), the branch starts from $Tr=0$ and extends to the direction of increasing $Tr$ up to $Tr=1200$ having many turnings on its way. The branch is entangled and experiences many smooth turnings and finally goes to the direction of decreasing $Tr$ up to $Tr=0$. Figures 7(b) and 7(c) show enlargements of Figure 7(a). Streamlines of secondary flow and axial flow and isotherms of temperature profiles on the first steady solution branch are shown in Figure 8 at several $Tr$. As seen in Figure 8, the branch is comprised with asymmetric two- and three-vortex solutions. Temperature distribution shows that the stream lines are consistent with secondary and axial flow distributions.

We take a curved rectangular duct of aspect ratio 2 with curvature $0.1$ and rotate it around the center of the curvature with an angular velocity ${{\Omega }_{T}}$ in the positive direction. Positive direction means that the rotational direction is in the same as the main flow direction. In the following, time evolutions of the unsteady flow characteristics are discussed, in detail, for $Dn=1000$, $Gr=100$ and $0\le Tr\le 2000$.

## 9.png

Figure 9. Unsteady solutions for $Dn=1000$, $Gr=100$ and $Tr=0$. (a) Time evolution of $\lambda$, (b) Phase plots in the $\lambda -\gamma$ plane, (c) Power spectra of the time evolution

## 10.png

Figure 10. Streamlines of secondary flow (top), axial flow (middle) and isotherms (bottom) for $Dn=1000$, $Gr=100$, $Tr=0$and $14.00\le t\le 18.00$

We performed time evolution calculation of $\lambda$ for $0\le Tr\le 1470$ at $Gr=100$. Figures 9(a), 11(a), 13(a) and 15(a) show the time evolution results for $Tr=0$, $Tr=1000$, $Tr=1200$ and $Tr=1470$ respectively. As seen in Figures 9(a), 11(a), 13(a) and 15(a), time-dependent solutions for $Tr=0$, $Tr=1000$, $Tr=1200$ and $Tr=1470$ oscillate irregularly that means the flow is chaotic and it is found that the flow is chaotic for all values of $Tr$ in the range $0\le Tr\le 1470$.

## 11.png

Figure 11. Unsteady solutions for $Dn=1000$, $Gr=100$ and $Tr=1000$. (a) Time evolution of $\lambda$, (b) Phase plots in the $\lambda -\gamma$ plane, (c) Power spectra of the time evolution

## 12.png

Figure 12. Streamlines of secondary flow (top), axial flow (middle) and isotherms (bottom) for , $Gr=100$, $Tr=1000$and $12.00\le t\le 17.00$

To clearly observe the unsteady flow characteristics, we also draw phase space of the time evolution results as shown in Figures 9(b), 11(b), 13(b) and 15(b) for $Tr=0$, $Tr=1000$, $Tr=1200$ and $Tr=1470$ respectively in the $\lambda -\gamma$ plane, where $\gamma =\int{\int{\psi dxdy}}$ and this quantity is zero at the cross-section. As seen in Figs. 9(b), 11(b), 13(b) and 15(b), most of the $\lambda -\gamma$ plane is covered with chaotic orbits, which shows that the unsteady flow for $Tr=0$, $Tr=1000$, $Tr=1200$ and $Tr=1470$ is chaotic. In order to justify the chaotic solution in more detail, the power spectra of the time change of $\lambda$ are calculated for $0\le Tr\le 1470$ as shown in Figures 9(c), 11(c), 13(c) and 15(c) for $Tr=0$, $Tr=1000$, $Tr=1200$ and $Tr=1470$ respectively, where continuous line spectrum with different frequencies are seen, which justify that the flow is chaotic for $Tr=0$, $Tr=1000$, $Tr=1200$ and $Tr=1470$. Then we draw streamlines of secondary and axial flow and isotherms of temperature profiles for the respective flow parameters as shown in Figures 10, 12, 14 and 16.

## 13.png

Figure 13. Unsteady solutions for $Dn=1000$, $Gr=100$ and $Tr=1200$. (a) Time evolution of $\lambda$, (b) Phase plots in the $\lambda -\gamma$ plane, (c) Power spectra of the time evolution

It is found that the streamlines of the secondary flow consist of two opposite vortices; one is an outward flow (anticlockwise direction) shown by solid line and the other one inward flow (clockwise direction) shown by dotted lines. The flow is accelerated due to combined action of the centrifugal, Coriolis and buoyancy forces; centrifugal force is created due to the motion through a curved channel, Coriolis force due to the rotation of the duct around the vertical axis while buoyancy forces because of the thermal gradient. It is found that the unsteady solutions at $Tr=0$ and $Tr=1000$ oscillates irregularly in the four- to seven-vortex solutions, the unsteady solutions at $Tr=1200$ oscillates in the four- to eight-vortex solutions, while the unsteady solution at $Tr=1470$ oscillates in the two- to six-vortex solutions The unsteady flow at $Tr=0,\text{ }Tr=\text{1000 }\!\!\And\!\!\text{ }Tr=\text{1200}$ are called strong chaos while that for $Tr=1470$ weak chaos [9]. It is found that maximum eight-vortex solution is attained at $Tr=1200$. Temperature distribution is found to be consistent with the secondary vortices and a strong interaction is observed between the heating-induced buoyancy force and the centrifugal instability, which stimulates fluid mixing and thus results in thermal enhancement in the flow. In this study, it is found that secondary flow enhances heat transfer in the flow particularly when Dean vortices emerge at the outer wall. It is also found that combined action of the centrifugal, Coriolis and buoyancy force help to increase the number of secondary vortices, and as the flow becomes chaotic, the number of secondary vortices increases and consequently heat is transferred substantially from the heated bottom wall to the fluid. If the Tr is increased a little, for example, $Tr=1480$, it is found that the chaotic flow turns into periodic.

## 14.png

Figure 14. Streamlines of secondary flow (top), axial flow (middle) and isotherms (bottom) for $Dn=1000$, $Gr=100$, $Tr=1200$and $13.00\le t\le 18.20$

In this study, contours of temperature profiles show that the streamlines of the heat flow is uniformly distributed to all parts of the contour transferring heat from bottom wall to the fluid, and the contribution of the rotation and pressure on secondary flows significantly change and increase the number of secondary vortices. It is clearly evident that heating the bottom wall causes the temperature contours to become asymmetrical in comparison to isothermal cases. This essentially arises from the interaction between the heating-induced buoyancy force and the centrifugal force that drives secondary vortices. In this regard, it should be noted that the centrifugal force due to the duct curvature creates two effects; it generates a positive radial fluid pressure field in the duct cross section and induces a lateral fluid motion driven from inner wall towards the outer wall. This lateral fluid motion occurs against the radial pressure field generated by the centrifugal effect and is superimposed on the axial flow to create the secondary vortex flow structure. As the flow through the curved duct is increased, the lateral fluid motion becomes stronger and the radial pressure field is intensified. In the vicinity of the outer wall, the combined action of adverse radial pressure field and viscous effects slows down the lateral fluid motion and forms a stagnant flow region. Beyond a certain critical value of $Dn$, the radial pressure gradient becomes sufficiently strong to reverse the flow direction of the lateral fluid flow. A weak local flow re-circulation is then established creating an additional pair of vortices in the stagnant region near the outer wall. This flow situation is known as Dean’s hydrodynamic instability while the vortices are termed as Dean vortices.

## 15.png

Figure 15. Unsteady solutions for $Dn=1000$, $Gr=100$ and $Tr=1470$. (a) Time evolution of $\lambda$, (b) Phase plots in the $\lambda -\gamma$ plane, (c) Power spectra of the time evolution.

In order to study the non-linear behavior of the unsteady solution, we investigate time evolution of $\lambda$ for $Tr=1480$ at $Gr=100$ and $Dn=1000$ as shown in Figure 17, where we see that the unsteady solution at $Tr=1480$ oscillates periodically. This periodic oscillation is well justified by drawing the orbits of the solution in the phase space as shown in Figure 17(b), where we see that periodic orbits overlap each other, which confirms that the flow at $Tr=1480$ is periodic. In order to investigate the transition from a chaotic solution to periodic oscillation in more detail, the power spectra of the time change of $\lambda$ is calculated for $Tr=1480$. The result is shown in Figure 17(c), in which not only the line spectrum of the fundamental frequency with large frequency but the other line spectrum with small frequency is seen, which shows that the oscillation presented in Figure 17(a) is periodic.

## 16.png

Figure 16. Streamlines of secondary flow (top), axial flow (middle) and isotherms (bottom) for $Dn=1000$, $Gr=100$, $Tr=1470$and $10.00\le t\le 15.65$

In fact, the periodic oscillation, which is observed in the present study, is a traveling wave solution advancing in the downstream direction which is well-justified in the recent investigation by Yanase et al. [34] for three-dimensional travelling wave solutions as an appearance of 2D periodic oscillation. Therefore, it is found that 2D calculations can accurately predict the existence of 3D traveling wave solutions by showing an appearance of 2D periodic oscillation. Then we draw streamlines of secondary and axial flow and isotherms of temperature profile for $Tr=1480$ as shown in Figure 18. It is found that the streamlines of secondary flow for the periodic oscillation at $Tr=1480$ consist of asymmetric two-vortex solution. It is found that the transition from chaotic state to the periodic oscillation occurs between $Tr=1470$and $Tr=1480$. If $Tr$ is increased further, for example $Tr=1490$, it is found that the periodic oscillation turns into steady-state solution as shown in Figure 19.

In order to search for the region of steady-state solutions, we investigate time evolution of $\lambda$for $1490\le Tr\le 2000$ and it is found that the flow is steady-state for all the values of $Tr$ in the range. Figures 19(a) and 20(a) show steady-state solutions for $Tr=1490$ and $Tr=2000$ respectively. Typical contours of secondary flow patterns and temperature profiles are shown in Figure 19(b) for $Tr=1490$ and in Figure 20(b) for $Tr=2000$ at time $t=20.00$. It is found that unsteady solution at $Tr=1490$ and $Tr=2000$ possess two-vortex solutions. In this study, it is found that if the $Tr$ is increased further, for example, $Tr>2000$, the flow is steady-state as well. Therefore, it is suggested that the transition from periodic solution to the steady-state occurs between $Tr=1480$ and $Tr=1490$. It is also found that temperature distribution is consistent with secondary vortices, and secondary flow enhances heat transfer in the flow through vortex generation.

## 17.png

Figure 17. Unsteady solutions for $Dn=1000$, $Gr=100$ and $Tr=1480$. (a) Time evolution of $\lambda$, (b) Phase plots in the $\lambda -\gamma$ plane, (c) Power spectra of the time evolution

## 18.png

Figure 18. Streamlines of secondary flow (top), axial flow (middle) and isotherms (bottom) for $Dn=1000$, $Gr=100$, $Tr=1480$and $17.10\le t\le 17.40$

## 19.png

Figure 19. (a) Time evolution of $\lambda$for $Dn=1000$, $Gr=100$ and $Tr=1490$, (b) Streamlines of secondary flow (top) and isotherm (bottom) for $Tr=1490$ at $t=20.00$

4.3 Validation of the numerical result

Here, we represent the validation of our numerical results with the experimental studies performed by some authors. By using visualization method, Yamamoto et al. [30] conducted experimental investigations (Figure 21(a)) of the flow through a rotating curved square duct of curvature $\delta =0.03$, where three of the duct walls, except the outer wall, rotate around the centerline of curvature at a constant revolution speed for the positive rotation at Tr =150.

In this study, however, we investigate flow features for the rotation of the whole system (not the three walls only), and compare our results with that of Yamamoto et al. [30] considering the same curvature and rotational speed. On the other hand, Figure 21(b) shows a comparative study of our numerical result with the experimental investigation obtained by Chandratilleke [35] for the flow through a curved rectangular duct of aspect ratio 2. We see that in both the cases our numerical results have a good agreement with the experimental data. Note that, till now no experimental studies have been found for rotating curved rectangular duct flow.

## 20.png

Figure 20. (a) Time evolution of $\lambda$ for $Dn=1000$, $Gr=100$ and $Tr=2000$, (b) Streamlines of secondary flow (top) and isotherm (bottom) for $Tr=2000$ at $t=20.00$

## 21.png

Figure 21. Experimental vs. Numerical results; left: experimental results, right: numerical results; (a) curved square duct, (b) curved rectangular duct

5. Conclusions

Due to ample applications from medical services to industrial activities, the study of flows and heat transfer through a rotating curved rectangular duct has attracted much attention to the researchers. In this study, a comprehensive numerical study on fluid flow and heat transfer through a rotating curved rectangular channel of aspect ratio 2 has been presented for a constant curvature $\delta =0.1$ over the Taylor number $0\le Tr\le 2000$ for $Dn=1000$. A temperature difference is applied between the horizontal walls for $Gr=100$, where the bottom wall is heated and cooling from the ceiling. After an extensive survey over the parametric ranges, three branches of asymmetric steady solutions with two- to ten-vortex solutions are obtained. It is found that there is no bifurcating relationship among the branches of steady solutions. We then investigated unsteady flow behavior by time evolution calculation justified by drawing the phase spaces and power spectra of the solutions. It is found that the unsteady flow undergoes in the scenario chaotic$\to$periodic$\to$steady-state’, it $Tr$ is increased in the positive direction. The study shows that combined action of the centrifugal and Coriolis force counteracts each other in a nonlinear manner which results in to turn the chaotic flow to steady-state flow. It is found that the unsteady flow is chaotic for $0\le Tr\le 1470$, periodic for $Tr=1480$and steady-state for $1490\le Tr\le 2000$. It is found that the transition from chaotic state to periodic flow occurs between $Tr=1470$ and $Tr=1480$; and from periodic to steady-state between $Tr=1480$ and $Tr=1490$. The present study elucidates the role of secondary vortices on convective heat transfer and it is found that convective heat transfer is significantly enhanced by the secondary flow; and the chaotic flow, which occurs at small $Tr$ but at large $Dn$, enhances heat transfer more effectively than the steady-state or periodic solutions. In this regard, it should be noted that, the occurrence of the chaotic state, as presented in the present study, is related with destabilization of the periodic or quasi-periodic solutions which reminds us the case of Lorenz attractor [36]. It may be possible that the transition in the present study is caused by a similar mechanism as that of Ruelle-Takens scenario [37] in the laminar flow.

Acknowledgement

Rabindra Nath Mondal, one of the authors, would gratefully acknowledge the financial support from the Bangladesh Ministry of Education (MoEdu) for advanced research in Science and Technology (37.20.0000.004.033.005.2014-1309/1(42)) to conduct this research work.

Nomenclature
 Dn Dean number Gr Grashof number Pr Prandtl number A Aspect ratio L Radius of the curvature X Horizontal axis Y Vertical axis Z Axis in the direction of the main flow U Velocity components in the $x-$direction V Velocity components in the $y-$direction W Velocity components in the $z-$direction T Temperature T Time Greek symbols $\delta$ Curvature of the duct $\rho$ Density $\lambda$ Resistance coefficient $\mu$ Viscosity $\kappa$ Thermal diffusivity $\upsilon$ Kinematic viscosity $\psi$ Sectional stream function
References

[1] Dierich, F., Nikrityuk, P.A. (2013). A numerical study of the impact of surface roughness on heat and fluid flow past a cylindrical particle. International Journal of Thermal Sciences, 65: 92-103. https://doi.org/10.1016/j.ijthermalsci.2012.08.009

[2] Chatterjee, S., Sugilal, G., Prabhu, S.V. (2018). Heat transfer in a partially filled rotating pipe with single phase flow. International Journal of Thermal Sciences, 125: 132-141. https://doi.org/10.1016/j.ijthermalsci.2017.11.024

[3] Khaled, A.M., Yongkai, Q., Wei, G., Guoqiang, X., Yan, Y. (2017). Contrastive study of flow and heat transfer characteristics in a helically coiled tube under uniform heating and one-side heating. Applied Thermal Engineering, 114: 77-84. https://doi.org/10.1016/j.applthermaleng.2016.11.168

[4] Lin, W.C., Feng, Y.M., Chieng, C.C. (2013). Numerical computations on flow and heat transfer characteristics of a helically coiled heat exchanger using different turbulence models. Nuclear Engineering and Design, 263: 77-86. https://doi.org/10.1016/j.nucengdes.2013.03.051

[5] Dean, W.R. (1927). Note on the motion of fluid in a curved pipe. Philos Mag, 4: 208-23. https://doi.org/10.1080/14786440708564324

[6] Berger, S.A., Talbot, L., Yao, L.S. (1983). Flow in curved pipes. Annual Review of Fluid Mechanics, 35: 461-512. https://doi.org/10.1146/annurev.fl.15.010183.002333

[7] Mondal, R.N., Islam, S., Uddin, K., Hossain, M.A. (2013). Effects of aspect ratio on unsteady solutions through curved duct flow. Applied Mathematics and Mechanics, 34: 1107-1122. https://doi.org/10.1007/s10483-013-1731-8

[8] Mondal, R.N., Watanabe, T., Hossain, M.A., Yanase, S. (2017). Vortex-structure and unsteady solutions with convective heat transfer through a curved duct. Journal of Thermophysics and Heat Transfer, 31(1): 243-254. https://doi.org/10.2514/1.T4913

[9] Chen, K.T., Yarn, K.F., Chen, H.Y., Tsai, C.C., Luo, W.J., Chen, C.N. (2015). Aspect ratio effect on multiple flow solutions in a two-sided parallel motion lid-driven cavity. Journal of Mechanics, 31(2): 153-160. https://doi.org/10.1017/jmech.2014.51

[10] Drazin, P.G., Reid, W.H. (2004). Hydrodynamic stability. Cambridge University Press, Cambridge. https://doi.org/10.1017/CBO9780511616938

[11] Sugiyama, S., Hayashi, T., Yamazaki, K. (1983). Flow characteristics in the curved rectangular channels: Visualization of secondary flow. Bulletin of JSME, 26(216): 964-969. https://doi.org/10.1299/jsme1958.26.964

[12] Winters, K.H. (1987). A bifurcation study of laminar flow in a curved tube rectangular cross section. Journal of Fluid Mechanics, 180: 343-369. https://doi.org/10.1017/S0022112087001848

[13] Yanase, S., Kaga, Y., Daikai, R. (2002). Laminar flow through a curved rectangular duct over a wide range of aspect ratio. Fluid Dynamics Research, 31(3): 151-183. https://doi.org/10.1016/S0169-5983(02)00103-X

[14] Hasan, M.S., Mondal, R.N., Lorenzini, G. (2019). Numerical prediction of non-isothermal flow with convective heat transfer through a rotating curved square channel with bottom wall heating and cooling from the ceiling. International Journal of Heat and Technology, 37(3): 710-726. https://doi.org/10.18280/ijht.370307

[15] Chandratilleke, T.T., Nursubyakto. (2003). Numerical prediction of secondary flow and convective heat transfer in externally heated curved rectangular ducts. International Journal Thermal Science, 42(2): 187-198. https://doi.org/10.1016/S1290-0729(02)00018-2

[16] Fellouah, H., Castelain, C., Moctar, A.O.E., Peerhossaini, H.A. (2006). Criterion for detection of the onset of Dean instability in Newtonian fluids. European Journal of Mechanics, B-Fluid, 25: 505-531. https://doi.org/10.1016/j.euromechflu.2005.11.002

[17] Zhang, J., Zhang, B., Jii, J. (2001). Fluid flow in a rotating curved rectangular duct. International Journal of Heat and Fluid Flow, 22(6): 583-592. https://doi.org/10.1016/S0142-727X(01)00126-6

[18] Mondal, R.N., Kaga, Y., Hyakutake, T., Yanase, S. (2007). Unsteady solutions and the bifurcation diagram for the flow through a curved square duct. Fluids Dynamics Research, 39(5): 413-446. https://doi.org/10.1016/j.fluiddyn.2006.10.001

[19] Wang, L., Yang, T. (2005). Periodic oscillation in curved duct flows. physica D, 200(3-4): 296-302. https://doi.org/10.1016/j.physd.2004.11.003

[20] Yanase, S., Nishiyama, K. (1988). On the bifurcation of laminar flows through a curved rectangular tube. Journal of the Physical Society of Japan, 57(11): 3790-3795. https://doi.org/10.1143/JPSJ.57.3790

[21] Hasan, M.S., Mondal, R.N., Lorenzini, G. (2019). Centrifugal instability with convective heat transfer through a tightly coiled square duct. Mathematical Modelling of Engineering Problems, 6(3): 397-408. https://doi.org/10.18280/mmep.060311

[22] Wang, L.Q., Liu, F. (2007). Forced convection tightly coiled ducts: Bifurcation in a high Dean number region. International Journal of Non-linear Mechanics, 42(8): 1018-1034. https://doi.org/10.1016/j.ijnonlinmec.2007.05.005

[23] Islam, M.Z., Mondal, R.N., Rashidi, M.M. (2017). Dean-Taylor flow with convective heat transfer through a coiled duct. Computers and Fluids, 149: 41-55. https://doi.org/10.1016/j.compfluid.2017.03.001

[24] Norouzi, M., Kayhani, M.H., Shu, C., Nobari, M.R.H. (2010). Flow of second-order fluid in a curved duct with square cross section. Journal of Non-Newtonian Fluid Mechanics, 165(7-8): 323-339. https://doi.org/10.1016/j.jnnfm.2010.01.007

[25] Chandratilleke, T.T., Nadim, N., Narayanaswamy, R. (2012). Vortex structure-based analysis of laminar flow behavior and thermal characteristics in curved ducts. International Journal of Thermal Sciences, 59: 75-86. https://doi.org/10.1016/j.ijthermalsci.2012.04.014

[26] Zhang, L.J., Zhang, W.P., Ming, P.J. (2014). Numerical simulation of secondary flow in a curved square duct. Applied Mathematics and Materials, 681: 33-40. https://doi.org/10.4028/www.scientific.net/AMM.681.33

[27] Hasan, M.S., Mondal, R.N., Kouchi, T., Yanase, S. (2017). Hydrodynamic instability with convective heat transfer through a curved channel with strong rotational speed. AIP Conference Proceedings, 2121(1): 030006. https://doi.org/10.1063/1.5115851

[28] Wu, X.Y., Lai, S.D., Yamamoto, K., Yanase, S. (2013). Numerical analysis of the flow in a curved duct. Advanced Materials Research, 706-708: 1450-1453. https://doi.org/10.4028/www.scientific.net/AMR.706-708.1450

[29] Mondal, R.N., Islam, M.Z., Pervin, R. (2014). Combined effects of centrifugal and coriolis instability of the flow through a rotating curved duct of small curvature. Procedia Engineering, 90: 261-267. https://doi.org/10.1016/j.proeng.2014.11.847

[30] Yamamoto, K., Wu, X., Kazuo, N., Yasutaka, H. (2006). Visualization of Taylor-Dean flow in a curved duct of square cross-section. Journal Fluid dynamics research, 38(1): 1-18. https://doi.org/10.1016/j.fluiddyn.2005.09.002

[31] Razavi, R., Soltanipour, S.E., Choupani, P. (2015). Second law analysis of laminar forced convection in a rotating curved duct. Thermal Science, 19(1): 95-107. https://doi.org/10.2298/tsci120606034r

[32] Li, Y., Wang, X., Yuan, S., Tan, S.K. (2016). Flow development in curved rectangular ducts with continuously varying curvature. Experimental Thermal and Fluid Science, 75: 1-15. https://doi.org/10.1016/j.expthermflusci.2016.01.012

[33] Gottlieb, D., Orazag, S.A. (1977). Numerical analysis of spectral methods. Society of Industrial and Applied Mathematics, Philadelphia, USA. https://doi.org/10.1137/1.9781611970425

[34] Yanase, S., Watanabe, T., Hyakutake, T. (2008). Traveling-wave solutions of the flow in a curved-square duct. Physics of Fluids, 20(12): 1-8. https://doi.org/10.1063/1.3029703

[35] Chandratilleke, T.T. (2001), Secondary flow characteristics and convective heat transfer in a curved rectangular duct with external heating. 5th World Conference on Experimental Heat Transfer, Fluid Mechanics and Thermodynamics [ExHFT-5], Thessaloniki, Greece.

[36] Lorenz, E.N. (1963). Deterministic non-periodic flow. J. Atmos. Sci, 20: 130-141. https://doi.org/10.1007/978-0-387-21830-4_2

[37] Ruelle, D., Takens, F. (1971). On the nature of turbulence. Communications in Mathematical Physics, 20: 167-192. https://doi.org/10.1007/BF01646553