Finite amplitude convection and heat transfer in inclined porous layer using a thermal non-equilibrium model

International Journal of Heat and Mass Transfer 113 (2017) 399–410

Contents lists available at ScienceDirect

International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt

Finite amplitude convection and heat transfer in inclined porous layer using a thermal non-equilibrium model M.N. Ouarzazi a,⇑, S.C. Hirata a, A. Barletta b, Michele Celli b a b

Laboratoire de Mécanique de Lille, CNRS, FRE 3723, Université de Lille I, Bd Paul Langevin, 59655 Villeneuve d’Asq, France Department of Industrial Engineering, Alma Mater Studiorum Università di Bologna, Viale Risorgimento 2, Bologna 40136, Italy

a r t i c l e

i n f o

Article history: Received 4 December 2016 Received in revised form 27 March 2017 Accepted 22 May 2017

Keywords: Porous medium Darcy’s law Nonlinear stability Inclined layer Local thermal non-equilibrium

a b s t r a c t Finite amplitude convection in a inclined porous layer heated from below is studied by using local thermal non-equilibrium (LTNE) as mathematical model which takes into account the heat transferred between the solid phase and the fluid phase. Consequently, in addition to Darcy-Rayleigh number Ra and the inclination angle /, two further non dimensional numbers are introduced: the inter-phase heat transfer parameter H and the porosity modified conductivity ratio c. In a recent paper (Barletta and Rees, 2015), the linear stability analysis of the basic monocellular flow indicated that the inclination angle promotes the appearance of longitudinal rolls as the preferred mode of convection. The current paper focuses on the nonlinear evolution of longitudinal rolls in a supercritical regime of convection. A weakly nonlinear analysis, using a derived amplitude equation, is adopted to determine the nonlinear effects of the parameters Ra; /; H and c. The results indicate that in inclined layers (i) the nonlinearity decelerates the mean flow; (ii) the heat transfer, determined by the evaluation of the Nusselt number (Nu) at the layer boundary, corresponds to the one obtained for horizontal layers by scaling Ra with cos /, i.e. Nu ¼ NuðRa cos /; H; cÞ; (iii) in accordance with existing laboratory experiments, the slope of Nu is less than 2, where 2 is the value predicted by the local thermal equilibrium model, and the slope represents the derivative of Nu with respect to the distance of the critical parameter from the threshold value for the onset of instability; (iv) increasing values of both H and c produce an enhancement of the heat transfer across the layer. Finally, the comparison between the LTNE theoretical predictions and existing experiments conducted with various combinations of solid matrix and fluids suggests a possible alternative way to determine the heat transfer coefficient H. Ó 2017 Published by Elsevier Ltd.

1. Introduction The problem of convective instability in a porous medium heated from below and saturated by a Newtonian fluid has been investigated extensively in the past. The work devoted to this area is well documented by the reviews of Nield and Bejan [1], Rees [2], Tyvand [3] and Barletta [4]. Among the experimental investigations aimed at visualizing the convective patterns and the temperature distributions we mention Elder [5], Combarnous [6], Close et al. [7], Shattuck et al. [8] and Howle et al. [9] for horizontal layers and Combarnous [10] when the porous layer is inclined to the horizontal.

⇑ Corresponding author. E-mail addresses: [email protected] (M.N. Ouarzazi), silvia.hirata@ univ-lille1.fr (S.C. Hirata), [email protected] (A. Barletta), michele.celli3@ unibo.it (M. Celli). http://dx.doi.org/10.1016/j.ijheatmasstransfer.2017.05.084 0017-9310/Ó 2017 Published by Elsevier Ltd.

For the horizontal configuration, the experimental results obtained with various combinations of solid particles and working fluids reveal that the heat transfer rate is not only a function of the Darcy-Rayleigh number Ra but it can also be significantly affected by the structure of the medium and the fluid properties as well. Particularly, the measured slope of the Nusselt number versus the relative distance to the critical Darcy-Rayleigh number was found to depend on the solid/fluid combination and was less than the theoretical prediction of 2 derived by Joseph [11]. The explanation for the difference between theory and experiments may lie with the effects of finite heat transfer coefficient between fluid and solid phases. Combarnous and Bories [12] proposed a local thermal non-equilibrium (LTNE) model with two-energy equations, which introduces a finite inter-phase heat transfer coefficient and the ratio between fluid and solid conductivities as additional parameters. We mention that the growing volume of work with the LTNE model is well documented in [1,12–27].

400

M.N. Ouarzazi et al. / International Journal of Heat and Mass Transfer 113 (2017) 399–410

Nomenclature a A ax ; az c d ^y ; e ^z ^x ; e e g; g h H J k km K L L; L0 LTE LTNE N Nu p P Ra Re; Im S; e S SLTE SLTNE t T0 T u

dimensionless wave number amplitude of convection components of the dimensionless wave vector specific heat diameter of the beads unit vectors in the ðx; y; zÞ-directions gravitational acceleration vector; modulus of g volumetric inter–phase heat transfer coefficient volumetric dimensionless inter–phase heat transfer parameter, Eq. (2) heat flux thermal conductivity effective thermal conductivity, vkf þ ð1  vÞks permeability layer thickness linear operators local thermal equilibrium local thermal non-equilibrium nonlinear expression Nusselt number dimensionless pressure disturbance amplitude, Eq. (10) dimensionless pressure disturbances, Eq. (8) Darcy–Rayleigh number, Eq. (2) real part; imaginary part transformed Darcy-Rayleigh numbers, Eqs. (11) and (14) slope of Nusselt number according to local thermal equilibrium slope of Nusselt number according to local thermal nonequilibrium dimensionless time, Eq. (2) temperature of the upper boundary dimensionless temperatures, Eq. (2) dimensionless velocity vector, ðu; v ; wÞ, Eq. (2)

On the other hand, in early experiments Bories and Combarnous [10] examined the secondary flow configurations of convection in a rectangular porous medium heated from below and inclined to the horizontal. The temperature recordings indicated that two main types of convective structures may be observed at the onset of convection. For small inclination angle /, the vortex patterns are polyhedral cells (i.e. oblique rolls). For higher values of /, the polyhedral cells are replaced by stationary longitudinal flow. The observed transition between the two types of convective patterns occurs at a critical angle /c ’ 15 . Some hysteresis effects associated to this transition were also observed. However, the heat transfer measured through the boundaries for a large range of slopes (0 –60 ) is found to be independent of the shape of the convective patterns, whether longitudinal rolls or polyhedral cells and a unique relation between Nusselt number and Ra cos / has been found. Caltagirone and Bories [28] examined the transition between the secondary flows in polyhedral cells and the longitudinal rolls both by three-dimensional numerical simulation of the problem and by performing a linear stability analysis of the basic state. A recent note by Nield [29] contains an interesting discussion on the results obtained by Caltagirone and Bories [28] and gives new insights into the question of the preferred patterns at the onset of the instability: rolls or polyhedral cells. Further results on the stability of an inclined porous layer were obtained by Storesletten and Tveitereid [30], Karimi-Fard et al. [31], Rees and Bassom [32], Rees et al. [33]. Karimi-Fard et al. [31] carried out

U ~ V x

dimensionless velocity disturbance vector, ðU; V; WÞ, Eq. (5) perturbations vector, ðHf ; Hs ; PÞT , Eq. (16) dimensionless position vector, ðx; y; zÞ, Eq. (2)

Greek symbols a thermal diffusivity am effective thermal diffusivity, km =ðqcÞf b thermal expansion coefficient d coefficient, Eq. (36) c dimensionless parameter, Eq. (2) DT reference temperature difference e dimensionless perturbation parameter, Eq. (5) h dimensionless temperature disturbance amplitudes, Eq. (10) H dimensionless temperature disturbances, Eq. (5) k dimensionless parameter, Eq. (2) m kinematic viscosity q density u arbitrary constant, Eq. (49) / inclination angle to the horizontal U transformed angle, Eq. (11) v porosity x complex dimensionless parameter, Eq. (10) s relaxation time, Eq. (41) Superscript, subscripts H dimensional quantity b basic solution c critical value f fluid phase s solid phase

