Onset of Thermal Instability in a Fluid-Saturated Porous Layer with Internal Heat Generation and Density Maximum.

Onset of Thermal Instability in a Fluid-Saturated Porous Layer with Internal Heat Generation and Density Maximum.

N. K. Enagi (1)* , Sridhar Kulkarni (2)

(1) KRCE Societys GGD Arts, BMP Commerce and SVS Science Degree College Bailhongal,Karnataka, India

(2) Government First Grade College, Paschapur-591122, Karnataka,India.

Abstract – The onset of thermal instability in a fluid-saturated porous layer subject to internal heat generation and a density maximum is investigated analytically under non-local thermal equilibrium (NLTE) conditions. A two-field model is employed to describe the thermal behavior of the porous medium, with separate energy equations formulated for the fluid and solid phases. Linear stability analysis is performed using the normal mode technique to determine the conditions governing the onset of thermal instability. Galerkins method is subsequently employed to obtain the Eigen values of the resulting system. The influence of internal heat generation and the thermal conductivity ratio on the stability characteristics of the porous layer is examined. The results demonstrate how these parameters modify the critical conditions for the onset of thermal instability and provide insight into the thermal behavior of porous media with density- dependent fluid properties.

Key Words: Convection, density maximum, internalheat generation, thermal non-equilibrium.

Nomenclature

kf Thermal conductivity tensor of the fluid phase.

ks Thermal conductivity tensor of the solid

phase.

(x, y, z) Space coordinates

Greek symbols

a Diffusive ratio

fJ Coefficient of thermal expansion

y

Ekf

Porosity-modified conductivity ratio(1-E)K

fJ1p0gd(pc)f(Tf – Tu)K/Eµkf

H Inter-phase heat transfer coefficient.

Qf Internal heat generator for fluid phase. Qs Internal heat generator for solid phase R2 Maximum density function.

RA Rayleigh number,

RD Maximum density function,

2

fJ2p0gd(pc)f(Tf – Tu) K/Eµkf

E Porosity

0 Non-dimensional temperature of fluid phase

µ Fluid viscosity

p Fluid density

<p Nondimensional temperature of solid phase

1/ Stream function

w Growth rate

Other symbols

K Permeability tensor K(ii + jj) + Kzkk

T Temperature

2 iJ2

V

1 iJx2

+ iJ2

iJy2

a Horizontal wave number, z2 + m2 C Specific heat

d Height of the porous layer

g Gravitational acceleration,

z, m Wave numbers in the x and y direction respectively

1. Pressure

2. Velocity vector, (u, v, w) t Time

V2 V2 + az

1 z

aZ

Subscript/ Superscripts

* Dimensionless quantity

Perturbed quantity

0 Reference

f Fluid phase

l Lower

s Solid

1. Introduction.

   The study of thermal convection in fluid-saturated porous media has received considerable attention because of its applications in engineering and natural systems. Early investigations mainly focused on the onset of convection under thermal non-equilibrium and non-Newtonian effects. Shivakumara et al. [1] studied convection in a visco-elastic fluid-saturated porous layer using the local thermal non-equilibrium (LTNE) model. Postelnicu [2] subsequently examined the effect of inertia on mixed convection under LTNE conditions. Sarvanan [3] and Shivakumara et al. [4] extended these studies by incorporating internal heat generation, density maximum, and variable viscosity.

   Subsequent studies considered additional physical effects. Israel C et al. [5] investigated thermal instability in a low-Prandtl-number fluid with internal heat generation, while Lee et al. [6] and Shivakumara et al. [7] examined ferro convection in porous media using thermal non-equilibrium models. Bhadauria et al. [8] studied the influence of internal heating on nonlinear convection under g-jitter. Later, Yadav et al. [9] investigated DarcyBrinkman convection in a rotating porous layer with internal heating.

   From 2019 onward, increasing attention was given to the combined effects of rotation, anisotropy, non- Newtonian fluids, and internal heating. Chavaraddi et al. [10] studied couple-stress-fluid convection in rotating anisotropic porous layers under LTNE conditions, while Lotten and Rees [11] considered an inclined anisotropic porous layer with internal heat generation. Recent developments have also included fractal porous structures and microscale transport effects [12,13]. Enagi et al. [14,15] investigated rotating porous layers with internal heating, density maximum, and couple-stress effects. More recently, Enagi et al.

   [16] examined the effect of anisotropy on DarcyBrinkman convection in a Maxwell-fluid-saturated porous layer.

   Thus, the research has progressively evolved from basic porous convection and LTNE effects to more complex combinations involving internal heating, rotation, anisotropy, density variation, and non-Newtonian fluid behavior. These developments motivate further investigation of convective instability in porous media involving the combined effects of these parameters.

