Proc. Indian Acad. Sci. (Math. Sci.), Vol. 96, No. 1, August 1987. pp. 1-25. '.(5"Printed in India.
Onset of convection in a porous mantle
G K PRADHAN and B PATRA Department of Mathematics, University College of Engineering. Burla 768 018, India MS received 27 October 1984
Abstract. The problem of thermal instability in a fluidsaturated porous spherical shellheated internally, due to uniform internal heat sources and in equilibrium under its own radial gravitational field is studied theoretically. A general disturbance is analysed into modes in terms of spherical harmonics of various orders, l, for different values of the thickness of the mantle and the criteria for the onset of convection for the first fifteen modes is obtained in four different cases when the outer and inner bounding surfaces are either impermeable or permeable. It is shown that as the thickness of the shell decreases, the pattern of convection which sets in at marginal stability shifts progressively to harmonics of higher order for all the three cases except when both the bounding surfaces are permeable, in which case the onset of convection occurs at a harmonic of order 1. A comparison of some representative results of these cases is made with that of continuous fluid shell with rigid or free boundary surfaces. The neutral stability plots for various thickness of the mantle, for five different models of the mantle, are plotted for the different types of boundary surfaces. Keywords. Porous mantle; mantle convection; Darcy-Rayleigh number; Darcy-Prandtl number; geophysical model.
1. Introduction The problem of convection in a spherical shell heated from within has attracted widespread a t t e n t i o n a m o n g geophysicists a n d astrophysicists interested in convective process. H e i s k a n e n a n d Vening Meinesz [9] suggested that sometime in the early history of the earth cooling, great current system existed in the mantle, chiefly due to temperature gradients caused by cooling at the surface. According to the presently available picture of the earth, its central core a b o u t 7,000 k m diameter, is composed of solid nickel and iron at temperature of nearly 5000~ a n d at pressure as high as 3.8 million kg/scm. W r a p p e d a r o u n d this core is the mantle, a 3000 km thick area of molten a n d semisolid rocks a n d a r o u n d this is the crust averaging 50 km of thickness and holding the oceans a n d continents. As a result, it has been suggested by Simpkins a n d Blythe [14] that thermal convection problems in porous media have application to convection in the mantle of the earth. The p r o b l e m of m a n t l e convection is discussed here.
2. Model formulation We consider a spherical shell confined between two concentric spheres of radii R 1 and R 2 c o n t a i n i n g fluid saturated porous matter, in e q u i l i b r i u m u n d e r its own
2
G K Pradhan and B Patra
gravitation and with a uniform distribution of heat sources such that in the absence of convection the temperature at each point will rise at a rate H,.. The solid incompressible substrate has a constant porosity ~ and a permeability K. The interstitial fluid is characterized by its thermal volume expansion a, density Pl, heat capacity (pc); and kinematic viscosity v. The ensemble formed by the solid substrate and the fluid is treated as a unique fictitious fluid of thermal conductivity k,. and heat capacity (p~),. = g(Pc)I + (l - g)(p~)~. This hypothesis has been extensively studied by Combarnous and Bories [7] and [8] for the case of a horizontal porous layer. In the steady state, the temperature distribution T(r) inside the spherical shell is given by VT = - 2/~(r)r
(1)
where r is the position vector with respect to the centre of the sphere,
fl(r) = flz + (flt/2r3),
f12 = -- (Hm/6krh'f)
and fll is a constant.
The equation of continuity, Darcy equation of motion and energy for a differential volume element of the fluid-porous medium are given by V-u=O
(2)
(Su/?t) + u'Vu = - V(6p/p,,) + 7(r)Or - (&'/K)u
(3)
C(80/?t) + gu-VO = k, KTV20 + 2gg(r)u.r.
(4)
The momentum equation (3) is similar to that used by Lapwood [11]. Here u is the disturbance velocity, 3p the perturbation in pressure, C is the ratio of heat capacities (Pc),./(Pc)f, kr the thermal conductivity ratio k,,/kr, k,, being the thermal conductivity of the solid-fluid mixture given by gkf + ( 1 - Qk~, f denoting the fluid value, s the solid value, m the solid-fluid mixture and Kf the thermal diffusivity of the fluid. In deriving the equation of motion and energy, we have assumed the D a r c y - O b e r b e c k Boussinesq approximation and have replaced pf occurring with the gravity term 9(r) by pf = p,,(1 - ctO) and pf occurring elsewhere as a constant equal to p,,, the density of the saturated medium; the perturbation density 6pl and that in the temperature O being related by 6pl = - pi~O and 7(r) = ~tg(r).
2.1. Eigenvalue problem for spherical shell To investigate the necessary criteria for the onset of thermal instability, we need only consider infinitesimal disturbances. Following Chandrasekhar [5], we apply the operator r.VXVX( ) to 13) and write the result in spherical polar coordinates (r, 0, ~b):
(K + ~)V2(u'r)=
-7(r)L20
(5)
where 2
L2 = r 2 1 -
__
~ . ~? ,, ~. s,n 0~:;,
slnuf't/
(o
1
~2
sin2 0 ~qq~2
(6)
Onset of convection in a porous mantle
3
and V 2 is usual Laplacian in spherical polar form. Restricting our attention to the linear problem describing small amplitude convection, we neglect the non-linear terms u-Vu and u.V| occurring in (3) and (4). The dependent variables u.r and | occurring in (4) and (5) are now expanded in a complete set of spherical harmonics u-r = W(r) YT'(O, qS)exp(pt) | = |
(7)
YT'(O, r exp(pt)
where r?(o, q~) = PT'(cos O)exp ( + Im r
(8)
and PT'(cos 0) are the associated Legendre polynomials, p being a constant. Substituting (7) and (8) in (4) and (5) with m = 0 (for axisymmetric disturbances) yield (1 + (K/gv)p)~ 1W(r) = - l(l + 1)(K/gv)7(r)|
(9)
and (k, x i ~ ~ - pC)O(r) = - 2~fl(r) W(r)
(1 O)
where
l(t+ l)
91 = (dZ/dr 2) + (2/r)(d/dr)
r2
(11)
We now measure r in units of a suitable radius R1, the outer radius of the shell and denote a = (K/&)p and rewrite (9) and (10) as K
2
(1 + t~)~l W = - l ( l + 1)7(r)c-Rl|
(12)
(cJ 1 - B C a ) O = - 2fl(r) k,~iR~ W
(13)
8!,'
and
where the Darcy Prandtl number B - l is defined by
()4)
K k, K: " 2.2. Boundary conditions
We 'use spherical polar coordinates (r, 0, r with origin at the centre of the core and 0 as colatitude and confine the shell between radii r = q and r = 1, the unit of length being chosen as the outer radius R 1 of the shell. The spherical boundaries are assumed to be either impermeable or else they are bounded by standing fluid identical with that saturating the medium. In both the cases the boundary conditions are |
forr=l
andr=q.
(15)
Further, on the impermeable boundaries, we have W=0
forr=l
and r = q
(16)
and on permeable boundaries 6p=O
forr=l
andr=q.
117)
4
G K Pradhan and B Patra
Considering the implication of imposing condition (17) on spherical boundaries from the 0- and @components of equation of motion, we deduce on making use of the incompressibility condition (2) that
(dW/dr) + (W/r) = 0
for r = 1 and r = r/
(18)
for permeable boundaries.
2.3. Exchange of stabilities for the case fl = constant and 7= constant With the restriction fl(r) = fl and 7(rl = i' and substituting K
F = -I(t + I)';~R2|
(19)
,gV
(12) and (13) reduce to (1 + o ) ~ W = F
(20)
( ~ -- BCajF = 1(1+ 1)ROW
(21)
2fl7 K ~4 Ro = ~t~ 1.
(22)
and where
The modified boundary conditions are now F=0
forr=l
andr=q
(23)
W = 0 for r = 1 and r = ~; on impermeable boundaries and
(dW/dr)+(W/r)=O
forr=l
andr=q
(24)
on permeable boundaries. We now multiply (21} by r2F *, superscript * denoting complex conjugate of the quantity, integrate the result over the r,~nge (r/, 1) and utilise (20) to obtain
dr
+ I(I + I)IFIZ + BCarelFI 2 dr
- d W ' I n - j n f l ~(r 2 d[~-r W z [ +I(t+ l)lWl z } dr ] = 0 . +l(l+ l)RD(l+tr*) [ rew-~-r (25)
In deriving (25) we have made use of the property of the operator ~
fnlr20~t~/dr=r2odl~ri- J,fl~r~ 2d~bd-~dr dr + l(l +
l)O~')dr,
(26)
where O(rl and elf) are any two bounded and continuous functions in (r/, 1) and the boundary conditions given by (23).
Onset of convection in a porous mantle
5
On the impermeable boundary W = 0 , the integrated part [r2W(dW*/dr)]~ occurring in (25) vanishes but if the boundary is permeable, then {d W*/dr) + (W*/r) = O, so that
rzwdW* l dr _], = [ - r l W l z ] ~ . Hence in either case the integrated part can be written in the form [-rlWl2]~.
127)
Vanishing of the imaginary part of (25) with the aid of (27) gives
lmla)IBCf.~rZlF'2dr+l(l+l)R~
r'~dW
~-I(1+ 1)1W,2) dr}
+ Jr[ B"IZ]~ ] 0.
(28)
For impermeable boundaries, the factor of Ira(a) in (28) is positive definite, and accordingly lmla) = 0
(29)
for impermeable boundaries. 2.4. Variational principle when fl and 7 are constants On the assumption of validity of P.E.S. when fl and 7 are constants, the equations governing the marginal state are given by c/~ W = F
(30)
~ F = 1(1+ I)RDW.
(31)
and These two equations together with boundary conditions (23) and (24) constitute a characteristic value problem for R e and the instability will set in with that value of I which leads to the lowest value of R o. The underlying characteristic value problem can be formulated in terms of a variational principle, the basis of which follows from (25), which, when a = 0 gives
1(l + I)R o
=
II~r2(dF) 2 } J" ( \ dr + 1(t + llF'- dr [rW2]~" J"+
(
\dr
/
.
I32)
+I(I+I)W 2 dr
It can be shown as in Chandrasekhar [5], that the lowest characteristic value for is the absolute minimum of the quantity on the right hand side of 132).
RD
3. Onset of thermal convection
Assuming the validity of principle of exchange of stabilities, the equations governing the marginal state (112) and (13)l are
6
G K Pradhan and B Patra ~ W = - l(l + 1)7(r)KR20
(33)
8V
and A
(34)
~ O = - 2fl(r)k~-R2 W where R 1 is the radius of outer spherical boundary. With the substitution ?(r) =)hc(r)
(35)
F = - l(l + 1)',,lc(r) ~K R20
(36)
fl(r) = fllb(r),
and
where fll and 71 are the values of fl(r) and ),(r) at r = 1, (33) and (34) are rewritten in the form cd~w = F
(37)
and b~t
=
(38)
l(l + 1 ) R O W
where the dimensionless parameter occurring in (38) is the Darcy-Rayleigh number Ro -
2fl t';1K krN f V
R~
(39)
which is defined in terms of the strength of the energy source H,,, the only means by which destabilising density differences are produced in the fluid. Solution of (37) and (38) must be sought which satisfy the boundary conditions (23) and (24). To solve the characteristic value problem by the variational method and determine the onset of instability as a spherical harmonic of order l, we expand F in a series of cylinder functions of order l + 89which vanish at r = 1 and r = )7. This expansion is possible because of the cylindrical symmetry (m = 0) of the problem. Let Ct+ 1 2.~(z) = J - ( t +
x2)(~xq)J,(z)
- Jr+ x 2(~q)J-,.(z)
(40)
where ~ is a constant not specified for the present. Now CI + i/2.1 § i / 2 ( ~ r ) = J - , i + I:2)(~xPI)JI + 1:2( ~r} - Jl + 1 2(~r/) J -(l + 12)(:~r)
(41)
clearly vanishes for r = r/and it will also vanish for r = I, provided J - a + I/2)(~
I/2(~x) - J r + i ; 2 ( ~ ( q ) J - a + i/2)(~) = O.
(42)
We now express F as a series in the form F = ZZAjCI+ /
J
1.2.J+ 1,.2(~f)
(43)
Onset of convection in a porous mantle
7
where A/s are constants. F being chosen in this form, the solution for W can be written in the form W = 2 AjWj,
(44)
J
where ~$.~is a solution of 1 cJ'tl4"j = r-~2 Cl + 123+
t,2(2j r}
(45)
which satisfies the necessary b o u n d a r y conditions at r = l and r = q. The general solution of (45) is given by 1 9~rl,2 Ct + 1,2,l+ 1,'2('~Xjr)
~Vj --
+
B(1j)rl +
B~ u
1,
(46)
where B u~ and B~ ~ are constants of integration to be determined by the b o u n d a r y conditions. Substitution of solutions for F and W in accordance with (43), (44) and (46) into (38), we get
= l(! + I)RD~Aj ~ C l + j
(~jr
1.2.~§ ~,2(:~jr) +
+
.
(47)
We now multiply (47) by r 3'2 Ct+ ~:2.~§ l/2{:qr) and integrate the result over the range of r from q to 1 and obtain
where p~:.~ =
_ f l ( r2C'+ ~,2"z+'/z(aar) C, c'+ ''za+ l'2(~jr) ) j
(49)
Ct+ t.2.z+ 1..z(ad')(B']'/+ 3,z + By,r-t+ 1...2)dr.
(50)
and
QRj =
fn Equation (48) now leads to the secular equation
Owing to the Hermitian character of the o p e r a t o r r2~t with respect to the function which vanish at the end of the intervals, it follows that P~:~ = U~:~. The matrix element
(52)
Qkj can be evaluated explicitly with the use of recurrence formulae
8
G K Pradhan and B Patra
satisfied by cylinder functions. Utilising the relevant recurrence formula, we have 1
QRj = - - [ C l + 1;2,1 + 3/2(0~k) - -
?]1+ 3 / 2 C l +
1/2.1+
3i2(~kq)] B~J}
9~k
l __ _ _ [ C l + ;(k
1~2,/_ 1 2(~k)
_ ~] - I ~- t " 2 C l + 1 , 2 , 1 -
1 2 ( ~ X k l l ) ] B ( ~ ).
(53)
Further, the use of recurrence formulae satisfied by the function Ct+ 1/2.~.(s) reduces the relation (53) to the form
Qkj
1
~Gt+t..2.kBl + Ht+t'2.kB2 ;"
(54)
Further reduction of QKj requires an explicit consideration of the nature of b o u n d i n g surfaces at r = 1 and r = q. We now consider four possible cases separately for further reduction of the matrix element Qkj.
Case I. Impermeable surfaces at r = 1 and r = q. We apply the conditions (24) for impermeable boundaries at r = 1 and r = q to the solution (46) and obtain the values of B~tjJ, and B~j~. O n inserting the values of B]j~ and B~ ~ so obtained in (54) for Qk.j, we get Qkj = 0.
(55)
Case 2. Impermeable surface at r -- 1 and permeable surface at r = r/. Condition (24) for impermeable b o u n d a r y at r = 1 and permeable b o u n d a r y at r -- r/ are applied to solution (46) to obtain the values of B]jl and B~ ~. The values of B]jl and B~ } so obtained, are then inserted into (54) to give
Qkj=
C'~+ ~.2(~?I)
9tk~Xj~(l+l)ql+l/2 +lq._ll+3,,2)t; [Glet 2k--Ht+l/2,k]
(56)
which is symmetric in k and j.
Case 3. Permeable surface at r = 1 and impermeable surface at r = r/. In a similar manner, we insert the values of B~/} and B~ 1 obtained on applying the relevant conditions from (24) to the solution (46), into (54) and reduce the value Of Qkj to
Qkj =
C'l + 1,,2(~)
i
~k~/lq(2,+ 1) + (l + 1)} {G,+ L,2.k -- '1'2t+ ~)Ht+ 1..2.kS.
(57)
We observe that the value of Qkj is symmetric in k and j.
Case 4. Permeable surfaces at r = 1 and r = q. We apply the b o u n d a r y conditions (24) for permeable boundaries to solution (46) and then insert the values of B]jl and B~jI so obtained into relation (54) and obtain
QI,j -
{IGl+l/z'jGl+liZ'k +(l + l)rl~21+t}Hl+l/2"jHl+t/2'k} ~t~Otk{l(1+ I)(1 -- qqzt+ 1~)}
where Qjk is again observed to be symmetric in k and j.
-
(58)
Onset of convection in a porous mantle
9
From the symmetric form of Q~j in k and j, the underlying characteristic value problem is seen to be self-adjoint and consequently, solving the secular equation (51) is equivalent to Afs as variational parameters and minimizing the expression (32) for R o.
4. Geophysical models We now consider a few models of geophysical importance.
4.1. Model of the mantle when both temperature distribution and gravity distribution are linear The fluid is subject to a spherically symmetric radial gravitational field g(r)r, where
g(r) is a function of r only. Hence for a spherical shell of density p,, overlying a core of radius R i and mass Mi, y(r)=fgJ{(Ali -3grcpmR i31)r3 + 4ltpm}
(59)
where ff is the constant of gravitation. For such a model the distribution of heat sources in the shell is uniform and the gravity varies as r and as such
b(r) = c(r) = 1.
(60)
Now (49) for specific condition (60) gives
pnl = __ ~1 r e Ct+ l/2.1+ t,2(ulr) (j C,. 1,2 l+ ,,,2(ukr) k.j j~ r~,2 ~ rU 2 dr = a~Nt+ 1.2.i6kj.
(61)
The secular equation (51) with the substitution of the value of P~,'] from (61) reduces to
Nt+l/2,k{l[l +~
ot21}~Skj+Qk j =0.
(62)
And in the first approximation (62) gives Nt+ t,,2,1~4
I[1 + I)R o = N1+12.1 -- Qllo~
(63)
where cr is a first zero of {42) given in Chandrasekhar and Elbert [6]. The above formula (63) together with the expressions for Q~s obtained in cases 1 4 earlier, has been used to determine the characteristic number R o for various values of l and r/ and for the set of four boundary conditions when the gravity varies as the radius. The results are tabulated in tables 1-4 and have been further illustrated in figures 1-4.
10
G K Pradhan and B P a t r a Table 1. Characteristic number R o for various values o f / and q when gravity varies as r. Impermeable surfaces at r = 1 and r = q.
1 l 2 3 4 5 6 7 8 9 10 1l 12 13 14 15
q = 0'2 { • 102}
q = 0-3 ( • 102)
q = 0.4 [ • 102j
q = 0-5 ( • 102J
q = 0"6 [ x 102)
q = 0"8 I • lOS1
2-41172 1.88564" 1.99335 2.24224 2.55400
323320 2-10474 2.06019* 2-26159 2.55922 2.90953
5-05558 2.72258 2.34118 2.39504 2.62086 2.93674 3-30909
9-32745 4.26398 3.15645 2.88668* 2-92913 3.12992 3-42796 3.79359 4-21033 4-66804 516056 5-68376 6-23549 6.81377 7.41780
2.10522 + 8.46827 5-44639 4-37503 3.98466 3.91102* 4-01674 4.23927 4-54633 4.91917 5.34624 5.81948 6.33336 6.88355 7-46689
3.10561~ 1.07681~ 5.70407 3-68690 2.68745 2.12608 1.78507 1.56793 1.42631 1.33384 1.27513 1.24072 1.22450 1.22246" 1.23175
* Denotes m i n i m u m value of Ro; t • 103; ~ • 104.
Table 2. Characteristic n u m b e r R o for various values o f / a n d rt when gravity varies as r. Impermeable surface at r = 1 and permeable surface at r = q.
I 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
q = 0"3 ( x 102)
q = 0"4 ( x 102)
q = 0"5 ( x 102)
q = 0'6 ( x 102)
q = 0"8 ( x 102)
2.13882 1.86146" 1.99294 2.24288 2.55425 2.90827
2.69862 2.03550* 2.06542 2-27407 2.56723 2.91347 3-29931
4.17776 2.59112 2.36535* 2.45463 2.67700 2.97951 3-33800 3.74029 4.17925 4.66685 5.15058 5.67826 6-23251 6-81219 7.41697
8.29264 4.18020 3.29216 3.08800* 3-14146 3.32867 3.60194 3.93869 4.32670 4.75834 5-22871 5.73403 6.27168 6.83938 7.43563
1.29763 t 3.924155+ 2.24362§ 1-58127~ 1.26053~ 1.08814~ 9-91470 9.38314 9.12473 9-05127 9-10774 9.27153 9.51255 9-82050 1.018545
* Denotes m i n i m u m value of Ro; t x 104; 3~ • 103.
O n s e t o f c o n v e c t i o n in a p o r o u s m a n t l e
11
Table 3. Characteristic number R n for various values of l and q when gravity varies as r. Permeable surface at r = I and impermeable surface at r = q,
l 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
q =0-2 ~• 10:}
q =0.3 (• 1021
q=0.4 (•
q =0-5 (• )
q =0.6 (x 10 z)
q =0-8 [• 102t
1-15110 1-11932" 1.32670 1-60126 1-91528
1.45243 1.20183" 1.35084 1.60777 1.91708 2.26231
2.12036 1.44640" 1.49245 1.66049 1.94035 2.27204 2.64179
3.65983 2.05567 1-80867" 1.85971 2.07514 2.35524 2.69162 3.07201 3-48913 3.93831 4.41663 4.92212 5.45385 6-01095 6.59244
7.79164 3.66510 2.75820 2-53222* 2.56233 2.72536 2.97413 3.28629 3.64973 4-05692 4.50298 4-98414 5.49776 6.04160 6.61396
1.06159t 3.85795~ 2.17674~ 1-51343~ 1.19156~ 1-01787~ 9-19768 8.65060 8.37564 8-28483* 8.32411 8.46849 8.69044 8.97900 932415
* Denotes minimum value ofRo; t • i0.~; ++ • 103.
Table 4. Characte.ristic number R o for various values o f / a n d r/when gravity varies as r. Permeable surfaces at r = 1 and r = ~l where 0"5=q~
1 1 2 3 4 5 6 7 8 9 10 II 12 13 14 15
q = 0"5 (• 102 )
rl = 0-6 (• )
q = 0-8 (• z)
1-04286" 1.17788 1-37402 1"60119 1.92173 2-25938 2"63192 3'03539 3.46714 3"92540 4-40927 4-91796 5.45156 6.00971 6-59179
1.32370" 1-43330 1"59519 1'80663 2"06454 2"36549 2.70653 3"08466 3.49720 3-94181 4-41662 4"91992 5"45048 6.00717 6-58932
3.84037* 3.91109 4.01680 4.15736 4-33227 4-54135 4"78374 5"04196 5"36741 5.70752 6.07600 6.48145 6.91417 7.37713 7.86950
* Denotes minimum value of R n.
12
G K Pradhan and B Patra
~D
II
Od
0
E~ ~
E
r o
~
od
61
o
___I
o
io
6
ZN 0
C]&l O~6o I
If
o
=~-
$. c~
d
I
I
0
~ a~l o~601
o
I
0
o
~
0
~.~
13
O n s e t o f c o n v e c t i o n in a p o r o u s m a n t l e
4.0
8
3.5
3.0 -o.6
? 0\ . 5 ~j -0.4 2.5-
0.3 0.2 2.0
1
I
I
I
I
I
I
I
0
2
4
6
8
10
12
14
16
L Figure 3. Neutral stability plots for the mantle at r = 1 a n d i m p e r m e a b l e
boundary
Ibir)
= I = ~_'(r))
with permeable boundary
a t r = q.
3-0
a~~ 2,5 o
2-0
0
I
I
I
l
I
I
I
2
4
6
8
10
12
14
16
L
Figure 4. Neutral stabihty plots for the mantle (bir) = l = c(r)) with permeable boundaries both at r = l andr=q(0'5~
14
G K Pradhan and B Patra
4.2. Model of the mantle with linear temperature and non-linear gravity distribution As suggested by Chandrasekhar [5] and Lyttkens (refer p. 243 of [5]), in the mantle of the earth, the value of gravity remains approximately constant. Here we consider the case of a model with uniform distribution of heat sources and constant gravity through the mantle, so that we have
b(r) = l
and
c(r) = r-1.
(64)
For such a specific case, (49) assumes the form pl,~(,) k,j =
~~f
i (rCl+
1 / 2 , / + l / 2 ( ' g t k r ) ' r C t + 1/2.1 + 1 , 2 ( ~ j r ) }
dr
= ~tk2(k [rlj).
(65)
The characteristic equation (51) with the aid of (65) becomes
~t2(klrlj) - l(l+ l)Ro[
"~6ki--Qkj]=O
(66)
and in the first approximation, (66) gives (llrll)~
l(l + l)Ro - Nl + 1.,2.l -- Q1 la~
(67)
This relation (67) together with the expressions for Q~j obtained for the four sets of boundary conditions has been used to determine Ro for various values of r and q, when the gravity is constant through the mantle. The values so obtained have been tabulated in tables 5-8 and have further been illustrated in figures 5-8. Table 5. Characteristic n u m b e r R o for various values of I and )7 when gravity is constant through the mantle. I m p e r m e a b l e surfaces at r = 1 and r = q.
I
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
t / = 0-2 ( X 10 2)
q = 0-4 ( • I0 2}
q = 05 ( • 10 2)
t / = 0.6 ( • 10 2)
q = 0.8 ( X 10 3}
1.49768 1.23321" 1-36943 1"60320 1.88485
3'55913 1-93575 t'68799" 1"75557 1.95778 2-23192 2.55793
7-00491 3"21504 2-39259 2'20254* 2-25543 2"43195 2-68752 3.00452 3"36826 3-77177 4,21102 4-67774 5-16922 5"68950 6.23837
1.68417t 6-78309 4.37345 3"52190 3"21562 3.16401* 3"26160 3.45924 3-72345 4"04848 4.42134 4.83599 5.29469 5-78218 6"30205
2-795055 9"69131 5-13367 3-31821 2.41870 1.91347 1-60835 1.41271 1"28510 1"20179 l'14089 1.11913 1-10450" 1.10266" 1-11104
* Denotes m i n i m u m value of Ro: ~" • 103: + • 10 "~.
O n s e t o f c o n v e c t i o n in a p o r o u s m a n t l e
15
Table 6. Characteristic n u m b e r R o for various values of I and r/when gravity is constant through the mantle. Impermeable surface at r = l and permeable surface at
r=q.
l
5 6 7 8 9 10 11 12 13 14 15
q = 0.2 I x 1021
q = 0.4 [ • 102~
q = 0.5 ( x 102)
~ = 0-6 I • 102}
q = 0.8 ( x 10zJ
1"23606 1-21219" 1'36488 1'60261 1'88504
1-89983 1"44724" 1"48917 1"66690 1"91772 2'21424 2'55037
3.13750 1.95371 1.79294" 1.87288 2.06129 2.31508 2.61699 2-96231 3.34340 , 3.77081 4.20287 4.67321 5.16675 5.68818 6-23767
6.63411 3.34834 2-64360 1.48580" 2.53516 2.69289 2-92477 3.21397 3.54357 3.91611 4-32415 4-76498 5.24312 5.74508 6.27567
1"16786 t 3-531745 2-01925~ 1.423145 1.134475 9.79329 8.93314 8.45421 8.22138 8.15519* 8.20607 8.36292 8.58032 8.85809 9-18722
* Denotes m i n i m u m value of Ro', t • 104; + x 103.
Table 7. Characteristic n u m b e r R o for various values ofl and q when gravity is constant t h r o u g h the mantle. Permeable surface at r = 1 and impermeable surface at r = r/.
l
5 6 7 8 9 10 11 12 13 14 15
q = 0"2 ( • 10 2 )
q = 0'4 ( • 10 z)
q = 0"5 ( x 10 2 )
q = 0"6 ( x 10 2 )
q = 0"8 ( • 10 2)
0.71483* 0.73204 0-91144 1"14490 1"41348
1.49273 1.02839" 1"07606 1.21714 1.44944 1-72675 2-04211
2.74853 1.54997 1.37097" 1.41896 1-59786 1-83002 2"11023 2.43303 2"79131 3"18215 3-60397 4-05091 4.52124 5.01914 5-54424
6-23331 2-93574 1.21483 2.03844* 2"06780 2.20482 2.41499 2-68161 2-98913 3.33885 3.72396 4.14182 4"59612 5.07494 5.58218
955427 t 347216 t 1"95907 t 1"36209' 1"07240 t 9"16086 8"28711 7"79419 7"54645 7"46463 7-50000 7"63858 7"83878 8"09906 8.41038
* Denotes m i n i m u m value of Ro: + x 10 3.
16
G K Pradhan and B Patra Table 8. Characteristic n u m b e r R o for various values of I a n d r/ when g r a v i t y is constant t h r o u g h the mantle. Permeable surfaces at r = 1 and r = r/.
l 1 2 3 4 5 6 7 8 9 10 1l 12 13 14 15
q = 0"2 ( x 102 )
q = 0'4 ( • 102')
~ = 0"5 ( x 102)
q = 0"6 ( x 102)
q = 0"8 (x 102)
0.58915* 0.71570 0-90904 1.44457 1.41345
0.65266* 0-77255 0.96335 1.16564 1.42499 1.71547 2-03710
0.78319* 0.88812 1.04150 1-22171 1.47973 1.75554 2-06342 2.40403 2.77371 3.17172 3.59797 4.04748 4-51934 5.01811 5-54369
1-05896" 1.14807 1-28094 1.45434 1.66608 1.91368 2.19770 2.51708 2.86421 3.24411 3-65255 4.08845 4-55660 5.04602 5.56139
3.45633* 3-51998 3.61512 3.74162 3.89904 4-08721 4-31015 4.54281 4.83603 5.14248 5.47445 584627 623658 6-65417 7.09829
* Denotes m i n i m u m value of R o. 4.5 8
4.0
o3. 5 1 r~" o 3.0
25i.
2.0
0
I
I
I
I
I
I
r
2
4
6
8
10
12
14
16 L Figure 5. Neutral stability plots for the mantle (h(rl = 1, c(r)= r -t) with impermeable b o u n d a r i e s both at r = I and r = tl.
17
Onset of convection in a porous mantle
'~. II
__
C',.I
_
0
tt
"C r .LE ,.~ E
~E ,o =
I
-
kO
b~
J
_6
i_
0
f
~
d
I
o
6
d
~
I
6
o
9
i1) 9
o
E
._~
a B Otflo I
i,,.. II II "" "C~ ',-Cb C',,J
tt
0 E
--00
_f I o
o
~b
j.h
I
I
~
o
d~
i0
O0
o o
aH ~
.~.E_
18
G K Pradhan and B Patra :3-0
2.5 cr o tK o
2~ 1"6
I 2
0
J 4
I 6
I 8
I 10
I 12
I 14
16
L Figure g. Neutral stability plots for the mantle (b(r)= I, c(r)= r -~) with permeable boundaries both at r = 1 and r = q.
Table 9. Characteristic number R o for various values of l at r/= 0 with impermeable and permeable bounding surface at r = 1 when gravity varies as r.
1
Impermeable ( x 102)
Permeable ( • 102)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
2.03832 1.83900 1.98707 2.24147 2-55390 2.90824 3.29732 3-71753 416672 4.64348 5.14693 5.67646 6.23149 6-81168 7-41669
0.98090 1-06266 1.27857 1-54919 1.85801 2.19920 2-56998 2.96887 3.39487 3.84727 4.32550 4-85730 5-38786 5.94338 6-52328
Values of R o for 1 = I to 11 have been obtained in the third approximation and the rest in second approximation.
Onset o f convection in a porous mantle
19
Table 10. Characteristic number R o for various values of 1 at q = 0 with impermeable and permeable bounding surfaces at r = 1 when gravity is constant through the mantle. Impermeable
Permeable
I
( • 10 2 )
( • 10 2)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
1.20261 1.19167 1.36313 1.60265 1.88478 2.19863 2.54553 2.91826 3-31671 3-74265 4.18960 4.66605 5.16591 5.68775 6.23743
0.57873 0-68860 0,87710 1.10767 1,37121 1,66259 1.98402 2.33056 2.70232 3.10090 3.52096 3.99270 4.46654 4-96272 5-48608
F u r t h e r , the values of characteristic n u m b e r R n for t / = 0 with either i m p e r m e a b l e or p e r m e a b l e b o u n d a r y for the first fifteen modes, when the gravity varies as r, o b t a i n e d by P r a d h a n a n d P a t r a [12] a n d when gravity is c o n s t a n t t h r o u g h the mantle have been presented in tables 9 a n d 10 respectively a n d p l o t t e d in figures 1 - 4 a n d 5-8. 4.3. Model o f thermal convection in the mantle with a massive core W e now c o n s i d e r a m o d e l of the m a n t l e overlying a massive core so that the c o n t r i b u t i o n of the m a n t l e to the gravity can be totally i g n o r e d a n d there is uniform d i s t r i b u t i o n of heat sources in the system. In the geophysical context, to a fair first a p p r o x i m a t i o n , the m a n t l e can be c o n s i d e r e d as a spherical shell of fluid-saturated p o r o u s m a t t e r overlying a massive solid core, so that c o n t r i b u t i o n of the mantle to the gravity can be c o n s i d e r e d negligible c o m p a r e d to that of the core. F o r this model, the radial g r a v i t a t i o n a l field for a core of mass Mi is given by M i _- (constant)r -2 ro(r)= ff -~-
(68)
E q u a t i o n (68) with the a i d of relation ~g(r) = ?(r) = ?lc(r)
(69)
c(r) = r - 3,
(70)
gives
where we have a s s u m e d that )'1 = ~ffMi. Since there is uniform d i s t r i b u t i o n of heat sources, for the m o d e l u n d e r discussion, we have b(r) = 1
and
c(r) = r -3
(71)
20
G K Pradhan and B Patra
The matrix element given in (49) with the aid of relation (71) reduces to pl,cr k,j
=~2
= ct2
1 rCl+t/2.t+l/2(~kr)~(~Ct+l/2,~+l/2(~qr)
f l r4{Ct+
1/2,l+
1,'2(~
1/2,1+ 1/2(~
}
dr
(72)
dr
Denoting p kl,.r~ _-- ~2(klc(r )- 11 J) = ct2(klr31j) j
(73)
the secular determinant (51) is obtained as ~t2(klr3lj) - l(l+ e ) R o [ ~ f k j - - Q k j ]
(74)
=0
and in the first approximation, we have (llr3l 1)~ I(l + I)R o - Nt + 1/2,1 - Qll ~x2"
(75)
The values of the characteristic number R o have been calculated from relation (75) for various values of l, for q = 0"5 and for four sets of boundary conditions considered. They are tabulated in table 11 and the corresponding neutral stability plots have been illustrated in figure 9. Table 11. Characteristic n u m b e r Ro for r / = 0"5 and for different laws of variation of gravity and temperature distributions (either b(r) = 1, c(r) = r - a or b(r)= r -3, c(r)= 1). Perm. at r = l r= q ( X 10 2) l 2 3 4 5 6 7 8 9 10 11 12 13 14 15
4"6144589 *t 5"2656253 t 6"2416493 t 7'4165087 t 9'1089848 t 1"0986330 1'3162442 1"5605619 1"8331306 2"1308233 2"4552142 2"8043270 3'1778860 3'5751278 3'9925100
Perm. at r = 1, Imp. at r = rt ( X 10 2)
Imp. at r = 1, Perm. at r = r/ ( X 10')
Imp. at r = 1, r = r/ ( X 10 2)
1"6201914 9.1897206 t 8.2160987 t 8.6139405 t 9"8361469 1"1452453 1"3461006 1'5793891 1'8447571 2'1378312 2'4593125 2'8066991 3"1792209 3"5758654 3"9929037
1"8494667 1'1583410 1"0744883" 1'1369534 1"2688958 1"4487992 1"6693605 1'9229668 2"2096342 2"5333043 2"8679979 3'2378665 3"6331264 4"0525166 4'4923044
4"1292005 1'9061807 1"4338549 1"3370735" 1"3884052 1"5219367 1"7143502 1'9503695 2"2260667 2"5339502 2'8735510 3"2410027 3"6348630 4"0534560 4'4928070
* Denotes the m i n i m u m value of R o. Perm. - Permeable: Imp. - Impermeable. txlO.
Onset of convection in a porous mantle 5.0 (i)
21 t b(r)=l, C(r] = " ~
(ii} b(r)= -~3 ,C(r) = 1
,\
o
s o 2.0
= 9 " ~ y
(a) (b) (c) (d)
(0) ~
o/1.5
I 2
0
I 4
I 6
Perm. at r=l, Perm. at r=1, ~rnp. at r=l, Imp. o t r = l ,
I 8 L
I 10
Perm.at r=~/ lrnp. atr='~ Perm. otr r~ Zmp. a t r = ~
I 12
I 14
16
Figure 9. Oarcy-Rayleigh number Ro for the onset of convection in the mantle as a spherical harmonic disturbance of order I for the shell thickness q = 0.5 for various types of boundaries when b(r)= |, c ( r ) = r -a or b ( r ) = r -3, r 1.
4.4. Model of the mantle when the gravity distribution is linear and all
the heat sources are confined to the core For such a geophysical model where all the heat sources are confined to the core and none exists in the mantle and the variation of gravity is same as in a homogeneous sphere, from the expressions for fl(r), constancy of fll and linear gravity distribution, we obtain
b(r) = r -3
and
c(r)=l.
(76)
By (52), the matrix element given in (49) can be written as pb(r), 1 k,j
pt,o(,)=~2(jlralk )
= --j,k
(77)
and as such the secular determinant (51) now becomes
92 ( j l r 3 l k ) - l ( l + 1)RoINt~k/Z'k 5kj--Qkj] = O.
(78)
Since Qkj is symmetric in k and j, the secular determinant (78) by transposing rows and columns is the same as the determinant (74) with c replaced by b. Hence Ro{for c(r)=f(r) anU b(r)= 1} = Ro{for c(r)= 1 and b(r)= f(r)}.
(79)
The values of the characteristic number R o for this model is therefore the same as for the model 4.3 presented in table 11 and illustrated in figure 9.
22
G K Pradhan and B Patra
4.5. Model of the mantle with constant density and with uniform distribution of heat sources As shown in (59) a constant value of gravity is strictly inconsistent with the model of the mantle of constant density overlying a core. A variation of gravity compatible with (59) is now given by rg(r)=constant
(6,) ~r+7~-
, (when 89
(80)
The distribution of heat sources being uniform, we also have
b(r) =
1.
(81)
F r o m (80) and (81) we obtain for the model under consideration
c(r)=
7
and
b(r)= 1.
(82)
By (82), the value of the matrix element (49) is given by
pl.~l~ =
rCl+ t,.2,1+t/2(oq,r)~) Cl+ l,,z,t+ l/z(~ir ) dr
=~t~fl{rCt+,.,z.t+l/z(.~kr)(6+~r-3)-~Ct+l.z.,+~/z(~jr)}dr.
Table 12. Characteristic n u m b e r R o for ~/= 0,5 and for laws of variation of gravity and temperature distribution: c(r) = (6 + r-3),'7 and b(r) = 1.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
Perm. at r = 1, r=q ( x 10 2 )
Perm. at r = 1, Imp. at r = t / (• 10 2 )
Imp. at r = 1, Perm. at r = t / ( • 10 2 )
Imp. at r = 1, r=q ( x 10 2 )
8.6091764 * t 9.7589595 t 1.1440644 1-3407913 16212123 1.9209903 2-2560936 2.6240633 30226790 3.4474762 3-9034517 4-3816692 4-8864267 5.4170035 59706851
3.0227843 1.7031616 1.5059715" 1.5572687 1.7506323 2-0024933 2.3072687 2,6557210 3.0418500 3.4588144 3.9099674 4.3853756 4.8884793 5.4181212 5,9712739
3-4505423 2-1467921 1-9694857" 2-0554379 22583742 2-5332658 28613486 32334421 3-6435019 4-0986535 4-5597209 5-0590605 5-5864199 6-1403391 6-7181134
77038319 3~5327884 26281874* 2.4172244* 2.4710764 2-6611487 2.9384627 3'2795194 3-6705976 4.0996986 4.5685560 5-0639607 55890910 61417632 6.7188652
* Denotes the m i n i m u m ,~alue of R o. Perm. -- Permeable: Imp. - Impermeable. 1 x 104.
(83)
23
O n s e t o f c o n v e c t i o n in a p o r o u s m a n t l e
3.0
2.5 g
(
..9.0
{b} Perm. at r - l , amp. at rffi'q
(O}
I t ) Imp. at r=l, Perrn.atr=~ (d) Imp. or r=l,'rmp, ot r=~ b(r)=l,C(r}=
1.5
I 2
0
I 4
I 6
[ 8
I 10
(e+r-3) 7
I 12
[ 14
16
I,
Figure 10. Darcy-Rayleigh number R o for the onset of convection in the mantle as a spherical harmonic disturbance of order I for the shell thickness r/= 0'5 for various types of boundaries when b(r) = 1, c(r) = (6 + r-3)/7. Denoting pl,c(r) kj = ~(klc(r)-
(84)
l lj)
the secular determinant (51) reduces to the form
(85) where
(klc(r)-'lJ)=f.z{(6r~Tf)C,+.2.t+,,'(=~r)C,+~/2.,+s/z(~ l
7r 4
(86) In the first approximation, we obtain
(11c(r)-x I1)~ l(l + I ) R o = Nt+1/2,1 - Q I I a ~
(87)
which gives the values of R o for r / = 0.5, for various values of I and for the set of four b o u n d a r y conditions. These values are presented in table 12 and further illustrated in figure 10.
5. Conclusion It is shown [figures 1-8] that in the first two models (when gravity varies as r and gravity remains constant through the mantle), the easiest mode for the onset of convection in the
24
G K Pradhan and B Patra
mantle gradually recedes to a spherical harmonic of higher order as the mantle thickness decreases. This is true in all the three cases of the shell having both impermeable, impermeable-permeable and permeable-impermeable bounding surfaces. But, in contrast, in the permeable-permeable case the onset of thermal convection occurs at spherical harmonic of order 1, irrespective of the thickness of the mantle. The first three cases are in harmony with the corresponding (ridig and/or free bounding) cases of the continuous fluid, where l increases with the decreasing mantle thickness, as established by Chandrasekhar [5]. But there is a difference in case both the bounding surfaces are permeable, as the lowest harmonic occurs at l = 1, in contrast to that of homogeneous fluid with free-free boundaries. The result for the third and fourth models [cf. w and w is the same because of the identity of the underlying characteristic value problem. From table 11 and figure 9 it is seen that for a particular value of the thickness of the shell, r/= 0.5, the order of the spherical harmonic shifts from 1 when both the surfaces are permeable, to 3 when one surface is permeable and the other impermeable, and finally to 4 when both surfaces are impermeable. For the fifth model of the mantle with constant density and with uniform distribution of heat sources for radius ratio r/= 0.5, the dependence of critical Darcy-Rayleigh number for the onset of convection on the order of the spherical harmonic is similar to models III and IV. The values of Ro for all the five models discussed here are comparatively lower than those for the continuous fluid shell and merit comparison with available data. It has been argued that the mantle with fluid magma and semisolid rocks behave similarly as a porous medium saturated with a fluid. The main difference between the continuous fluid model discussed in Chandrasekhar [3] and the porous fluid model discussed in this paper seem to be the sensitivity to kinematic boundary conditions that the porous model exhibits. This brings into sharp focus a number of qualitative differences between the two models. For example, consider the manner of the onset of thermal instability in a porous mantle. In the case of a porous mantle occupying half the radius of the sphere, we have seen that if both the boundaries are permeable, the current system must belong principally to the mode l = 1, instead to the modes l = 3, 4 and 5 as in a homogeneous fluid shell or to the modes l = 3, 4 and 6 as in a porous fluid shell with impermeable boundaries or to the modes l = 2, 3 and 4 as in a porous fluid shell with one of the boundaries permeable and the other impermeable. It is further observed from the figures that in all the five models discussed in this paper the neutral stability plots for different mantle thickness do not cross but lie symmetrically one above the other. The critical values of the Rayleigh number for the models discussed here are comparatively lower than those for the continuous fluid shell and merit comparison with available data.
Acknowledgement This work was supported by the University Grants Commission, New Delhi, under Faculty Improvement Fellowship awarded to the junior author (B P). This paper constitutes a part of the PhD Thesis in Applied Mathematics submitted to Sambalpur University.
O n s e t o f c o n v e c t i o n in a p o r o u s m a n t l e
25
References [I] Abramowitz M and Stegun I A, Handbook t~['mathematicalfimctions (New York: Dover Pub.) Ch. 10 (1965) [2] Busse F, Pattern of convection in spherical shells J. Fluid Mech. 22 [ 1975) 67-85 [3] Chandrasekhar S, The onset of convection by thermal instability in spherical shells Philos. Ma.q. Ser. 7 44 (1953) 233 241 [4] Chandrasekhar S, A correction Philos. May. Set. 7 44 (1953) 1129-1130 [5] Chandrasekhar S, Hydrodynamic and hydromagnetic stability IOxford: Clarendon Press) Ch. 6 (1961) [6] Chandrasekhar S and Elbert D. The roots of J-it. ~ ~ 2~(,~/lJ1 ~ i 2(2) - J l . ~~2(2r/)J -~1 + t,2~0-) = 0 Proc. Camb. Philos. Soc. 49 (1953J 446-448 [7] Combarnous M A and Bories S A, Modelisation de La Convection naturelle au Sein d'une Couche poreuse horizontale a" Solide-fluide Int. J. Heat Mass TransJbr 17(1974) 505-515 [8] Combarnous M A and Bories S A, Hydrothermal convection in saturated porous media Advances in hydroscience 10 {1975) 231-307 [9] Heiskanen W A and Vening Meinesz F A, The earth and its gravity field {New York: McGraw Hill) Ch. [1 (1958) [10] Joseph D D, Stability t~/[luid motion I and II [Berlin: Springer Verlag) (1976) [ l l ] Lapwood E R, Convection of a fluid in a porous medium Proc. Camb. Philos. Soc. 44 (1948) 508-521 [12] Pradhan G K and Patra B, Porous convection in a globe with internal heat sources Indian J. Pure Appl. Math. 14 (1983) 661 670 [13] Relton F E, Applied Besselfimctions ILondon: Blackie and Sons Ltd) Ch. 3, 4 [1949) [14] Simpkins P G and Blythe P A, Convection in a porous layer Int. J. Heat Mass Transfer 23 (1980) 881-887