an investigation of oscillatory instability for the case of doublediffusion. Rees and Bassom [32] defined a Squire-like transformation allowing a general study of normal modes with an arbitrary orientation. Storesletten and Tveitereid [30] included in the stability analysis the effect of anisotropy in the porous medium, while Rees et al. [33] extended this analysis by considering an arbitrary orientation of the principal axes of anisotropy. Recently Barletta and Rees [34] revisited the topic of instability in an inclined porous layer and performed a linear stability analysis in the framework of LTNE model by assuming an infinite extent of the porous cavity in the transverse and longitudinal directions. These authors showed that the longitudinal rolls are the preferred mode of instability at the onset of convection. The neutral stability for the longitudinal rolls is found to correspond to the one obtained for a horizontal layer, by scaling the Darcy-Rayleigh number with cosine of the inclination angle. At the pore level, the interface between solid and fluid cannot display any discontinuity of temperature. However, when average temperatures are evaluated over a reference elementary volume (REV) for the solid and the fluid, these temperatures may well be different. This effect can arise in an unsteady regime, but also under steady conditions in cases where the thermal conductivities of the fluid and the solid are markedly different [1]. The LTNE model of heat transfer in porous media provides a description of how the different temperature fields of the solid phase and of the fluid phase interact.

401

M.N. Ouarzazi et al. / International Journal of Heat and Mass Transfer 113 (2017) 399–410

In the present paper, we propose a weakly nonlinear stability analysis by using a LTNE model motivated by a desire to shed new light on: (i) the nonlinear effects of the Darcy-Rayleigh number, the finite inter-phase heat transfer coefficient, the porosity modified conductivity ratio and the inclination angle on the properties of the convective longitudinal rolls and on the heat transport estimated by the Nusselt number; (ii) the experimental results for both horizontal and inclined porous layers and the determination of inter-phase heat transfer coefficient; (iii) the limits of Darcy’s law in a porous medium in which the pore scale is not microscopically small. This task will be accomplished by reducing the full dynamic equations describing the system in the vicinity of the critical conditions to time-dependent equation describing the nonlinear evolution of the longitudinal rolls amplitude. Therefore this work is a natural extension of [34]. 2. Mathematical model and the basic solution Let us consider an inclined porous layer saturated by a fluid. We 



denote as / 2 ½0 ; 90  the inclination angle to the horizontal. The boundary planes, yH ¼ 0; L, are assumed to be impermeable and isothermal with different temperatures: T 0 þ DT is the temperature of the lower boundary, while T 0 is the temperature of the upper boundary, with DT > 0. A sketch of the layer is given in Fig. 1. We assume that the saturated porous medium is isotropic and homogeneous, that the effect of viscous dissipation can be neglected, and that the local thermal non-equilibrium (LTNE) can be described by a two-temperature model [1]. Thus, according to the Oberbeck-Boussinesq approximation, we can write the local mass, momentum and energy balance equations in a dimensionless form as

$  u ¼ 0; $u ¼ k

1þc

c

ð1aÞ    ^x þ cos / e ^y ; Ra $ T f sin / e

@T s ¼ r2 T s þ HcðT f  T s Þ; @t

@T f þ u  $T f ¼ r2 T f þ HðT s  T f Þ: @t

ð1bÞ

The dimensionless quantities employed in Eqs. (1) are defined as

ðx; y; zÞ ¼ ðxH ; yH ; zH Þ

ð1dÞ

t ¼ tH

u ¼ ðu; v ; wÞ ¼ ðuH ; v H ; wH Þ

L

vaf

af L2

;

¼ uH

L

;

T s;f ¼

TH s;f  T 0 ; DT

vaf 2 gbDTKL vkf af hL Ra ¼ ; c¼ ; k¼ ; H¼ : am m ð1  vÞks as vkf

ð2Þ

Here, the stars denote the dimensional time, coordinates and fields, while the subscripts ‘‘s” and ‘‘f” denote the solid phase and the fluid phase, respectively. In Eqs. (2), L is the layer thickness, a is the thermal diffusivity, k is the thermal conductivity, b is the thermal expansion coefficient of the fluid, m is the kinematic viscosity of the fluid, v is the porosity and K is the permeability of the porous medium, g is the modulus of the gravitational acceleration vector g; am ¼ km =ðqcÞf is the effective thermal diffusivity, namely the ratio between the effective thermal conductivity of the saturated porous medium, km , and the product between the density, q, and the specific heat, c, of the fluid. The inter-phase heat transfer coefficient h serves to model the heat exchange between the fluid phase and the solid phase. We point out that h describes a local volumetric heat transfer and, as a consequence, it is expressed in W=ðm3 KÞ instead of W=ðm2 KÞ as it usually happens when modelling surface heat transfer. The four dimensionless parameters, appearing in the governing Eqs. (1) and defined by Eq. (2), are the Darcy-Rayleigh number, Ra, the conductivity ratio, c, the diffusivity ratio, k, and the inter-phase heat transfer parameter, H. We note that the definition of DarcyRayleigh number given in Eq. (2) is consistent with that declared under local thermal equilibrium (LTE) conditions. The same definition has been recently employed also by Barletta and Rees [27]. An alternative definition has been proposed, for instance, in Banu and Rees [17]. These authors defined a Darcy-Rayleigh number that coincides with Ra multiplied by ð1 þ cÞ=c. Eq. (1b) expresses the local momentum balance including, on the right hand side, the buoyancy force contribution. This equation is obtained by taking the curl of Darcy’s law, so that the pressure force term does not appear explicitly. The boundary conditions on the velocity field and the temperature field are expressed in dimensionless form as

y¼0: ð1cÞ

1 ; L

y¼1:

v ¼ 0; v ¼ 0;

T s ¼ T f ¼ 1;

ð3Þ

T s ¼ T f ¼ 0:

A basic solution of the problem defined in Eqs. (1) and (3) can be written as a combination of a thermal stratification in the y direction for both the solid and the fluid phase and a non-uniform flow in the x direction with a vanishing mass flow rate,

ub ¼

  1þc 1  y Ra sin /; c 2

v b ¼ 0;

wb ¼ 0;

ð4Þ

T sb ¼ T fb ¼ 1  y; where the subscript ‘‘b” stands for ‘‘basic solution”. This basic flow is induced by the buoyancy force, and describes a single-cell circulation having an infinite width along the x-direction. 3. Governing equations for finite amplitude disturbances The perturbation of the basic velocity and temperature fields Eq. (4) can be expressed as

u ¼ ub þ U;

T s;f ¼ T sb;fb þ Hs;f :

ð5Þ

Thus, Eqs. (1) yield the nonlinear set of equations Fig. 1. The fluid saturated porous layer and the thermal boundary conditions.

$  U ¼ 0;

ð6aÞ

402

M.N. Ouarzazi et al. / International Journal of Heat and Mass Transfer 113 (2017) 399–410

$U

   ^x  cos / e ^y ¼ 0; Ra $  Hf sin / e

1þc

c

@ Hs  r2 Hs  HcðHf  Hs Þ ¼ 0; @t

k

ð6bÞ

ð6cÞ

ð6dÞ

The boundary conditions for the disturbances ðU; HÞ are

y ¼ 0; 1 :

V ¼ 0;

Hs ¼ Hf ¼ 0:

ð7Þ

A pressure representation of the velocity disturbance U can be written as follows

U¼

1þc

c

  ^x þ cos / e ^y  $; Ra Hf sin / e