2. Mathematical formulation

   Consider a porous matrix of height (h), extending horizontally between two free surfaces. The lower and upper surfaces are maintained at temperatures Tf and Tu respectively, with the lower surface being heated and the upper surface being cooled. Thus, Tf > Tu , A uniform temperature difference is maintained between the two surfaces, producing a constant vertical temperature gradient, VT = Tf Tu. The mathematical model consists of a momentum equation that includes the time-derivative term, together with two separate temperature equations. The physical configuration considered in this study is illustrated schematically in Figure 1.

   Z

   = (0, 0, )

   z=d

   Tu

   z=0

   Tl

   X

   Heated from below

   Y

   Figure 1: Geometrical Representation of physical problem

   Internal heat source

Porous layer

g = (0, 0, -g)

Height of the fluid layer

Thus, the corresponding governing equations are

V ·q=0 (1)

g –

1 aq = – 1 Vp + p V q (2)

E at po po k

aTf

E(p ) + (p ) (q V)T

= Ek V2T + (T – T )+ Eq

(3)

c f at c f

f f f s f f

(1 – E)(p ) aTs = (1 – E)k V2T – (T – T ) + (1 – E)q

(4)

c s at

s s s f s

2

p = p0[1 – fJ1(Tf – Tu) – fJ2(Tf – Tu ) (5)

To remove the pressure term from the momentum equation and making equations (3) & (4) dimensionless by using following transformations

(x, y, z) = d(x*, y*,z*) (u, v, w) = Ekf

(pc)fd

(u*,v*,w*), p = kfµ

(pc )f K

p*,

( ) *

( ) *

(pc)fdz

* (6)

Tf =

Tl – Tu Tf + Tu , Ts =

Tl – Tu Ts

+ Tu , t = t

k

f

1 a 2

2 2 2 2

azw

(V w) = R1V1 Tf + R2V1 Tf – V1 w – (7)

PrD at aZz

aTf azTf 2 azTf (

(8)

at + w aZz = V1 Tf + aZz + H Ts – Tf ) + Qf

a aTs = V 2T + azTs – yH(T – T ) + Q

(9)

at 1 s

aZz

s f s

RA = fJ1p0gd(pc )f(Tl – Tu)K/Eµkf , the Rayleigh number corresponding to the properties of fluid phase.

RD = fJ2p0gd(pc )f(Tl – Tu )2K/Eµkf , as a measure of the density maximum.

H = d2/(1 – E)ks the non-dimensional interphase heat transfer coefficent.

y = Ekf/(1 – E)ks, the conductivity ratio a = (pc)skf/(pc)fks , the diffusivity ratio.

Qf = qf/(Tl – Tu )kf , fluid phase internal heat generator parameters,

Qs = d2qs/(Tl – Tu)ks solid phase internal heat generator parameters

2.2 Basic State.

The basic state is assumed to be quiescent and is given by

u = v = w = 0 Tf = Tfb (z) Ts = Tsb(z) (10)

The basic state temperature of fluid phase and solid phase satisfy the equations

dzTfb

dZz + H(Tsb – Tfb ) = 0 ,

dzTsb – yH(T – T ) = 0 (11)

z sb fb

dZ

, with the boundary conditions Tfb = Tsb = 1 at z = 0 , Tfb = Tsb = 0 at z = 1 (12) So that the conduction state solutions are given by

T Qf 2 Qf

Qs 2 Qs

fb = – 2 z

+ ( 2 – 1) z + 1, Tsb = – 2 z

+ ( – 1) z + 1 (13)

2

2.3. The perturbed state:

The basic state is perturbed and quantities in the perturbed state are given by

(u, v, w) = (u1, v1, w1), q = q1, Tf = Tfb + 0, Ts = Tsb + <p (14)

z

Substituting equation (14) into (7) to (9) and using the basic state solutions, we obtained following linear equations for perturbed quantities (after neglecting the primes)

1 a (V2w) = R V20 + 2R

{-Q Z Qf

} 2 2

azw

(15)

PrD at

1 1 2

f 2 + ( 2 – 1) z + 1

V1 0 – V1 w – aZz

ae Qf

2 aze + H(<p- 0) (16)

aZ

at + w (-Qfz + 2 – 1) = V1 0 + z

a a<p

2 az<p – yH(<p – 0) (17)

at = V1 <p+ aZz

Since the fluid and solid phases are not in thermal equilibrium, the use of appropriate thermal boundary condition may pose a difficulty. However, the assumption that the solid and fluid phases have equal temperatures at the boundary surfaces made at the beginning of this section helps in overcoming this difficulty. Accordingly, equations (15 to 18) are solved impermeable isothermal boundaries. Hence the boundary conditions are

w = aw = 0 = <p = 0 at z = 0, 1 (18)

aZ

2.4. Linearity Properties.

To study the linear stability theory, we use the linearized version of equations (16) to (18).The principle of exchange of stabilities holds in the presence of isotropy and non-LTE effects (there is only one destabilizing agency) so that the onset of convection is stationary (ie w =0).We seek the solutions to the linearized equations in the form

[w, 0, <p] = [W(z), 0(z), <P(z)]ei(lx+my) sin nZ+wt (19) Here l, m is the wave numbers in horizontal plane and wis growth rate. Infinitesimal perturbation of the rest state may be either damp or grow depending on the values of the parameter w. Substituting the equations

(19) in the equations (16 to 17) we get the following equations

w 2 2 2 2

[ (a – D ) – ((a – D )] W – [R

+ 2R

{- Qf 2 Qf

}] 2

(20)

PrD

1 2 z

+ ( – 1) z + 1

2

2

a 0 = 0

(-Q z + Qf + 1) W + (w + a2 + H – D2)0 – H<P = 0 (21)

f 2

yH0 + (aw + a2 + yH – D2)<P = 0 (22)

We use Galerkins technique to solve the Eigen value problem. The trail function is

W1 = z3(1 – z)20 = <P = z(1 – z)and integrating the equations so obtained over the layer from 0 to 1 we get

[az + 2] 1

Q +14

[ f 2

– +R { }] a B

= 0 (23)

66 9 1 2

2 36 1

(18 – Qf)11 – 188(a2 + H + 10)B1 + 188H11 = 0 (24)

-yHB1 + (a2 + yH + 10)11 = 0 (25)

Now to solve RA the above equation can be put in the form of the following matrix, we get

az + 2

Qf Qf

– [R + R {- 2

}] 2 0

11 0

66 9

1 2 z

+ ( + 1) z + 1 a

2

2

18 + Qf -188(a2 + H + 10) 188H

0 -yH a2 + yH + 10

B1 = 0 I I

l J

(26)

l11J 0

By setting the determinant of the coefficient matrix to zero we get

9

R = 448(az+10) [az + 2] [1 + H ] 2

A az(18+Qf) 22 3

az+yH+10

– R2(Qf + 14) (27)

1. Results and discussion

   The Rayleigh number given by Equation (27) is evaluated numerically for various values of the governing physical parameters to investigate their effects on the onset of convection. In particular, the influence of internal heat generation and the maximum density function in a saturated porous layer is examined using linear stability theory. Since oscillatory convection is considered highly improbable under the conditions investigated, the present study focuses exclusively on the onset of convection in the stationary mode.

   Figures 1(a)1(c) illustrate the neutral stability curves for different values of the inter-phase heat transfer parameter (H), the maximum density function parameter (R2), and the conductivity ratio ( ),

   respectively. The neutral curves are topologically connected, allowing the stability criterion to be expressed in terms of a critical Rayleigh number. The system remains stable when the Rayleigh number is below its critical value and becomes unstable when it exceeds this threshold.

   Figure 1(a) shows that the minimum critical Rayleigh number increases with increasing (H), indicating that a higher inter-phase heat transfer coefficient delays the onset of convection and thus has a stabilizing effect on the system. In contrast, Figures 1(b) and 1(c) demonstrate the effects of the maximum density function parameter (R2) and the conductivity ratio ( ) on the neutral stability curves.

   The results indicate that the minimum critical Rayleigh number decreases as the internal heat generation, maximum density function, and conductivity ratio increase. This reduction in the critical Rayleigh number implies that these parameters promote the onset of convection and therefore exert a destabilizing influence on the system.

   H=100

   H=50

   H=20

   H=10 H=5

   2

   f

   400

   350

   300

   R

   A

   250

   200

   480

   Q =1

   f

   Q =5

   f

   Q =10

   f

   Q =15

   f

   Q =20

   f

   2

   400

   R

   A

   320

   240

   160

   150	y=0.1, R =10,	Q =5			80	y=0.1, R =10, H=50
   1000	4 8 a	12	16	20	0	4 8 12 a	16

   Figure 1(a):Neutral stability curves Vs ‘a’ for different values of H

   Figure 1(b)::Neutral stability curves Vs ‘a’ for different values of Qf

   R =5

   R =15

   2

   2

   R =30

   2

   R =50

   2

   R =100

   2

   y=0.1, Q =5, H=50

   f

   400

   320

   R

   A

   240

   160

   80

   0 4 8 12 16

   a

   2

   Figure 1(c):Neutral stability curves Vs ‘a’ for different values of R

   Figures 2(a) and 2(b) illustrate the variation of the critical Rayleigh number with the inter-phase heat transfer coefficient (H) for different values of the porosity-modified conductivity ratio ( y ) and the internal heat generation parameter (Qf), respectively, while keeping the other governing parameters fixed.

   Figure 2(a) shows the effect of the porosity-modified conductivity ratio ( y ) on the critical Rayleigh number as (H) varies. It is observed that, as (H) approaches zero, the critical Rayleigh number becomes independent of ( y ) and approaches the corresponding value for the local thermal equilibrium (LTE) case. Furthermore, the critical Rayleigh number decreases with increasing ( y ), indicating that an increase in the porosity- modified conductivity ratio has a destabilizing effect on the system and promotes the onset of convection. It is also observed that, for sufficiently large values of ( y ), the critical Rayleigh number becomes independent of (H), suggesting tht the influence of inter-phase heat transfer on the onset of convection becomes negligible in this parameter regime.

   Figure 2(b) depicts the variation of the critical Rayleigh number with the inter-phase heat transfer coefficient (H) for different values of the internal heat generation parameter (Qf). The results indicate that the critical Rayleigh number decreases as (Qf) increases. This behavior implies that internal heat generation enhances the destabilization of the system, thereby facilitating the onset of convection.

   Q =10, R =20

   f D

   y=0.1

   y=0.3 y=0.5

   y=0.7 y=1

   2.0k

   1.5k

   R

   A

   1.0k

   500.0

   30k

   R =20, y=0.01

   2

   Qf=0 Qf=2

   Qf=5 Qf=10

   Qf=20

   R

   A

   20k

   10k

   0.0

   0 1 2 3 4 5

   00 1 2 3 4 5 6

   H

   10

   log H

   log

   10

   Figure 2 (a): Variation of R

   Figure 2 (b): Variation of RA

   against

   f

   A

   against H for different values of y

   H for different values of Q

   Figure 3(a) illustrates the variation of the critical wave number with the inter-phase heat transfer coefficient (H) for different values of the conductivity ratio ( ). The results reveal that the critical wave

   number is independent of Figure 3(a) illustrates the variation of the critical wave number with the inter- phase heat transfer coefficient (H) for different values of the conductivity ratio (\gamma). The results reveal that the critical wave number is independent of ( ) for both small and large values of (H). However, at

   intermediate values of (H), the critical wave number increases as ( ) decreases, reaching a maximum

   before subsequently decreasing. These findings indicate that the influence of the conductivity ratio on the critical wave number is most pronounced at intermediate values of the inter-phase heat transfer coefficient. for both small and large values of (H). However, at intermediate values of (H), the critical wave number increases as ( ) decreases, reaching a maximum before subsequently decreasing. These findings indicate

   that the influence of the conductivity ratio on the critical wave number is most pronounced at intermediate values of the inter-phase heat transfer coefficient.

   y =0 .1

   y =0.3

   y=0 .5

   y=0 .7

   y=1

   Q = 10, R =2 0

   f

   D

   14

   12

   a

   c

   10

   8

   6

   0 2 4 6

   lo g 10 H

   F ig 3 ( a ) :V a r ia tion o f ‘a ‘ a g ain s t ‘H ‘ f o r d i f f ere n t v a lu es o f y .

   Figures 4(a)4(c) illustrate the variation of the critical Rayleigh number with the internal heat generation parameter (Qf) for different values of the conductivity ratio ( ) the inter-phase heat transfer

   coefficient (H), and the maximum density function parameter (R2), respectively.

   Figure 4(a) demonstrates the influence of the conductivity ratio ( ) on the critical Rayleigh

   number. The results indicate that the critical Rayleigh number decreases with increasing (\gamma), suggesting that a higher conductivity ratio has a destabilizing effect on the system and promotes the onset of convection.

   Figure 4(b) depicts the variation of the critical Rayleigh number with (Qf) for different values of the inter-phase heat transfer coefficient (H). It is observed that the critical Rayleigh number decreases as the internal heat generation parameter increases and reaches zero at a finite value of (Qf).This behavior indicates that internal heat generation progressively destabilizes the system by reducing the threshold required for the onset of convection. When the critical Rayleigh number reaches zero, convection can arise without an additional externally imposed buoyancy force, under the assumptions of the model.

   2

   Figure 4(c) illustrates the effect of the maximum density function parameter (R2) on the critical Rayleigh number. The results show that the critical Rayleigh number decreases with increasing (R2), indicating that an increase in the maximum density function has a destabilizing influence on the system and facilitates the onset of convection.

   R =20, H=100

   2

   y=0.01

   y=0.1

   y=1 y=2 y=5

   600

   450

   R

   A

   300

   150

   600

   R

   A

   400

   200

   H=100

   H=50

   H=20 H=10 H=0

   R =20, y=0.01

   0-1 0 1 2 3 0

   Q

   10

   f

   log

   -1 0 1 2

   Q

   10

   f

   Q

   f

   log

   A

   f

   Fig 4 (a): Plots of R versus Q for different values of y

   Fig 4 (b): Plots of RA versus for different values of H

   R =0

   600 2 R2=5

   2

   R =20

   400

   R

   A

   2

   R =50

   R =10

   2

   y=0.01,H=100

   200

   Q

   0 0 2 4

   10

   f

   log

   Fig 4 (c) Plots of RAversus Qf for different R2

2. Conclusion.

The present study investigates the linear stability porous layer saturated with a, incorporating internal heat generation and a density maximum, while accounting for the temperature difference between the fluid and solid phases. The eigen value problem governing the onset of convection is solved using the Galerkin method. Particular emphasis is placed on examining the effects of internal heat generation, the conductivity ratio, and the maximum density function on the stability characteristics of the system. The results reveal that internal heat generation significantly reduces the critical Rayleigh number, thereby lowering the threshold for the onset of convection and promoting thermal instability. Similarly, an increase in the conductivity ratio and the maximum density function parameter decreases the critical Rayleigh number, indicating their destabilizing influence on the system. These findings demonstrate that internal heat generation, the conductivity ratio, and the maximum density function play important roles in governing the onset of convection in rotating couple-stress fluid-saturated porous media. The study provides valuable insight into the thermal stability of such systems and highlights the importance of these parameters in predicting and controlling convective instability.

References

[1]. I.S. Shivakumara, M.S. Malashetty, K.B. Chavaraddi, Onset of convection in a visco-elastic-fluid- saturated sparsely packed porous layer using a thermal non equilibrium model, Can. J. Phys. 84(11) (2006) 973990, https://doi.org/10.1139/p06-085

[2]. Postelnicu A.: Effect of Inertia on the onset of mixed convection in a porous layer using the thermal non-equilibrium model. J Porous Media 10 (5), 515-24 (2007). https://doi.org/10.1615/JPorMedia.v10.i5.80

[3]. S. Sarvanan, Thermal non-equilibrium porous medium with heat generation and density maximum, Transp. Porous Media 76 (2009) 3543. https://doi.org/10.1007/s11242-008-9232-x

[4]. I.S. Shivakumara, A.L. Mamatha, M. Ravisha, Effects of variable viscosity and density maximum on the onset of Darcy-Benard convection using a thermal non equilibrium model, J. Porous Media 13(7) (2010) 613622., DOI:10.1615/JPorMedia.v13.i7.30

[5]. Israel C. Cookey et al.:Onset of thermal instability in a low Prandtl number fluid with an internal heat source in porous media. American jou of scientific & Industrial Research 1 (1),18-24 (2010). https://www.researchgate.net/publication/265694080_

[6]. Lee, I.S. Shivakumara, M. Ravisha, Effect of thermal non-equilibrium on convective instability in a ferromagnetic fluid saturated porous medium, Transp. Porous Media 86 (2011) 103124. Springer Nature Link

[7]. I.S. Shivakumara, Jinho Lee, M. Ravisha, R. Gangadhara Reddy, The onset of Brinkman ferro convection using a thermal non-equilibrium model, Int. J. Heat Mass Transf. 54 (2011) 21162125.. https://doi.org/10.1016/j.ijheatmasstranser.2010.12.016

[8]. B.S. Bhadauria, I. Hashim, P.G. Siddeshwar, Effect of internal heating on weakly nonlinear stability analysis of Rayleigh-Benard convection under g-jitter, Int. J. Non-Linear Mech. 54 (2013) 3542 https://doi.org/10.1016/j.ijnonlinmec.2013.03.001

[9]. D. Yadav, J. Wang, Jinho Lee, Onset of DarcyBrinkman convection in a rotating porous layer induced by purely internal heating, J. Porous Media 20(8) (2017) 691706.

DOI:10.1615/JPorMedia.v20.i8.20

[10]. K.B. Chavaraddi, N.K. Enagi, S. Kulkarni, Onset of convection in a couple stress fluid saturated rotating anisotropic porous layers using thermal non-equilibrium model, J. P. J. Heat Mass Transf. 16(1) (2019) 125142. http://dx.doi.org/10.17654/HM016010125.

[11]. L.S. Lotten, D.A.S. Rees, Onset of convection in an inclined anisotropic porous layer with internal heat generation, Fluids 4 (2019) 75. https://doi.org/10.3390/fluids4020075

[12]. C.Q. Wu, H.J. Xu, C.Y. Zhao, A new fractal model on fluid flow/heat/mass transport in complex porous structures, Int. J. Heat Mass Transf. 162 (2020) 120292.

[13]. H.J. Xu, Thermal transport in micro channels partially filled with micro-porous media involving flow inertia, flow/thermal slips, thermal non-equilibrium and thermal asymmetry, Int. Commun. Heat Mass Transf. 110 (2020) 104404. https://doi.org/10.1016/j.icheatmasstransfer.2019.104404

[14]. N. K. Enagi, K. B. Chavaraddi, S. Kulkarni and G. K. Ramesh, Effect of maximum density and internal heating on the stability of rotating fluid saturated porous layer using LTNE model, Heliyon 8(6) (2022), e09620., https://doi.org/10.1016/j.heliyon.2022.e09620.

[15]. N. K. Enagi, Krishna B. Chavaraddi, Sridhar Kulkarni, Couple stress fluid saturated rotating porous layer with internal heat generation and density maximum, JP Journal of Heat and Mass Transfer 35(1) (2023), 1-19. http://dx.doi.org/10.17654/0973576323039.

[16]. N. K. Enagi, Krishna B. Chavaraddi and Sridhar Kulkarni, The effect of anisotropy on Darcy- Brinkman convection in a Maxwell fluid saturated porous layer, Advances and Applications in Fluid Mechanics 31(1) (2024), 1-22. https://doi.org/10.17654/0973468624001.