This paper deals with convective instability in a fluid-saturated, rotating porous layer subject to alternating direction of the centrifugal body force, when the layer fails to exhibit thermal equilibrium. The Darcy model is used to describe the flow, and a two-field model is used to take care of the energy exchange. The normal mode approach of the linear theory and the energy approach of the nonlinear theory are used to find the stability characteristics. Unconditional and sharp nonlinear thresholds are found. The study brings out the failure of the linear theory in describing the instability in most parts of the parameter space of interest where possible subcritical instabilities may arise. The stability boundaries are presented graphically and it is found that the inter-phase heat transfer coefficient has a significant effect in the thermal non-equilibrium regime.
The dispersion phenomena through porous media coupled with heat transfer exhibit their presence in a wide variety of applications. The extensive work carried out so far in this area has been well reviewed in the book of Nield & Bejan . This paper considers heat transfer in rotating porous media mainly motivated by applications such as food and chemical processing, solidification and centrifugal casting of metals, rotating turbo machineries, etc. In such applications, one important observation is the following: at high rotation speeds in the terrestrial environment, the resulting centrifugal acceleration can induce buoyancy that could dominate the gravity-induced one. Convection caused by this type of buoyancy, referred to as centrifugal convection, is the only possibility if a system is rotating in zero-gravity condition. Unlike gravity, the centrifugal acceleration can act in any direction and varies with position. Centrifugal convection has received less attention when compared with the gravitational counterpart, though there are many promising applications. The initial work of Vadasz  dealt with the onset of centrifugal convection in a Darcian porous medium using linear stability analysis. The conditions imposed at the boundaries were in such a way that the resulting temperature gradient remains collinear with the centrifugal body force. Later, Vadasz  investigated the effect of including alternating direction of the centrifugal body force. These studies accounted for the effect of the centrifugal body force while neglecting gravity.
Although linear theory does provide a useful result, in order to have a complete understanding, one has to perform a nonlinear analysis that provides a threshold for global stability. The nonlinear result, obtained through the energy approach , is certainly much stronger when the stability obtained is unconditional, i.e. for all initial data or for at least finite (non-vanishing) initial data. By introducing a coupling parameter in the energy approach and by selecting it optimally, it is possible to sharpen the stability bound in many physical problems. The energy approach has been applied to rotating fluid problems, and only conditional results have been derived. However, unconditional nonlinear results have been obtained for rotating porous systems. Straughan  derived nonlinear energy results for thermal convection in a Darcian porous layer which is rotating about an axis orthogonal to the planes containing the layer. Recently, we  obtained nonlinear results for the problem discussed by Vadasz [1,2].
In most of the investigations involving porous media, the state of local thermal equilibrium (LTE) has been profoundly used. Nevertheless, in many practical applications involving sudden and high speed flows, the hot fluid stream penetrates well into the relatively cold porous structures, and hence, in a representative elementary volume, its temperature becomes sufficiently higher than that of the adjacent solid phase. This non-local thermal non-equilibrium (NLTE) situation can be handled by considering separate energy equations for solid and fluid phases with a coupling in between them to represent the energy exchange. A brief review of various two equation models accounting for the NLTE effect in porous convection and many recent studies including free and external forced convection boundary layers were given by Rees & Pop . The onset criterion for convection in a porous medium exhibiting NLTE was first determined by Banu & Rees . They predicted the emergence of tall thin convection cells under circumstances that depend on the conductivity ratio of the two phases. This was followed by similar studies in the Brinkman porous medium confined between stress-free [9,10] and rigid  boundaries. These studies on NLTE were performed based on linear stability analysis, whereas Straughan  carried out a nonlinear analysis and concluded that the linear results are important as they coincide with the nonlinear limits.
The later studies on the NLTE situation have included some additional effects. These include double diffusive convection [13–15], ferroconvection [16,17], convection in complex fluids [18,19] and non-uniformly heated fluids  and thermovibrational convection . Recent articles have also discussed this situation in vertical porous media [22,23]. Apart from the earlier-mentioned convective instability problems, recent studies on NLTE have also addressed free convection from a vertical plate , convection through a cylindrical annulus , free convection in a vertical layer heated through a heat flux , an effective porosity concept applicable to NLTE models , pore pressure and thermal stress distributions in a porous medium subjected to convective cooling/heating on its boundary under transient thermal loads  and mixed convection in a vertical channel adopting the non-Darcy–Brinkman–Forchheimer extended model . Motivated by the growing volume of studies on thermal non-equilibrium in recent times, this paper aims to carry out linear as well as nonlinear analyses for the onset of convecton in a vertical rotating porous layer that lacks thermal equilibrium.
2. Problem formulation
We consider a fluid-saturated, infinitely tall vertical porous layer −L/2<x<L/2 subject to a constant rotation rate ω about an axis in a zero-gravity environment (figure 1). The porous layer is saturated with a permeable fluid. The layer is heated on its right boundary (Th) and cooled on the left boundary (Tc), and as a result of these imposed thermal boundary conditions, a uniform temperature gradient β is acting across the layer. The axis of rotation is parallel to the z-axis and is placed anywhere within the boundaries of the porous domain at a dimensionless distance x0 from the centre of the layer. The offset distance is presented in a dimensionless form, representing the ratio between the dimensional offset distance and the thickness of the porous layer in the form . Free convection occurs as a result of the centrifugal body force that can have alternating directions within the layer. The Coriolis force is considered small, because most of the inertial forces are neglected when using Darcy's law. The only inertial effect considered is the centrifugal acceleration, as far as the changes in density are concerned. These assumptions are compatible with Darcy's law, which is applicable as long as the Reynolds number is not too large. The thickness of the layer W is assumed to obey W≪L, which enables us to neglect the y component of the centrifugal acceleration and hence to use a Cartesian coordinate system. The fluid is assumed to be viscous and incompressible, and the Boussinesq approximation is used to account for the effects of the density variations.
The conservation equations for the above set up in the rotating frame of reference take the form 2.1where is the filtration velocity, Kp is the permeability of the porous material, p is the pressure, ρ is the fluid density, μ is the dynamic viscosity and is the unit vector in the x-direction. Here, (2.1)1 represents the mass balance for an incompressible fluid. Equation (2.1)2 represents the momentum balance expressed in terms of an extended form of Darcy's law, including the centrifugal rotation effect. One should note that it has taken into account only the x-component of the centrifugal force. A two-field model that represents the temperature fields of the fluid and solid phases separately is used for the energy balance: 2.2Here, T is the temperature, ϕ is the porosity, t is the time, CP is the specific heat capacity at constant pressure, k is the thermal diffusivity and h is the inter-phase heat transfer coefficient. The subscripts ‘f’ and ‘s’ refer to the fluid and the solid phases, respectively. In this model, the energy equations are coupled through the last term, which account for the heat lost to or gained from the other phase. The inter-phase heat transfer coefficient, h, depends on the nature of the porous matrix and the saturating fluid, and the value of this coefficient has been the subject of intense experimental interest. When h=0, the energy equations become decoupled and take the form of a single equation as reported in Saravanan & Brindha  under LTE assumption. The equation of state is given by 2.3where ρ0 is the reference density, T0=(Tc+Th)/2 is the average temperature and α is the coefficient of thermal expansion of the fluid.
Under the earlier-mentioned set-up, the governing equations admit an equilibrium state described by 2.4We now consider where the primed quantities denote the disturbances on the corresponding terms. Upon substitution of (2.4) into (2.1) and (2.2), the governing equations for the disturbances may be written as 2.5where the primes have been omitted for convenience. These equations are then non-dimensionalized, using the scales L for length, L2/kf for time, ϕkf/L for velocity, βLTf and βLTs for temperatures of fluid and solid phases, respectively, and μϕkf/Kp for pressure. Then, the perturbed non-dimensional equations are 2.6where Rc=αβρ0ω2KpL3/(ϕμkf), H=hL2/(ϕkf(ρCP)f), γ=ϕkf(ρCP)f/((1−ϕ)ks(ρCP)s), χ=kf/ks are the centrifugal Rayleigh number based on the fluid properties, the scaled inter-phase heat transfer coefficient, the porosity-modified conductivity ratio and the diffusivity ratio, respectively.
Because the fluid and solid phases are not in thermal equilibrium, the use of appropriate thermal boundary conditions may pose a difficulty. However, the assumption that the solid and fluid phases share the same temperatures at the boundary helps in overcoming this difficulty. Accordingly, the boundary conditions for the perturbation variables are given by 2.7These correspond to no-penetration and isothermal conditions on the physical boundaries.
3. Linear stability analysis
In order to find linear instability threshold of the equilibrium state, we shall follow the normal mode technique  whereby disturbances are resolved into a complete set of normal modes, each of which may be treated separately, to analyse the development of an arbitrary initial disturbance. Accordingly, we assume a temporal growth of disturbances in the form 3.1where σ is a complex constant representing the eigenvalue of each mode. The sign of the real part of σ will decide the stability of the flow. It is important to note that the linear analysis approach assumes the disturbances to be small enough so that terms of quadratic and higher order can be neglected. Hence, the resulting system, obtained by substituting (3.1) into (2.6), is 3.2We shall now prove that over-stable motion does not arise for this model. Applying curl twice on (3.2)2 and decomposing the small disturbances in the resulting equation via the normal mode , the linear instability equations (3.2) then reduce to 3.3where k is the wavenumber. We follow the moment method of Mikaelian  to prove our assertion. The choice of the moment method is solely due to the presence of the spatial variable x in (3.3). This method is superior in the sense that it does not suffer from the ambiguities of satisfying some higher-order boundary conditions. Accordingly, we multiply (3.3)1 by Um, (3.3)2 by Θmf (3.3)3 by Θms and integrate over the layer to find 3.4where 〈⋅〉 is the usual integration with respect to x from to . Here, m=0 in (3.4) corresponds to moment equations. Equations (3.4) can be combined into a single equation 3.5where F=〈(x−x0)Θf〉/〈Θf〉 is real. Now letting and equating the real and imaginary parts separately, we arrive at 3.6which, in turn, imply 3.7a negative definite quantity. This shows the validity of the principle of exchange of stabilities, and hence the possibility of over-stable motions is ruled out.
We now analyse the marginal state governed by (3.3) with σ=0. Let us denote the lowest eigenvalue by RcL. Thus, the relevant equations are 3.8subjected to the conditions 3.9Equations (3.8) and (3.9) constitute an eigenvalue system of equations for the linear centrifugal Rayleigh number RcL. In the numerical calculations, we determine the critical centrifugal Rayleigh number of linear instability theory as 3.10
4. Nonlinear stability analysis
A linear analysis can show only where the system is definitely unstable. It may be noted that the solution can become unstable well before the threshold predicted by the linear theory. Hence, it is highly useful to perform a nonlinear stability analysis that may yield results close to those of the linear theory. We now use the energy approach, successfully applied to many convection models , to investigate the nonlinear stability of the diffusive equilibrium state.
Let Ω be a typical period cell for the disturbance, and let denote the integration over Ω and let ∥⋅∥ denote the norm on L2(Ω). We multiply (2.6)2 by and integrate over Ω to find 4.1Similarly, multiplying (2.6)3 by Tf and integrating over Ω, one can find 4.2Multiplying (2.6)4 by Ts and integrating over Ω, we obtain 4.3In order to study the nonlinear stability of the equilibrium state, we now introduce positive coupling parameters λ1 and λ2 and define an energy functional E(t) by 4.4This definition is a generalization of that introduced in Saravanan & Brindha  through an additional coupling parameter λ2. The idea behind the introduction of the coupling parameters is to pick these in an optimal way to derive as sharp a stability boundary as possible. Now, using (4.2)–(4.4), the evolution of the energy functional is derived as 4.5where the production term and the dissipation term are given by 4.6One may note that the expression for agrees with that in Saravanan & Brindha  for H=λ2=0. Here, in order to guarantee the positive definiteness of , the last term on the right-hand side of must be positive. We use the Young inequality with some δ>0 to show that 4.7The positive definiteness of can be guaranteed by choosing the free parameter δ such that λ1−(λ1+λ2)/2δ≥0 and λ2−(λ1+λ2)δ/2≥0.
To establish a nonlinear stability result, we must investigate the maximum problem 4.8where is the space of admissible functions over which we seek a maximum such that , Tf,Ts∈H1(Ω) and on . This equation leads to the following relation 4.9Thus, on integration it follows, with the aid of the Poincare inequality, that if Rc<RcN then as at least exponentially in time. We note that the decay of E is for all initial disturbances, regardless of how large they may be and, hence the nonlinear stability is unconditional.
The nonlinear stability threshold is now given by the variational problem (4.8). The Euler–Lagrange equations arising from (4.8) are given by 4.10This leads to the Euler–Lagrange equations 4.11with the appropriate boundary conditions where π is a Lagrange multiplier.
In order to actually find the sharp limits, we solve these Euler–Lagrange equations for RcN. After removing π by operating curl, (4.11) may then be reduced to 4.12where ux(x,z)=U(x)eikz, , Ts(x,z)=Θs(x)eikz in which k is the wavenumber. The relevant boundary conditions are 4.13The critical centrifugal Rayleigh number requires the determination of the lowest eigenvalue RcN(k,λ1,λ2). At the same time, we exploit the coupling parameters to make RcN as large as possible. This is achieved by performing the optimization 4.14The number RcN,cr so found is the critical centrifugal Rayleigh number predicted by the nonlinear theory.
5. Numerical methodology
The compound matrix method [30,33], which was used by us in our earlier works [6,17], was used to solve the eigenvalue problems obtained in §§2 and 3. We shall now outline the procedure to solve the eigenvalue problem constituted by (4.12) and (4.13).
We first rewrite the eigenvalue equations as a system of first-order equations by defining . Suppose U1, U2 and U3 are respectively the three independent solutions of the system (4.12) and (4.13) with the initial values (0,1,0,0,0,0)T, (0,0,0,1,0,0)T and (0,0,0,0,0,1)T at . We then define 6C3 new variables y1–y20 as the 3×3 minors of the 6×3 solution matrix whose first, second and third columns are U1, U2 and U3, respectively. For example, we define The idea is to define y2–y20 similarly and then obtain differential equations for the yi by differentiation. By differentiating each yi in turn and simplifying, we arrive at the following differential equations for the yis: 5.1
From the initial conditions on Ui, we see that the system (5.1) has to be integrated numerically from 0 to 1 with the initial condition 5.2The appropriate final condition that satisfies (4.13) is found to be 5.3The eigenvalue RcN is varied until (5.3) is satisfied.
6. Results and conclusion
Thermal instability in a porous medium in the form of centrifugal convection is investigated when the medium moves away from the state of thermal equilibrium. The compound matrix method was used to solve the resulting eigenvalue problems. The stability limits obtained through the two different analyses discussed earlier are presented graphically as functions of the location of the axis of rotation x0, the inter-phase heat transfer coefficient H and the porosity modified conductivity ratio γ.
The linear and nonlinear marginal stability curves are displayed in figure 2. They clearly show that both the analyses predict different results. These curves clearly partition the Rc−k plane into two disjoint stable and unstable regions. Each of these curves has a single minimum called the critical Rayleigh number Rccr, accompanied by a corresponding critical wavenumber kcr, at which the instability in the form of a series of counter-rotating convective cells near the hot wall sets in Vadasz . This is anticipated as the centrifugal acceleration has a destabilizing effect to the right of the axis of rotation. The effect of H on the marginal stability curves is shown in figure 2a when the axis of rotation lies at the midplane of the porous layer (i.e. x0=0). It is observed that the marginal curves move upward indicating that the effect of H is to stabilize the system. It is found from figure 2b that an increase in the values of γ leads to a destabilizing effect.
The effect of H on the stability characteristics when the axis of rotation is placed at the left boundary (cold wall) is shown in figure 3. An important observation is the coincidence of the nonlinear boundaries almost with the linear ones. From this, one can conclude that the linear analysis itself captures the physics of the problem and hence becomes important. The result also demonstrates that the coupling parameters can be optimized to make the nonlinear limit sharp enough. We may also infer from figure 3a that the centrifugal acceleration can also be a physical mechanism that can cause the proximity of RcN,cr to RcL,cr such as gravity and thermocapillarity .
We observe from figure 3 that the stability characteristics change significantly only for intermediate values of H, which in turn delays the onset of convection. This is because a strong non-equilibrium effect is achieved between the two extremes and . We observe that for very small values of H, Rccr is independent of γ. This is because the conductivities do not play any role when there is no significant exchange of heat between the two phases. As H starts increasing the two phases start interacting, depending on the conductivities of both phases of the medium. Hence, Rccr is found to be a function of γ for higher values of H. Moreover, we see that for moderate and large values of H, Rccr decreases with increasing values of γ. Therefore, the effect of γ is to reduce the stabilizing effect of H.
From figure 3, one can also observe that kcr increases with H, reaches a maximum value and then decreases back for a further increase in H. We find that kcr approaches π, the LTE limit [2,6] when and . This is because when , the solid and fluid phases do not interact at all, whereas when , the solid and fluid phases have identical temperatures, and hence may be treated as a single phase. Thus, a strong non-equilibrium effect is seen only in the intermediate H range and the instability manifests in the form of flat circulation patterns. Figure 3c depicts the variation of Rccr(γ/(1+γ)) against H. It should be noted that Rccr(γ/(1+γ)) is the critical centrifugal Rayleigh number based on the mean properties of the porous medium which is used when assuming the LTE assumption. It is found to increase monotonically depending on γ. It is of interest to note that as , Rccr(γ/(1+γ)) of both linear and nonlinear analyses reach the well-known LTE limit 77.10 [2,6].
The influence of location of the rotating axis is displayed in figure 4. These results correspond to γ=0.01 and H=100. It is observed that RcL,cr increases monotonically against x0 and becomes unbounded as , representing unconditional stability from the linear theoretical point of view. Nevertheless, RcN,cr agrees satisfactorily well with the linear one for x0<0.2. It starts deviating from the linear limit as x0 increases beyond 0.2. This is natural because the linear theory guarantees stability for laminar flows Rc<RcL,cr only in the conditional sense, for small disturbances. In fact, RcN,cr always remains finite for all values of x0. Thus, a clear instability region of subcritical nature can be seen. This indicates that one should use the nonlinear theory to predict reasonable results as the axis of rotation moves towards the hot wall. The corresponding change in kcr as a function of the moving axis is shown in figure 4b. The nonlinear results show that convection always sets in with a finite wavenumber in contrast to the linear theory, which predicts the onset of convection with an unbounded wavenumber as .
In general, though the linear and nonlinear limits suffer qualitatively same changes against the thermal non-equilibrium parameters H and γ, the linear analysis is found to over predict the onset condition compared with the nonlinear analysis. Thus, instabilities may arise before the linear threshold is reached. We note that this result is unconditional and delimits the region of subcritical bifurcations. This is quite understandable, because the linear stability theory gives sufficient conditions for instability, whereas the nonlinear one gives the sufficient conditions for stability. Thus, we can conclude that the energy method produces practically useful optimal results that cannot be determined by the linear theory. We would like to add a final remark that the entire analysis is also applicable for any orientation of the porous layer as there is no gravity.
The authors thank the University Grants Commission, India for its support through DRS Special Assistance Programme in Fluid Dynamics. This work was carried out as a part of a research project (No. SR/FTP/MS-02/2007) awarded by SERC, Department of Science and Technology, India.
- Received November 2, 2012.
- Accepted March 20, 2013.
- © 2013 The Author(s) Published by the Royal Society. All rights reserved.