ð8Þ

allows one to satisfy identically Eq. (6b), and to express the other governing Eqs. (6) and boundary conditions (7) as

  @ Hf @ Hf 1 þ c 1 1þc @P  y Ra sin /  Ra Hf cos / þ þ c 2 c @y @t @x  r2 Hf  HðHs  Hf Þ   @ Hf 1 þ c @P Ra Hf sin /  ¼ c @x @x   @ Hf 1 þ c @P @ Hf @P þ  ; Ra Hf cos /  c @y @y @z @z k

@ Hs  r2 Hs  HcðHf  Hs Þ ¼ 0; @t

r2 P 

1þc

c

y ¼ 0; 1 :

  @ Hf @ Hf Ra sin / þ cos / ¼ 0; @x @y @P ¼ 0; @y

Hs ¼ Hf ¼ 0:



Ra cos / ¼ S cos U;

ð11Þ

ax Ra sin / ¼ aS sin U:

  @ Hf @ Hf 1 þ c 1  y Ra sin /  V  r2 Hf  HðHs  Hf Þ þ c 2 @t @x ¼ U  $Hf :

replaced with a transformed angle U, and the Darcy-Rayleigh number Ra has been replaced with a transformed Darcy-Rayleigh number S,

ð9aÞ

The minimum of function SðaÞ defines the critical values ac and Sc . These critical values were computed in Barletta and Rees [34] for different values of transformed angle U. The most unstable normal modes are found to be longitudinal rolls defined by

U ¼ 0 ;

S ¼ Ra cos /:

ð12Þ 

This instability exists for every tilt angle / smaller than 90 . At neutral stability, the dependance of S on H and c is given by

S¼

2   a2 þ p2 a2 þ p2 þ H þ Hc : a2 ða2 þ p2 þ HcÞ 1þc

c



ð13Þ

We note that, on account of Eq. (12), the marginal stability condition expressed in terms of the Darcy-Rayleigh number Ra exploits a dependence on the inclination angle / through a factor 1= cos /. As a consequence, an increasing inclination to the horizontal implies a stabilizing effect on the basic flow, viz. a monotonic increment of the lowest eigenvalue Ra corresponding to assigned values of ða; H; cÞ, evaluated through Eqs. (12) and (13). We mention that the special case of c ! 0 is of great importance for nonlinear stability analysis. In this limit, Barletta and Rees [34] introduced a parameter e S defined as

ð9bÞ

1þc e S: S¼

ð9cÞ

The minimum of e SðaÞ defines the critical values for longitudinal rolls and are given in the limiting case of c ! 0, with H  Oð1Þ by

ð9dÞ

ð14Þ

c

ac ¼

pffiffiffiffi

p H þ p2

1=4

;

e S c ¼ H þ 2p



pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi

p þ H þ p2 :

ð15Þ

For c – 0, the variation of e S c with inter-phase heat transfer coefficient H for different values of the porosity modified conductivity ratio c is shown in Fig. 2. This figure indicates that e S c increases from

4. Review of linear stability analysis The linear stability analysis is the starting point for treating the nonlinear problem. Since the latter is the main objective of the present work, we shall briefly sketch the linear analysis, by emphasizing the most important points. The linear problem is governed by Eqs. (9) without the right hand side in Eqs. (9a). As usual, we study the normal modes expressed as

e S c ¼ 4p2 for H ¼ 0 independently of c to the asymptotic LTE condition e S c ¼ 1þc 4p2 or equivalently Sc ¼ 4p2 . The LTE condition is c

ð10Þ

attained with large values of H for c  Oð1Þ. The variation of the critical wave number ac with H for different values c is depicted in Fig. 3. We observe from this figure that for small values of c, the critical wave number increases from ac ¼ p to reach its maximum at a well defined value of H before decreasing to the LTE condition ac ¼ p when H ! 1. This figure also indicates that the

pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ^ x þ az e ^z is the wave vector, a ¼ a2x þ a2z is the wave numwhere ax e ber, and x is a complex parameter. The real part of x; ReðxÞ, is the angular frequency of the wave. The imaginary part, ImðxÞ, is the exponential growth parameter. If ImðxÞ < 0, the normal mode describes a stability condition. On the other hand, if ImðxÞ > 0, the normal mode yields instability. The breakup of the single-cell flow regime occurs at the marginal stability threshold, namely when ImðxÞ ¼ 0. The marginal stability condition is defined by ImðxÞ ¼ 0. General modes with 0 6 ax 6 a and 0 6 az 6 a are termed oblique rolls. The modes with ax ¼ a and az ¼ 0 are called transverse rolls, while the modes with ax ¼ 0 and az ¼ a are the longitudinal rolls. In Barletta and Rees [34], the authors used a Squire-like transformation, to map all the oblique rolls onto corresponding transverse rolls. Accordingly the tilt angle of the layer, /, has been

Fig. 2. Critical value of e S as a function of H. The arrow indicates increasing values of c, namely c ¼ 0:05; 0:2; 0:4; 0:6; 1; 10.



Pðx; y; z; tÞ

Hs;f ðx; y; z; tÞ



 ¼

pðyÞ hs;f ðyÞ

exp ½iðax x þ az z  xtÞ;

403

M.N. Ouarzazi et al. / International Journal of Heat and Mass Transfer 113 (2017) 399–410

Above the linear threshold, we introduce a small parameter e which measures the distance to criticality by setting e S 2 is of order unity. This fixes the temporal S 2 , where e S¼e S c þ e2 e scaling to

t2 ¼ e2 t:

ð20Þ

Temporal derivative is then replaced by

@ @ @ ! þ e2 : @t @t @t 2

Fig. 3. Critical value of a as a function of H. The arrow indicates increasing values of c, namely c ¼ 0:05; 0:2; 0:4; 0:6; 1; 10.

critical wave number ac is not significantly influenced by H when c  Oð1Þ. 5. Nonlinear stability analysis Once the onset of convective instability has been identified, a nonlinear stability approach is needed to determine the evolution of the unstable longitudinal rolls and the associated heat transfer as functions of the non dimensional parameters ðRa; /; c; HÞ. In this part, the control parameter used is e S ¼ Ra cos /ð1 þ cÞ=c instead of

Ra. We restrict our analysis to the case where the parameter e S is close to the threshold e S c of the onset of longitudinal rolls. In this case, the weakly nonlinear dynamics of the linearly unstable longitudinal rolls can be described by an amplitude equation. 5.1. Derivation of the amplitude equation

e ¼ N; ðL0 @ t þ LÞ V

ð16Þ

e ¼ ðHf ; Hs ; PÞT contains the temperature perturwhere the vector V bation of the fluid, the temperature perturbation of the solid and the pressure perturbation. The expressions of the linear operators L0 and L and the nonlinear expression N are the following:

1 0 0

B L0 ¼ @ 0

1

C 0 A;

k

B B L¼B @

@2  @y 2



@2 @z2

þHe S

H

@ @y

2

0

2

Hc

@ @  @y 2  @z2 þ H c

@ e S @y

0

  T @P @ Hf @P @ Hf e S Hf þ þ ; 0; 0 : @y @y @z @z

ð22Þ

e i depend also on the slow variable t 2 . By substituting The vectors V the expressions (22) and (21) in the system (16), and then collecting coefficients of OðeÞ, a set of equations is obtained. 5.1.1. Solution of the first order problem To first order, the set of linear homogeneous equations for e 1 ¼ ðHf1 ; Hs1 ; P1 ÞT at criticality is given by V

e 1 ¼ 0: Lc V

ð23Þ

with Lc is the linear operator L for e S¼e S c . The boundary conditions e for V 1 are

@P1 ¼ 0; @y

y ¼ 0; 1 :

Hs1 ¼ Hf1 ¼ 0:

ð24Þ

Hf1 ¼ Aðt2 Þeiac z sinðpyÞ þ cc; Hs1 ¼

@2 @y2

2

@ þ @z 2

C C C; A

ð25Þ

Hc Aðt 2 Þeiac z sinðpyÞ þ cc; Hc þ p2 þ a2c

P1 ¼ 

ð26Þ

p e S Aðt Þeiac z cosðpyÞ þ cc; p2 þ a2c c 2

ð27Þ

where cc stands for a complex conjugate. The quantities e S c and ac are already obtained by linear stability analysis [34]. 5.1.2. Adjoint problem Proceeding with the higher orders of the expansion requires the determination of the adjoint mode. In the present case, the appropriate scalar product is

D

E e i; V ej ¼ V

1 2p=ac

Z 2p=ac Z

1 0

0

ei  V e  dz dy; V j

ð28Þ

e  is the complex conjugate of V e j. where V j With this definition, the adjoint operator of L is

0 B B Ly ¼ B @

1

@ @ e  @y 2  @z2 þ H  S 2

Hc

H

@ @  @y 2  @z2 þ H c

2

2

@  @y

ð18Þ

@ e S @y

2

0

0 @2 @y2

þ

@2 @z2

1 C C C: A

ð29Þ

e y of the adjoint At criticality, one easily finds the eigenfunctions V problem,

and

N¼

e ¼ eV e 1 þ e2 V e2 þ    V

ð17Þ

0 0 0 0

~ in The evolution equations are obtained by expanding the vector V power series of e:

The solutions of Eq. (23) satisfying the boundary conditions (24) are

In order to derive the amplitude equation, we use the perturbation expansion near the critical threshold. In solving weakly nonlinear hydrodynamic problems with this method for a non-selfadjoint linear operator, we must also use the corresponding adjoint operator. Then we perform inner product using both, the solutions of the homogeneous linear operator and of the adjoint homogeneous linear operator. The corresponding amplitude equation is then obtained from a solvability condition, known as the Fredholm alternative, which states that the inhomogeneity must be orthogonal to the solution of the adjoint problem. Since the corresponding linear operator in Eqs. (9) is non-self-adjoint we used the above generalization of a perturbation expansion. The governing equations may be write in the compact notation

0

ð21Þ

ð19Þ

ey ¼ V

1

Hf1 ; Hs1 ;

c

1 P1 e Sc

!T :

ð30Þ

404

M.N. Ouarzazi et al. / International Journal of Heat and Mass Transfer 113 (2017) 399–410

5.2. Average heat transfer

5.1.3. Solution of the second order problem e 2 ¼ ðHf2 ; Hs2 ; P 2 ÞT At Oðe2 Þ, the non homogenous problem for V is given at criticality by,

e 2 ¼ N2 ; Lc V

ð31Þ

1

c

with boundary conditions

@P2 ¼ 0; @y

The stationary form of Eqs. (1c) and (1d) governing the temperature fields is

r2 T s ¼ HðT s  T f Þ;

ð42aÞ

r2 T f  u  $T f ¼ HðT f  T s Þ:31

Hs2 ¼ Hf2 ¼ 0 at y ¼ 0; 1:

ð32Þ

The nonlinear term N 2 is the nonlinear expression N evaluated to Oðe2 Þ. e yi ¼ 0 Calculations showed that the solvability condition hN 2 ; V

ð42bÞ

By eliminating HðT f  T s Þ from Eqs. (42a) and (42b) we obtain

1 r2 ðT f þ T s Þ þ u  $T f ¼ 0:

ð43Þ

c

is trivially satisfied. The equation at this order has then a solution

The conductive state supports a uniform heat flux Jcond ¼ Jf;cond þ Js;cond flowing vertically through the layer with the fluid and the solid contributions,

Hf2 ¼ djAðt 2 Þj2 sinð2pyÞ;

ð33Þ

J f;cond ¼ 

dT f;b ¼ 1; dy

ð44Þ

ð34Þ

J s;cond ¼ 

1 dT s;b 1 ¼ ; c dy c

ð45Þ

Hs2 ¼ 

P2 ¼

Hc djAðt 2 Þj2 sinð2pyÞ; 4p 2 þ H c

1 e S c djAðt2 Þj2 cosð2pyÞ; 2p

ð35Þ

with

d¼

2p



ðp2 þ a2c þ H þ HcÞð4p2 þ HcÞ   : p2 þ a2c p2 þ a2c þ Hc ð4p2 þ Hc þ HÞ

ð36Þ

The above expressions means either that the heat flux through the solid phase is negligible if c ! 1 or that the heat flux through the solid phase is the dominant mechanism of heat transfer through the layer if c ! 0. In the presence of a convective velocity field, the total stationary vertical flux J is independent of y and it reads

JðzÞ ¼ J cond þ wðT f;b þ Hf Þ  5.1.4. Third order solvability condition At Oðe3 Þ, the following non homogenous e V 3 ¼ ðHf3 ; Hs3 ; P3 ÞT is obtained

e 3 ¼ L0 @ t V e 1 þ N3 ; Lc V 2

problem

for

ð37Þ

with boundary conditions

@P3 ¼ 0; @y

Hs3 ¼ Hf3 ¼ 0 at y ¼ 0; 1:

ð38Þ

The nonlinear term N 3 is the nonlinear expression N evaluated to Oðe3 Þ. Similarly, applying the solvability condition e y i ¼ 0 at Oðe3 Þ yields, the behavior of A at the e 1 þ N3 ; V hL0 @ t V 2

 d Hf þ 1c Hs dy

:

ð46Þ

For longitudinal rolls with axes along the x direction, we define the lateral average heat flux by

hJðzÞi ¼

1 2p=ac

Z 2p=ac

ð47Þ

JðzÞdz: 0

The relative contributions of conduction and convection on heat transport is expressed using dimensionless Nusselt number Nu defined by

+ *  d Hf þ 1c Hs hJðzÞi hJðz; y ¼ 0Þi c jy¼0 : Nu ¼ ¼ ¼1 dy J cond J cond 1þc

ð48Þ

6. Results of the weakly nonlinear approach and discussion

time scale t2 6.1. Finite perturbation fields

S @A e s ¼ e2 A  pdAjAj2 : @t 2 S c

ð39Þ

One can get rid of the parameter e by re-introducing the original variables t ¼ t 2 =e2 , e S 2 ¼ ðe Se S c Þ=e2 , and defining A ¼ eA. This yields an equation for A in terms of the original time variable of the landau type:

@A e Se S s ¼ e c A  pdAjAj2 : @t Sc

ð40Þ

In addition to d, the coefficient s appearing in (40) depends on k; c and H. Its explicit expression is:

s¼

p2 þ a2c a2c e Sc

" 1 þ k

2

H c Hc þ p2 þ a2c

# 2 :

The amplitude Eq. (40) possesses a stationary non trivial nonlinear solution,

e Se Sc A¼ e pd S c

!1=2 e iu ;

for any arbitrary constant phase u. Up to the first nonlinear correction and by re-defining ðHf;s ; PÞ ¼ eðHf;s1 ; P1 Þ þ e2 ðHf;s2 ; P2 Þ, the perturbation temperature field for the fluid phase and for the solid phase are respectively

Hf ¼

e Se Sc pdeS c

!1=2



e  Se Sc eiu eiac z þ cc sinðpyÞ  sinð2pyÞ; peS c

ð41Þ

According to the expression (36), the cubic Landau constant pd is positive, indicating that the bifurcation from the conductive state to the convection state is supercritical, independently of the parameters c; H and /.

ð49Þ

Hs ¼

Hc

Hc þ p þ

þ

2

a2c

e Se Sc pdeS c

!1=2



 eiu eiac z þ cc sinðpyÞ

e Hc Se Sc sinð2pyÞ; 4p þ H c p e Sc 2

ð50Þ

ð51Þ

M.N. Ouarzazi et al. / International Journal of Heat and Mass Transfer 113 (2017) 399–410

P¼

1 p2 þ a2c

"

#1=2

peS c ðeS  eS c Þ

ðeiu eiac z þ ccÞ cosðpyÞ þ

d

e Se Sc cosð2pyÞ: 2p2 ð52Þ

It is interesting to note that in the limit H ! 1 with c  Oð1Þ or in the limit c ! 1 with H  Oð1Þ, Eq. (51) yields Hs ¼ Hf . These two limiting cases correspond to the local thermal equilibrium (LTE) between the two phases. We also note that, assumed d  Oð1Þ and e S c  Oð1Þ, in the limit

c ! 0 with H  Oð1Þ or in the limit H ! 0 with c  Oð1Þ, Eq. (51) yields the temperature field Hs ! 0. Substitution of Eqs. (50) and (52) into Eq. (8) leads to the expressions of the components of the perturbation velocity field,

"

e Se Sc Þ S c ðe U ¼ tan / pd  tan /

V¼

e Se Sc

p

#1=2 iu iac z

ðe e

þ ccÞ sinðpyÞ ð53Þ

sinð2pyÞ;

" #1=2 e S c ðe a2c Se ScÞ ðeiu eiac z þ ccÞ sinðpyÞ; p2 þ a2c pd

ac W ¼ i 2 p þ a2c

"

peS c ðeS  eS c Þ

ð54Þ

#1=2

d

ðeiu eiac z Þ cosðpyÞ þ cc:

ð55Þ

405

We also note from Eq. (58) that when c ! 0 with H  Oð1Þ, i.e. when conduction is the leading heat transfer mechanism, the slope SLTNE tends to 0 and therefore Nu ! 1. We mention that this limit may have a special interest when studying the behavior of porous media such as the metallic foams. In fact, the limit c ! 0 describes an asymptotic condition where kf  ks . For c – 0, the Nusselt number may be written as a function of the parameter S ¼ Ra cos / instead of e S,

Nu ¼ 1 þ SLTNE

  Ra cos /  Sc : Sc

ð59Þ

According to Eq. (58) the slope SLTNE ! 2 as H ! 1 with

c  Oð1Þ or c ! 1 with H  Oð1Þ. This result is exactly the same as that drawn for the case of LTE and derived, for instance, in [11]. The slope SLTNE of Nusselt number relationship in Eq. (58) is plotted both in Fig. 4 as a function of cfor different values of H and in Fig. 5 as a function of H for different values of c. As it can be noticed from these figures, the effect of both H and c is to increase the slope of Nusselt number. Therefore the effect of each of these parameters is to enhance the heat transfer, which reaches its maximum value with SLTNE ¼ 2 in the LTE conditions. The physical reason for this is that when H ! 1 with c  Oð1Þ or c ! 1 with H  Oð1Þ), the fluid and solid phase have almost equal temperatures and therefore may be treated as a single phase. 6.3. Relation to experiments and determination of H

The expression (53) states that a non-vanishing inclination angle / induces a component of the perturbation velocity field in the x direction, transforming a two dimensional problem when / ¼ 0 to a three dimensional one. Moreover, the inclination angle and the nonlinearity of the problem act together to produce a mean flow correction  tan / sinð2pyÞðe Se S c Þp1 . With the nonlin-

The main objective of this section is the determination of the dimensionless inter–phase heat transfer parameter H by using the theoretical results obtained in the current work in connection with experiments. In the following, results from the representative experiments [6,7,10] will be considered in detail.

ear contribution, the total mean flow including the basic flow given by Eq. (4) may be written as

6.3.1. Horizontal porous layer with very small pore size (d=L  1) Experimental informations in the case of natural convection in horizontal porous layer are available in Combarnous [6] where they are presented in the form of tables. We collected from [6] the data corresponding to the measured Nusselt number for different values of Darcy-Rayleigh number. In experiments [6], an extensive pattern of thermo-couples was set in the porous medium with a view to observe the cellular organization of convection and to measure the heat transfer. The fluid saturated porous media employed are different, each one crafted by using different solid beads (glass, quartz) characterized by various diameters d. The saturating fluids used are different (water, oil) depending on the given experimental series. The height of the experimental cell is L ¼ 53:5 mm for all the series so that (d=L  1). Selected results

umean

"  #  e e 1 S  Sc e y  ¼ S sinð2pyÞ tan /: 2 p

ð56Þ

The effect of the nonlinear modification of the mean flow is that the fluid can be decelerated. With increasing the inclination angle and/or the distance to the critical threshold e S c , the deceleration becomes stronger, specially along the surfaces y ¼ 1=4 and y ¼ 3=4. 6.2. Effects of H and c on heat transfer Substituting the solutions of the temperature fields given by Eqs. (50) and (51) into Eq. (48) yields the following relationship

Nu ¼ 1 þ SLTNE

e Se Sc ; e Sc

ð57Þ

where SLTNE is the slope of the curve Nu as a function of the relative distance to criticality, as determined by local thermal nonequilibrium model

SLTNE ¼

  2c H : 1þ 2 4p þ H c 1þc

ð58Þ

Since SLTNE and e S c are independent on the inclination angle /, Eq. (57) states that the dependence of the Nusselt number on / is only through the parameter e S ¼ c1 ð1 þ cÞRa cos /. Consequently, the average heat transfer for the longitudinal rolls is the same as that for the Darcy-Bénard problem in a horizontal layer with LTNE conditions, provided that Ra is replaced by Ra cos /.

Fig. 4. The Slope SLTNE of Nusselt number as a function of c with different values of H. The arrow indicates increasing values of H, namely H ¼ 0:1; 1; 10; 100; 1000.

406

M.N. Ouarzazi et al. / International Journal of Heat and Mass Transfer 113 (2017) 399–410

Fig. 5. The Slope SLTNE of Nusselt number as a function of H with different values of c. The arrow indicates increasing values of c, namely c ¼ 0:05; 0:1; 0:2; 0:4; 0:6; 1.

from experiments [6] will be presented in the following to illustrate the dependence of Nusselt number with respect to DarcyRayleigh number up to Ra ¼ 140 (i.e. more that three times the critical Rayleigh number). Figs. 6–8 present the heat transfer results in the ðRa; NuÞ plane for different series. The continuous lines represent a quadratic fit of the experimental data. The slopes of measured Nu as a function of the relative distance to criticality ðRa  Rac Þ=Rac are reported in Table 1 for the selected eight series, together with the ratio of the diameter of the beads to the height of the medium d=L and the corresponding porosity modified conductivity ratio c as it has been evaluated from [6]. As can be seen from Table 1, the discrepancy between the experimental slopes and the theoretical ones stemming from the limit of LTE model SLTE ¼ 2 is very important. We emphasize that this discrepancy is probably due to the non validity of the assumption of an infinite heat transfer coefficient between the fluid and the solid phases. Now, we address the question of the determination of the dimensionless inter–phase heat transfer parameter H. In the case of horizontal porous layer, the Nusselt number predicted by the weakly nonlinear analysis is given by Eq. (59) with / ¼ 0,

Nu ¼ 1 þ SLTNE

  Ra  Rac : Rac

ð60Þ

The set of Eqs. (58)–(60) is sufficient for comparison with experiments and the measured Nusselt number as a function of the distance to criticality allows the determination of the dimen-

Fig. 6. The measured Nusselt number as a function of Ra for glass beads/water combinations (from Combarnous [6]). The symbols are experimental data for series 11 (diamonds), series 6 (squares) and series 2 (circles) respectively. The lines are the interpolation of the experimental results.

Fig. 7. The measured Nusselt number as a function of Ra for quartz/water combinations (from [6]). The symbols are experimental data for series 14 (squares), series 5 (diamonds) and series 13 (circles) respectively. The lines are the interpolation of the experimental results.

Fig. 8. The measured Nusselt number as a function of Ra for glass/oil combinations (from [6]). The symbols are experimental data for series 8 (squares) and series 7 (circles) respectively. The lines are the interpolation of the experimental results.

sionless inter–phase heat transfer parameter H. Indeed, Eq. (58) yields the following relationship

H¼

  4p2 SLTNE 1 :  1þc ð1  SLTNE =2Þ 2c

ð61Þ

In order to use the experimental slope for the calculation of H, one has to know the physical parameters concerning the system, namely the porosity of the porous medium, the thermal conductivity of the fluid and the thermal conductivity of the solid phase. The accuracy on the final value of H is mainly determined by the accuracy on these thermophysical properties which allow the evaluation of the parameter c. They can either be found in literature or measured experimentally. For the problem under consideration, these physical properties are documented in [6]. For the selected series presented in this paper, the determined values of H are reported in Table 1. We have also used the average value of c for series with the same solid/fluid combination. Precisely, we used c ¼ 0:244 for glass/water (series 2, 6 and 11), c ¼ 0:051 for quartz/water (series 5; 13 and 14) and c ¼ 0:055 for glass/oil (series 7 and 8). With these prescribed values of c, the Nusselt number is represented as a function of the relative distance to the critical Rayleigh number for glass/water, quartz/water and glass/oil combinations in Figs. 9–11 respectively. In these figures the symbols are experimental data [6] and the continuous lines are the theoretical results as predicted by the set of Eqs. (58)–(60) represented for

407

M.N. Ouarzazi et al. / International Journal of Heat and Mass Transfer 113 (2017) 399–410

Table 1 Determination of H by combining theoretical and experimental results for different combinations of solid matrix and fluids and different d=L ([6,7] for horizontal cavity and [10] for inclined cavity) References

Solid/liquid

d L

c

Slope of Nu (experimental)

H

hd vkf

[6] (series2) [6] (series6) [6] (series11) [6] (series5) [6] (series13) [6] (series14) [6] (series7) [6] (series8) [7] (1 layer) [7] (2 layers) [7] (3 layers) [7] (4 layers) [7] (5 layers) [10]

Glass/water Glass/water Glass/water Quartz/water Quartz/water Quartz/water Glass/oil Glass/oil Polyprop./air Polyprop./air Polyprop./air Polyprop./air Polyprop./air Glass/water

0.033 0.075 0.056 0.037–0.047 0.034–0.037 0.037–0.047 0.037 0.003–0.019 0.95 0.56 0.41 0.32 0.26 0.088

0.2504 0.2359 0.2462 0.0504 0.0538 0.0479 0.0541 0.0560 0.22 0.15 0.127 0.112 0.110 0.235

1.04 0.93 0.9 0.64 0.50 0.55 1.20 0.92 0.4 0.52 0.80 1.17 1.18 1

105.90 86.18 73.68 316.81 194.96 263.42 996.81 527.38 4.41 46.99 148.85 408.27 425.88 98.33

0.11 0.48 0.23 0.56 0.24 0.46 1.36 0.15 3.98 14.74 24.41 40.77 28.35 0.076

different values of H. A close inspection of these figures allows one to draw some remarks. As can be seen from Fig. 9, the value of H which captures the experimental results seems to depend on the relative distance to the critical Rayleigh number. Two regions are approximately identified. When the distance to criticality is less than 0:5, the value H ¼ 10 describes very well the behavior of Nusselt number. By increasing the distance to criticality beyond 0:5, the value of H ¼ 10 becomes inadequate and is replaced by H ’ 87:86. The same observation may be made for results presented in Fig. 10 when a jump in the appropriate value of H from H ’ 100 to H ’ 255 is observed. We can conclude from these observations that the influence of the filtration velocity induced by convection on the heat transfer coefficient H must be taken into account. Nevertheless, the influence of convective movements on H becomes less significant for a well developed convection (i.e., sufficiently far from criticality) as can be seen from Fig. 9 with H ’ 87:86, from Fig. 10 with H ’ 255 and from Fig. 11 with H ’ 725. In that case, the determination of heat transfer coefficient H requires solely the knowledge of the porosity modified conductivity ratio c.

2

Fig. 10. The Nusselt number as a function of the relative distance to the critical Rayleigh number for quartz/water combinations. The symbols are experimental data [6] for series 14 (squares), 5 (diamonds) and 13 (circles) while continuous lines are the theoretical predictions using different values of H.

6.3.2. Horizontal porous layer with d=L ¼ Oð1Þ Experiments on convection in horizontal porous layers conducted by Close et al. [7] revealed that the ratio of the particle diameter to the height of the porous medium d=L had a strong influence on Nusselt number. Between one and five layers of polypropylene (polyprop. in Table 1) spheres with diameter d ¼ 18:4 mm were used in these experiments with air as a working fluid. Depending on the number of layers, the height of the porous

Fig. 11. The Nusselt number as a function of the relative distance to the critical Rayleigh number for glass/oil combinations. The symbols are experimental data [6] for series 8 (squares) and 7 (circles), while continuous lines are the theoretical predictions using different values of H.

Fig. 9. The Nusselt number as a function of the relative distance to the critical Rayleigh number for glass/water combinations. The symbols are experimental data [6] for series 11 (diamonds), 6 (squares) and 2 (circles), while continuous lines are the theoretical predictions using different values of H.

cavity has been adjusted to L ¼ 19:3 mm, L ¼ 32:4 mm; L ¼ 45:4 mm, L ¼ 58:2 mm and L ¼ 71:3 mm respectively. Consequently the ratio of particle diameter to the height of the porous medium varies from d=L ¼ 0:95 for one layer configuration to d=L ¼ 0:26 for five layers packing, as it is indicated in Table 1. Plots of the variations of Nusselt number as a function of Darcy-Rayleigh number are illustrated in Fig. 5 of [7] for one to five layers. The slops of the curves of Nusselt number as a function of

408

M.N. Ouarzazi et al. / International Journal of Heat and Mass Transfer 113 (2017) 399–410

the relative distance to criticality is reported in Table 1. These data show a consistent trend of increasing the slope of Nusselt number with number of layers or equivalently with the ratio of the particle diameter to the height of the porous medium d=L. It is interesting to note that these slopes lie between 0:4 and 1:18, which does not agree with the predicted value of 2 by local thermal equilibrium model. From [7], we have also evaluated the corresponding porosity modified conductivity ratio c for layers between one and five. As c and the slope of the measured Nusselt number are determined, Eq. (61) yields the values of the dimensionless inter–phase heat transfer parameter H. The values of H are reported in Table 1 for different number of layers. The above results provide a test of the usual assumption of homogeneous porous media which is necessary to use Darcy’s law. In order to understand the influence of the ratio of the particle diameter to the height of the porous medium d=L on convection, let us write the dimensionless inter-phase heat transfer parameter H as the product of two dimensionless parameters, 2

H¼

2

hL hd L2 ¼ : vkf vkf d2

ð62Þ

The effect of local heat transfer is described by the dimensionless 2

parameter hd =vkf , while the influence of the pore size on heat 2

transfer depends on L2 =d . For homogeneous porous media, the 2

dimensionless parameter hd =vkf is affected solely by the thermal characteristics of the solid and the fluid phases as well as the texture of the medium. We can therefore identify a more or less

sion in both x and z directions. We believe, as Nield [29] suggested, that the presence of vertical end walls in real experiments may play a key role in the pattern selection mechanism. We expect that the lateral and transversal confinement of the porous cavity may stabilize the longitudinal rolls and consequently may promote the appearance of oblique rolls for small inclination angle. From the phenomenological point of view, there is a close analogy between the problem under consideration and the mixed convection in a porous medium heated from below and subjected to an imposed horizontal through flow [36,37]. In these papers, the authors showed that the lateral confinement induces a transition of the convective structures from oblique rolls to longitudinal rolls. However, in [10] the heat transfer measured through the boundaries is found to be independent of the shape of the convective cells, whether longitudinal or oblique rolls. Bories and Combarnous [10] presented in Fig. 7 of their paper the Nusselt number Nu as a function of Ra cos / for different values of the inclination angle /, namely 7:5 15 22:5 30 45 60 . The slope of the curve Nu as a function of Ra cos / is found to be independent of /, in agreement with our theoretical prediction. We also note that the experimental slope of Nu as a function of the relative distance to criticality is approximately equal to 1, and it is close to the slopes found in the case of horizontal cavities for glass/water (series 2, 6 and 11) as shown in Table 1. For the determination of H, we have used the average value of c ¼ 0:235 and the average value of d=L ¼ 0:088. In these conditions, we found H ¼ 98:3 which is comparable to its values for series 2, 6 and 11 in the case of horizontal porous cavity.

2

heterogeneous porous medium by the influence of L=d on hd =vkf . 2

When the variation of hd =vkf as a function of L=d is significant, the porous medium is more heterogeneous and the usual theoretical hypothesis of spatial uniformity may not be justified. In Table 1, 2

we report the values of the dimensionless parameter hd =vkf for experiments [6,7] for different combinations of solid matrix and fluids and different d=L. A close inspection of the last column in Table 1 2

reveals that the variation of hd =vkf is weak for experiments [6] where d=L  1, while its variation becomes significant for d=L  Oð1Þ in experiments [7]. These observations state that the structure of the porous medium may play a non negligible role when the pore scale is not too smaller with respect to the height of the porous layer. This reasoning allows one to conclude that the modification of the existing theory based on Darcy’s law is an important issue. One approach that can be used to solve the issue is assuming that the medium properties, such as the permeability, the thermal conductivity or the Darcy-Rayleigh number, are functions of the coordinates. In particular, theoretical and numerical investigation performed by Ouarzazi et al. [35] showed that slowly spatial variations in Rayleigh number may lead to the localization of the convective patterns, as it was observed experimentally by Shattuck et al. [8] and Howle et al. [9]. 6.3.3. Inclined porous layer In early experiments Bories and Combarnous [10] examined the secondary flow configurations of convection in a rectangular porous medium heated from below and inclined to the horizontal. The porous media were made of spherical glass beads with diameters of 5:25; 4:7 and 3:25 mm. De-aerated water has been used as the working fluid. The temperature measurements indicated that two main types of convective structures may be observed at the onset of convection. For small inclination angle /, the vortex patterns are polyhedral cells. For higher values of /, the polyhedral cells are replaced by stationary longitudinal flow. The observed transition between the two types of convective patterns occurs at some defined critical angle. This transition is not predicted in the current work where the porous layer is assumed of infinite exten-

6.3.4. Range of validity of the weakly nonlinear analysis In the present paper, we compare the results of convective heat transfer given by the amplitude equation with the data corresponding to measured Nusselt number for different DarcyRayleigh numbers. The range of validity of the weakly nonlinear analysis based on the cubic amplitude equation has to be discussed. The rigorous derivation of such equation requires the theoretical assumption that Darcy-Rayleigh number is near its critical value at the onset of convection (i.e.   1). Moreover, the range in which higher orders in  resulting from higher modes are negligible is not obvious. In practical applications of the weakly nonlinear analysis, the range of validity of amplitude equations depends on the physical system under consideration. Generalis and Fujimura [38] examined the range of validity of amplitude equations in the context of Rayleigh-Bénard problem in pure fluids without porous matrix. They derived amplitude equation to the higher order (from the 3th order up to the 25th order) and evaluate their steady solutions for various Prandtl numbers. They also carried out fully numerical analysis on bifurcation of steady solutions. Comparing two results, they showed that the amplitude equation has a wide range of validity for Prandtl number larger than Oð1Þ. For fluids with small Prandtl numbers, on the other hand, they showed unexpectedly a very narrow range of validity. In Figs. 7–9, we present results from experiments for DarcyRayleigh number up to Ra ’ 140 or similarly to 2 ¼ 2. We remark from these figures a linear dependence of Nusselt number on Darcy-Rayleigh number. The data from experiments for Ra > 140, not reported in these figures showed a deviation from this linear dependence. The linear variation of measured Nusselt number as a function of Darcy-Rayleigh number encouraged us to push the comparison up to Ra ’ 140, although we are aware of the limitation of the prediction of the cubic amplitude equation up to this extreme value of Rayleigh number. In order to assess the validity of the weakly nonlinear analysis up to this range of DarcyRayleigh number, numerical simulations of the full problem are

M.N. Ouarzazi et al. / International Journal of Heat and Mass Transfer 113 (2017) 399–410

needed. This task is out of the scope of the present investigation and is postponed to a future work.

7. Conclusions The present work is focussed on the weakly nonlinear stability analysis of a basic stationary single-cell flow in an inclined porous layer saturated by a fluid. The mathematical formulation is based on local thermal non-equilibrium (LTNE) model and it introduces two dimensionless parameters: the finite inter-phase heat transfer coefficient H and the porosity modified conductivity ratio c. A perturbation analysis of the basic stationary single-cell flow showed in linear theory that convection starts to develop in the form of stationary longitudinal rolls. A weakly nonlinear theory has been applied in order to derive an amplitude equation for the marginally unstable longitudinal rolls. In the framework of the derived amplitude equation, the finite perturbation fields are obtained and a simple analytical expression of Nusselt number Nu is determined as a function of Darcy-Rayleigh number Ra, the inclination angle /; H and c. The main results that have been obtained can be summarized as follows. (a) The bifurcation from the basic stationary single-cell flow to longitudinal rolls is a supercritical pitchfork bifurcation independent of the parameters /; H and c. (b) The nonlinearity introduces a modification of the main flow that tends to decelerate the fluid. The deceleration is more significant for increasing values of both the inclination angle and the departure of the critical parameter from the threshold value for the onset of instability. (c) Nu depends on / through a cosine factor, thus allowing a simple correspondence to the Nusselt number relative to the horizontal layer. It is worth nothing that this result is in line with existing laboratory experiments. (d) Nu increases with increasing values of H and c meaning that these dimensionless parameters promote the heat transfer. We mention that in the limit c ! 0, where the solid phase is much more conductive than the fluid phase, Nu ! 1 meaning that the heat transferred by convection is negligible compared to the heat transferred by conduction. On the other hand, the heat transfer associated to local thermal equilibrium (LTE) is recovered in both limits H ! 1 with c ¼ Oð1Þ or c ! 1 with H ¼ Oð1Þ. (e) A quantitative comparisons of the Nu predicted theoretically in this study with existing experiments is performed. This comparison demonstrates that the LTE model is not appropriate to describe the heat transfer due to convection in both horizontal or inclined porous cavity. We have shown that the comparison between the measured and the theoretically predicted heat transfer in the nonlinear regime is a possible way to deduce the unknown inter-phase heat transfer coefficient H in LTNE model. In addition, the comparison of the theoretically predicted results with existing experiments performed for porous media characterized by a pore scale that is not microscopically small, showed that the assumption of homogeneous system fails. Consequently, Darcy’s law may be modified to take into account the heterogeneous character of this kind of porous media. Work in this direction is in progress. The current weakly nonlinear stability analysis was carried out by assuming that the porous cavity is infinite in the horizontal plane. However, since experiments can only take place in finite containers, the comparison between theoretical and experimental results is not always an easy task as sidewalls play an important

409

role. An interesting extension of this work may deal with the analysis of the influence of the horizontal confinement on the linear and the nonlinear pattern selection by using the LTNE model. This task is specially challenging for understanding the transition observed in the experiments [10] between oblique rolls and transversal rolls at some critical inclination angle.

Conflict of interest All authors declare no conffiicts of interest. References [1] D.A. Nield, A. Bejan, Convection in Porous Media, fifth ed., Springer, New York, 2017. [2] D.A.S. Rees, The stability of Darcy-Bénard convection, in: K. Vafai, H.A. Hadim (Eds.), Handbook of Porous Media, CRC Press, New York, 2000, pp. 521–558 (Chapter 12). [3] P.A. Tyvand, Onset of Rayleigh-Bénard convection in porous bodies, in: D.B. Ingham, I. Pop (Eds.), Transport Phenomena in Porous Media II, Pergamon, New York, 2002, pp. 82–112 (Chapter 4). [4] A. Barletta, Thermal instabilities in a fluid saturated porous medium, in: A. Öchsner, G.E. Murch (Eds.), Heat Transfer in Multi–Phase Materials, Springer, New York, 2011, pp. 381–414. [5] J.W. Elder, Steady free convection in a porous medium heated from below, J. Fluid Mech. 27 (1967) 29–48. [6] M. Combarnous, Convection naturelle et convection mixte en milieu poreux, ph.d. thesis, University of Paris, 1970. [7] D.J. Close, J.G. Symons, R.F. White, Convective heat transfer in shallow, gasfilled porous media: experimental investigation, Int. J. Heat Mass Transf. 28 (1985) 2371–2378. [8] M.D. Shattuck, R.P. Behringer, G.A. Johnson, J.G. Georgiadisi, Convection and flow in porous media. Part 1. Visualization by magnetic resonance imaging, J. Fluids Mech. 332 (1997) 215–245. [9] L.E. Howle, R.P. Behringer, J.G. Georgiadisi, Convection and flow in porous media. Part 2. Visualization by shadowgraph, J. Fluids Mech. 332 (1997) 247– 262. [10] S.A. Bories, M.A. Combarnous, Natural convection in a sloping porous layer, J. Fluid Mech. 57 (1973) 63–79. [11] D.D. Joseph, Stability of Fluid Motions, Springer, New York, 1976. [12] M. Combarnous, S. Bories, Modelisation de la convection naturelle au sein d’une couche poreuse horizontale a l’aide d’un coefficient de transfert solide– fluide, Int. J. Heat Mass Transf. 17 (1974) 505–515. [13] A. Amiri, K. Vafai, T.M. Kuzay, Effect of boundary conditions on non-Darcian heat transfer through porous media and experimental comparisons, Numer. Heat Transf. A 27 (1995) 651–664. [14] D.A. Nield, Effects of local thermal non-equilibrium in steady convective processes in a saturated porous medium: forced convection in a channel, J. Porous Med. 1 (1998) 181–186. [15] A.V. Kuznetsov, Thermal nonequilibrium forced convection in porous media, in: D.B. Ingham, I. Pop (Eds.), Transport Phenomena in Porous Media, Pergamon, Oxford, 1998, pp. 103–129. [16] B. Alazmi, K. Vafai, Constant wall heat flux boundary conditions in porous media under local thermal non–equilibrium conditions, Int. J. Heat Mass Transf. 45 (2002) 3071–3087. [17] N. Banu, D.A.S. Rees, Onset of Darcy-Bénard convection using a thermal non– equilibrium model, Int. J. Heat Mass Transf. 45 (2002) 2221–2228. [18] A. Postelnicu, D.A.S. Rees, The onset of Darcy-Brinkman convection in a porous layer using a thermal nonequlibrium model – Part I: Stress–free boundaries, Int. J. Energy Res. 27 (2003) 961–973. [19] D.A.S. Rees, I. Pop, Local thermal non-equilibrium in porous medium convection, in: D.B. Ingham, I. Pop (Eds.), Transport Phenomena in Porous Media III, Pergamon, Oxford, 2005, pp. 147–173. [20] A. Postelnicu, Effect of inertia on the onset of mixed convection in a porous layer using a thermal nonequilibrium model, J. Porous Med. 10 (2007) 515– 524. [21] A. Nouri-Borujerdi, A.R. Noghrehabadi, D.A.S. Rees, The effect of local thermal non–equilibrium on conduction in porous channels with a uniform heat source, Transp. Porous Med. 69 (2007) 281–288. [22] A. Nouri-Borujerdi, A.R. Noghrehabadi, D.A.S. Rees, Onset of convection in a horizontal porous channel with uniform heat generation using a thermal nonequilibrium model, Transp. Porous Med. 69 (2007) 343–357. [23] A. Postelnicu, The onset of a Darcy-Brinkman convection using a thermal nonequilibrium model. Part II, Int. J. Therm. Sci. 47 (2008) 1587–1594. [24] K. Yang, K. Vafai, Analysis of temperature gradient bifurcation in porous media – an exact solution, Int. J. Heat Mass Transf. 53 (2010) 4316–4325. [25] A. Postelnicu, The effect of a horizontal pressure gradient on the onset of a Darcy–Bénard convection in thermal non-equilibrium conditions, Int. J. Heat Mass Transf. 53 (2010) 68–75. [26] D.A.S. Rees, Microscopic modeling of the two-temperature model for conduction in heterogeneous media, J. Porous Med. 13 (2010) 125–143.

410

M.N. Ouarzazi et al. / International Journal of Heat and Mass Transfer 113 (2017) 399–410

[27] A. Barletta, D.A.S. Rees, Local thermal non-equilibrium effects in the DarcyBénard instability with isoflux boundary conditions, Int. J. Heat Mass Transf. 55 (2012) 384–394. [28] J.P. Caltagirone, S. Bories, Solutions and stability criteria of natural convective flow in an inclined porous layer, J. Fluid Mech. 155 (1985) 267–287. [29] D.A. Nield, A note on convection patterns in an inclined porous layer, Transp. Porous Med. 86 (2011) 23–25. [30] L. Storesletten, M. Tveitereid, Onset of convection in an inclined porous layer with anisotropic permeability, Appl. Mech. Eng. 4 (1999) 575– 587. [31] M. Karimi-Fard, M.C. Charrier-Mojtabi, A. Mojtabi, Onset of stationary and oscillatory convection in a tilted porous cavity saturated with a binary fluid: linear stability analysis, Phys. Fluids 11 (1999) 1346–1358. [32] D.A.S. Rees, A.P. Bassom, Onset of Darcy-Bénard convection in an inclined layer heated from below, Acta Mech. 144 (2000) 103–118.

[33] D.A.S. Rees, L. Storesletten, A. Postelnicu, The onset of convection in an inclined anisotropic porous layer with oblique principle axes, Transp. Porous Med. 62 (2006) 139–156. [34] A. Barletta, D.A.S. Rees, Local thermal non-equilibrium analysis of the thermoconvective instability in an inclined porous layer, Int. J. Heat Mass Transf. 83 (2015) 327–336. [35] M.N. Ouarzazi, F. Mejni, A. Delache, G. Labrosse, Nonlinear global modes in inhomogeneous mixed convection flows in porous media, J. Fluid Mech. 595 (2008) 367–377. [36] A. Delache, M. Ouarzazi, M. Combarnous, Int. J. Heat Mass Transf. 1485. [37] A. Delache, M. Ouarzazi, Weakly nonlinear interaction of mixed convection patterns in porous media heated from below, Int. J. Therm. Sci. 47 (2008) 709– 722. [38] S.C. Generalis, K. Fujimura, Range of validity of weakly nonlinear theory in the Rayleigh-Bénard problem, J. Phys. Soc. Jpn. 78 (2009) 8–084401.