Transp Porous Med (2010) 85:397–414 DOI 10.1007/s11242-010-9569-9

Natural Convection Boundary-Layer Flow in a Porous Medium with Temperature-Dependent Boundary Conditions J. H. Merkin · I. Pop

Received: 5 November 2009 / Accepted: 12 March 2010 / Published online: 31 March 2010 © Springer Science+Business Media B.V. 2010

Abstract The natural convection boundary-layer flow on a surface embedded in a fluidsaturated porous medium is discussed in the case when the wall heat flux is related to the wall temperature through a power-law variation. The flow within the porous medium is assumed to be described by Darcy’s law and the Boussinesq approximation is assumed for the density variations. Two cases are discussed, (i) stagnation-point flow and (ii) flow along a vertical surface. The possible steady states are considered first with the governing partial equations reduced to ordinary differential equations by similarity transformations and these latter equations further transformed to previously studied free-convection problems. This identifies values of the exponent N in the power-law wall temperature variation, N = 3/2 for stagnation-point flows and 3/2 ≤ N ≤ 3 for the vertical surface, where similarity solutions do not exist. Time development for stagnation-point flows is seen to depend on N , for N < 3/2 the steady state is approached at large times, for N ≥ 3/2 a singularity develops at finite time leading to thermal runaway. Numerical solutions for the vertical surface, where the temperature-dependent boundary condition becomes more significant as the solution develops, show that, for N < 3/2, the corresponding similarity solution is approached, whereas for N > 3/2 the solution breaks down at a finite distance along the surface. Keywords Convective flow · Porous media · Boundary-layer flow · Temperature-dependent wall conditions

J. H. Merkin (B) Department of Applied Mathematics, University of Leeds, Leeds, LS2 9JT, UK e-mail: [email protected] I. Pop Faculty of Mathematics, University of Cluj, 3400 Cluj, Romania e-mail: [email protected]

123

398

J. H. Merkin, I. Pop

1 Introduction Boundary-layer flows in fluid-saturated porous media driven by convection have already received considerable attention, see Ingham and Pop (2005); Vafai (2005); Nield and Bejan (2006) for examples. A wide range of geometrical configurations and heating conditions has been treated for both natural convection as well as cases where there is an interaction between the buoyancy-induced convection and some external imposed flow (mixed convection). The majority of these previous studies have assumed either some prescribed wall temperature or some prescribed wall heat flux distribution. An alternative way in which a convective flow can arise is if the wall heat flux is dependent on the local wall temperature. This situation can arise, for example, through Newtonian heating, where the wall heat flux is proportional to the wall temperature Lesnic et al. (1992), or when there is an exothermic catalytic reaction on the wall. This latter case has been examined by Mahmood and Merkin (1998); Merkin and Mahmood (1998) taking an Arrhenius temperature dependence for the reaction rate. A different way of modelling a surface reaction and which also gives a relation between the wall heat flux and the wall temperature is through some power-law variation in the temperature. This is the case we treat here. We take k

N  ∂T = −γ T − T∞ ∂y

on y = 0

(1)

where T is the wall temperature and T∞ is the (constant) ambient temperature, with k being the thermal conductivity, y measuring distance normal to the wall and γ is a constant. A related form of local heating, now assumed to take place within the boundary-layer flow, has been discussed by Magyari et al. (2007); Mealey and Merkin (2008); Merkin (2008, 2009), where it was seen that the boundary-layer development depended critically on the exponent N . The case of boundary-layer flow in a viscous fluid with a wall heat flux—temperature relation given by (1) has been treated by Merkin (1994) for Newtonian heating, effectively taking N = 1 in (1) and by Merkin and Chaudhary (1996) for a general temperature dependence of the wall heat flux. Here we start by considering steady natural convection flows set up on a bounding surface on which the boundary condition (1) is applied. We treat the separate cases of flow near a stagnation point and flow on a vertical surface. We find that we are able to transform these problems into standard free convection problems that have been treated previously. This identifies values of the exponent N for which solutions do not exist. We then consider the time evolution of the stagnation-point flow starting with the fluid temperature everywhere at the ambient temperature. We find that the development of this flow depends critically on the exponent N in (1); for N < 3/2 the solution approaches the corresponding steady state at large times, whereas for N > 3/2 a singularity develops at a finite time leading to thermal runaway. For the vertical surface case we consider a modified version of boundary condition (1) through which the flow starts from a prescribed wall heat flux condition with the wall temperature dependence in (1) coming increasingly into play as the flow develops from the leading edge. We find that the similarity solution corresponding to (1) is attained at large distances along the wall if N < 3/2. However, for N > 3/2 a singularity arises in the numerical solutions to this modified problem at a finite distance along the wall, the nature of which is discussed, as is the critical case N = 3/2. We start by describing our model.

123

Natural Convection Boundary-Layer Flow

399

2 Model The basic boundary-layer equations for our model, which is based on Darcy’s law, are given in Ingham and Pop (2005); Vafai (2005); Nield and Bejan (2006), for example. For our specific model we make these equations dimensionless following Mahmood and Merkin (1998); Merkin and Mahmood (1998) by introducing a reference temperature Tref and a velocity scale U0 defined by  Tref =

k γ



gβ K ναm

1/2 2/2N −3 ,

U0 =

gβ K Tref ν

 N =

3 2

 (2)

gβ K Tref and where K is the permeability αm ν of the porous media, ν the kinematic viscosity, β the coefficient of thermal expansion, g the acceleration due to gravity, αm the effective thermal diffusivity and  a length scale. Expressions 2 already indicate that having N = 3/2 is a special case that will require further consideration. We make the equations dimensionless by putting with a corresponding Rayleigh number Ra =

u = U0 u, v = U0 Ra −1/2 v, x =

x , 

y=

y Ra −1/2 , T − T∞ = Tref T 

(3)

where (u, v) are the velocity components and x measures distance along the wall. This leads to the steady dimensionless version of our model as ∂u ∂v ∂x + ∂y = 0 u = S(x)T

(4) (5)

u ∂∂Tx + v ∂∂Ty = ∂ T2 ∂y 2

(6)

where S(x) = x for stagnation point flows and S(x) = 1 for flow on a vertical surface. Equations 3, 4, 5 are subject to the boundary conditions v = 0,

∂T = −T N on y = 0, ∂y

T → 0 as y → ∞

(7)

For the time-dependent stagnation-point flow problem we include a time derivative in Eq. 6 which now becomes ∂T ∂T ∂T ∂2T +u +v = ∂t ∂x ∂y ∂ y2

(8)

together with the initial condition T = e−y

at t = 0,

(0 ≤ y < ∞)

(9)

designed to satisfy boundary conditions (7) and to input heat initially into the system and thus generate a nontrivial solution. In (8) we have used σ /U0 as a time scale. We start by considering stagnation point flows and then treat the flow on a vertical surface, both for different values of the exponent N .

123

400

J. H. Merkin, I. Pop

3 Stagnation Point Flow 3.1 Steady Flow Here S(x) = x and we look for a solution of the steady problem in the form ψ(x, y) = x f (y),

T = T (y)

(10)

where ψ is the stream function, defined in the usual way. Equation 5 gives 7) then giving f  + f f  = 0,

f  (0) = − f  (0) N ,

f (0) = 0,

f

= T , with (6,

f  → 0 as y → ∞

(11)

where primes denote differentiation with respect to y. We can construct a solution to the problem given by (11) following the method used previously by Merkin and Mahmood (1998); Merkin and Chaudhary (1996); Chaudhary et al. (1995). Suppose that f  (0) = θ0 , where θ0 is to be found. We then make the transformation N /3

f = θ0

f˜,

N /3

y˜ = θ0

y

(12)

This leads to the problem f˜  + f˜ f˜  = 0,

f˜(0) = 0,

f˜  (0) = −1,

f˜  → 0 as y˜ → ∞

(13)

where primes now denote differentiation with respect to y˜ . The problem given by (13) has a solution for which f˜  (0) = a0 , with a0 = 1.36472. To solve our original problem (11) we note that   3 2N /3 ˜  2N /3 3/(3−2N ) θ0 = f  (0) = θ0 a0 , giving θ0 = a0 (14) N = f (0) = θ0 2 From (14) we can see that the wall temperature θ0 → 1 as N → ∞ (from below) and as N → −∞ (from above) and that for N = 1, θ0 = 2.5417. A graph of θ0 against N is shown in Fig. 1. This graph shows that θ0 → 1 as N → ±∞ and that the solution is singular at N = 3/2. If we put N=

3 + δ, 2

|δ|  1

(15)

we have, from (14), that θ0 ∼ e−3 log a0 /2δ

as δ → 0

(16)

Expression (16) shows that θ0 is exponentially large as N → 3/2 from below (δ < 0) and is exponentially small as N → 3/2 from above (δ > 0), as can be seen in Fig. 1. We also note that, when N = 3/2, problem (11) is invariant to the transformation f = α f˜,

y˜ = αy

for any α > 0

(17)

3.2 Time Development We now write ψ(x, y, t) = x f (y, t), T = T (y, t), following from (10), in Eq. 8 and solve the resulting equation, subject to boundary conditions (7) and initial condition (9), numerically using an approach based on the Crank-Nicolson method with the resulting set of nonlinear finite-difference equations being solved using Newton-Raphson iteration. This technique has

123

Natural Convection Boundary-Layer Flow

10

401

0

8

6

4

2

0 -3

-2

-1

0

1

2

3

N 5

4

Fig. 1 The wall temperature θ0 plotted against N obtained from (14) for the steady stagnation point flow

1.0

2.6 0

2.4 2.2 2.0 1.8 1.6

0

1.4

-1.0 1.2

t

1.0 0

1

2

3

4

5

6

Fig. 2 The wall temperature θ0 plotted against t for N = 1, 0, −1 obtained from the numerical solution of Eq. 8 for the unsteady stagnation point flow, showing that the solution approaches its respective steady state for t large

been applied previously in several related problems, see Mahmood and Merkin (1998); Mealey and Merkin (2008); Merkin and Pop (2000); Chaudhary and Merkin (1996) for example. In Fig. 2 we plot the values of θ0 = T (0, t) against t for representative values of N < 3/2. In each case we see that the numerical solution approaches the steady state solution given by (14). The same behaviour was seen for all the other cases of N < 3/2 tried and we can 3/(3−2N ) conclude that, for N < 3/2, the steady state θ0 = a0 is approached at large times. For N > 3/2 we found that the numerical solution broke down at a finite time ts , the value of ts depending on N , with θ0 becoming large as t → ts . We illustrate this case in Fig. 3 with

123

402

J. H. Merkin, I. Pop

0

14

2.0

3.0 12 10 8

5.0 6 4 2

t 0.0

0.05

0.1

0.15

0.2

0.25

0.3

0.35

0.4

0.45

Fig. 3 The wall temperature θ0 plotted against t for N = 2, 3, 5 obtained from the numerical solution of Eq. 8 for the unsteady stagnation point flow, showing the development of a singularity at finite time

plots of θ0 against t for N = 2, 3, 5. In each case the values of θ0 become large close to ts , with the breakdown in solution occurring at earlier times as N is increased. To discuss the nature of this singularity in the solution at ts we put τ = ts − t and look for a solution of Eq. 8 valid for τ > 0 and τ  1. The singularity arises in an inner layer close to the wall where we put f = τ (N −2)/2(N −1) φ(ζ, τ ),

ζ =

y τ 1/2

(18)

with Eq. 8 then becoming ∂ 3φ ∂ζ 3

−

2 ∂φ ζ ∂ 2φ 1 ∂ 2φ (2N −3)/2(N −1) ∂ φ − φ =0 + τ + τ 2 ∂ζ 2 2(N − 1) ∂ζ ∂τ ∂ζ ∂ζ 2

(19)

together with the boundary conditions φ = 0,

∂ 2φ ∂ζ 2



∂φ =− ∂ζ

N on ζ = 0

(20)

(the outer boundary condition is relaxed at this stage). Then, with N > 3/2, we look for a solution by expanding φ(ζ, τ ) = φ0 (ζ ) + τ (2N −3)/2(N −1) φ1 (ζ ) + · · ·

(21)

The leading-order problem is φ0 −

123

ζ  1 φ − φ = 0 2 0 2(N − 1) 0

(22)

Natural Convection Boundary-Layer Flow

403

where primes now denote differentiation with respect to ζ . The solution to Eq. 22 which satisfies the boundary conditions in (20) can be expressed in terms of Confluent Hypergeometric Functions Slater (1960) as   1 1 ζ2  φ0 (ζ ) = A0 U ; ; (23) 2(N − 1) 2 4 where A0 is a constant determined from the boundary condition on ζ = 0 in terms of Gamma Functions Copson (1935) by

 N   2(NN−1)

(24) A0N −1 = π (N −1)/2  2(N1−1) From (23) we have φ0 (ζ ) =

A0 (N − 1) ζU (N − 2)



1 3 ζ2 ; ; 2(N − 1) 2 4

 (N  = 2)

(25)

with A0 given by (24). Evaluating expression 25 for ζ large and then using (18) we have that f ∼

A0 (N − 1)21/(N −1) (N −2)/(N −1) + ··· y (N − 2)

(N  = 2)

(26)

We note that, for N = 2, expression 26 involves a term in log y. This leads to an outer region in which y is left unscaled. We are unable to determine the leading-order term in this outer region explicitly, knowing only that it satisfies (26) for y  1 and the outer boundary condition for y 1. This indeterminacy in the outer region arises from the asymptotic nature of this solution, as described in Stewartson (1955) and has arisen in similar situations, Merkin and Mahmood (1998); Merkin (1994) for example. From (18) and (23) we have that √   ∂f A0 π −1/2(N −1)  −1/2(N −1)

θ0 = ∼τ φ0 (0) = τ as τ → 0 (27) ∂ y y=0  2(NN−1) showing that θ0 becomes large as t → ts . To check our analysis we should, from (27), get −2(N −1) straight line behaviour for t close to ts if we plot θ0 against t. We illustrate this in −2 Fig. 4 for N = 2 where we plot θ0 against t. We do see the required straight line form as t approaches ts .√ Further analysis of the data for this case gives a slope of approximately 2.757 estimating A0 π, from (27), as 0.602 giving A0 0.340. This is in reasonable agreement with its value of π −1 from (24). This approach does not work at the critical value of N = 3/2. Our numerical results for the time-dependent problem (7– 9) still show that the solution approaches a singularity at a finite time as t increases, as can be seen in Fig. 5a where we plot θ0 against t for N = 3/2. This figure shows that the numerical solution again breaks down at a finite value of t, with θ0 again becoming large. To obtain the nature of the solution near this singularity we can still apply transformation (18), now putting y f = τ −1/2 φ(ζ, τ ), ζ = 1/2 (28) τ The leading-order problem is now φ0 −

ζ  φ − φ0 + φ0 φ0 = 0, 2 0

φ0 (0) = 0, φ0 (0) = −[φ0 (0)]3/2

(29)

123

404

J. H. Merkin, I. Pop 0.8

-2 0

0.7 0.6 0.5 0.4 0.3 0.2 0.1

t

0.0 0.1

0.15

0.2

0.25

0.3

0.35

0.4

Fig. 4 A plot of θ0−2 against t for N = 2 obtained from the numerical solution of Eq. 8 for the unsteady stagnation point flow. Here • corresponds to values of t used in the numerical integration

A numerical solution to (29) also subject to the extra condition that φ0 → 0 as ζ → ∞ gives φ0 (0) = 1.1232. To check the validity of this approach we plot θ0−1 against t in Fig. 5b. This gives a straight line behaviour for t  0.4, with a slope we estimate as 0.884, giving φ0 (0) 1.131 in good agreement with the numerical solution of (29). 4 Vertical Surface 4.1 Similarity Solutions Here S(x) = 1 in Eq. 5 and to consider this case we make the transformation ψ = x (N −2)/(2N −3) F(η),

η=

y

x

, (N −1)/(2N −3)

 N =

3 2

 (30)

This results in the problem F  +

N −2 1 F F  + F 2 = 0 2N − 3 2N − 3

(31)

subject to the boundary conditions F(0) = 0,

F  (0) = −F  (0) N ,

F  → 0 as η → ∞

(32)

where primes now denote differentiation with respect to η. Before considering Eq. 31 for general values of N , we note that for N = 1 it has the solution F(η) = 1 − e−η

123

(33)

Natural Convection Boundary-Layer Flow

405

(a) 100 0

80

60

40

20

t

0 0.0

(b)

0.2

0.4

0.6

0.8

1.0

1.2

-1 0

0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1

t

0.0 0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1.0

1.1

Fig. 5 a The wall temperature θ0 plotted against t for N = 3/2 obtained from the numerical solution of Eq. 8 for the unsteady stagnation point flow, showing the development of a singularity at finite time. b A plot of θ0−1 against t for N = 3/2. Here • corresponds to values of t used in the numerical integration

Following from the previous case, we can again take F  (0) = θ0 , with θ0 to be found and make a transformation equivalent to (12). This results in the problem F  +

N −2 1  F F  + F 2 = 0, 2N − 3 2N − 3

F(0) = 0, F  (0) = − 1,

F  → 0 as η → ∞ (34)

primes denoting differentiation with respect to η. The problem given in (34) has been considered previously by Merkin and Zhang (1992); Wright et al. (1996) and has been, in the 3 present notation, shown not to have a solution for N in the range ≤ N ≤ 3. A graph of 2

123

406

J. H. Merkin, I. Pop

(a) _ F’(0)

6

5

4

3

2

1

N -10

-5

0

5

10

(b) 0

1.2

1.0

0.8

0.6

0.4

0.2

N 0.0 -10

-5

0

5

10

15

20

Fig. 6 a A plot of F  (0) obtained from a numerical solution of (34) plotted against N . b The wall temperature

θ0 for the flow along a vertical surface calculated from (36) and the values of b0 (N ) obtained from the solution to (34) plotted against N

F  (0) plotted against N obtained from the numerical solution of (34) is shown in Fig. 6a, with this graph showing, as reported in Merkin and Zhang (1992); Wright et al. (1996), the solution becoming singular as N → 3/2 (from below) and as N → 3 (from above). We then have, following from (14), that 2N /3

θ0 ≡ F  (0) = θ0



2N /3

F (0) = θ0

b0 (N )

(35)



where b0 = F (0), now depending on N , is shown in Fig. 6a. From (35)   3 N = θ0 = b0 (N )3/(3−2N ) 2

123

(36)

Natural Convection Boundary-Layer Flow

407 

as before. We can use expression 36 and the values of b0 = F (0) obtained from (34) to calculate the wall temperature θ0 for a given value of N . The results are shown in Fig. 6b with a plot of θ0 against N . From expression 36 and the results shown in Fig. 6a, θ0 → 0 as N → 3/2 and as N → 3, the limits of existence of solutions to (34), as can be seen in Fig. 6b. Also θ0 achieves a maximum value of 1.3067 at N 0.31, with again θ0 → 1 as N → ±∞. The nature of the solution to (34) as N → 3/2 and as N → 3 in the present terminology was considered in Merkin and Zhang (1992). From these results we have 1/3  3 3 −N (from below) F  (0) ∼ 1.3452 + ··· as N → 2 2 (37) F  (0) ∼ 0.7211 (N − 3)−2/3 + · · ·

as N → 3

(from above).

The behaviour given in (37) can be seen in Fig. 6a. Expressions 37 then show that   1 3 (log(3 − 2N ) + 0.1956) + · · · as N → (from below) θ0 ∼ exp (3 − 2N ) 2 (38) θ0 ∼ 1.3867 (N − 3)2/3 + · · ·

as N → 3

(from above).

The rapid exponential fall in θ0 as N → 3/2 and the somewhat slower algebraic decrease in θ0 as N → 3 can be seen in Fig. 6b. 4.2 Flow Development Questions then arise as to how the similarity solution given by (31, 32) for N < 3/2 and N > 3 might appear in a more general boundary-layer flow and also what might happen in such a situation when N is in the range 3/2 ≤ N ≤ 3. Clearly in the first case any solution to (31, 32) could arise at the start of a boundary-layer flow, i.e. as the leading term in an expansion for small x, with the solution developing from this as x increases. A more interesting question is, can this similarity solution arise as an asymptotic (large x) solution in a boundary-layer flow. To get some idea about this we used transformation (30), now with F = F(x, η), leading to the equation

  ∂3 F ∂ F ∂2 F ∂F 2 N −2 1 ∂2 F ∂ F ∂2 F + + =x F − (39) 2N − 3 ∂η2 2N − 3 ∂η ∂η ∂ x∂η ∂ x ∂η2 ∂η3 subject to the boundary condition     ∂F N ∂2 F −βx −βx =− e + (1 − e ) ∂η ∂η2

on η = 0

∂2 F

(40)

= −1 for x small ∂η2 η=0 and boundary condition (32) for x large. Thus the numerical solution of (39, 40) should change from the results shown in Fig. 6a for x small to those shown in Fig. 6b for x large. We find that this is the case when N < 3/2, as can be seen in Fig. 7 where we plot θ0 = Fη (x, 0) against x for representative values of N . The approach to the similarity solution (31, 32) becomes slower as N becomes closer to the critical value of N = 3/2, as might be expected. where β is a constant. This form of boundary condition gives

123

408

J. H. Merkin, I. Pop 1.4

0.5 -1.0

1.2

1.0 0

0.8

0.6

1.25

0.4

0

2

4

6

8

10

12

14

16

18

x 20

Fig. 7 The wall temperature θ0 plotted against x for values of N < 3/2, labelled on the figure, obtained from the numerical solution of (39, 40) with β = 1.0 in (40) for the flow along a vertical surface, showing the solution approaching the similarity solution (31, 32) for x large

However, for N > 3 we find that the numerical solution breaks down at a finite value x s of x with θ0 becoming large. We illustrate this situation in Fig. 8a with a plot of θ0 against x for the case N = 4, where we can see that θ0 increases very rapidly as the solution develops. We were unable to continue our numerical solution to larger values of x than those shown in the figure. We also plot δ = F(x, ∞) in Fig. 8a and we see that δ remains virtually constant, increasing only very slightly as x increases. We found similar behaviour, namely θ0 becoming large and δ remaining almost constant, in our numerical integrations for all the other values of N > 3 tried. We found that the breakdown in solution, i.e. the value of x s , decreased as N was increased. To determine the nature of this singularity that arises when N > 3 we put ξ = x s − x and look for a solution valid for ξ small. Expressions (30) suggest the transformation y Y = (N −1)/(2N −3) (41) ψ = ξ (N −2)/(2N −3) G(ξ, Y ), ξ When (41) is substituted into Eqs. 5, 6 and a solution then sought by expanding in powers of ξ , we find the the problem for the leading-order term G 0 (Y ) is G  0 −

N −2 1 G 0 G 0 − G 2 = 0 2N − 3 2N − 3 0

(42)

subject to the boundary conditions G 0 (0) = 0,

G 0 (0) = −G 0 (0) N

(43)

the outer boundary condition being relaxed at this stage. The solution to Eq. 43 has the form G 0 (Y ) ∼ B0 Y (N −2)/(N −1) + · · ·

as Y → ∞

(44)

for some constant B0 which depends on N . Hence we should regard (41 – 43) as the solution in an inner region. From (41, 44) we then have ψ ∼ B0 y (N −2)/(N −1) + · · ·

123

(45)

Natural Convection Boundary-Layer Flow

409

(a) 6.0 5.5 0

5.0 4.5 4.0 3.5 3.0 2.5

x

2.0 0.0

0.02

0.04

0.06

0.04

0.06

0.08

0.1

(b) 0.008

-5

0

0.007 0.006 0.005 0.004 0.003 0.002 0.001

x

0.0 0.0

0.02

0.08

0.1

Fig. 8 a The wall temperature θ0 and δ = F(x, ∞) plotted against x for N = 4 obtained from the numerical solution of (39, 40) with β = 0.2 in (40) for the flow along a vertical surface, showing the solution approaching a singularity at a finite value of x. b A plot of θ0−5 against x for N = 4. Here • corresponds to values of x used in the numerical integration

at the outer edge of the inner region. This, and further consideration of the solution in the inner region, leads to a regular expansion ψ(ξ, y) = ψ0 (y) + ξ ψ1 (y) + · · ·

(46)

in the outer region. Here the function ψ0 (y) cannot be determined explicitly and depends, in part, on how the solution has developed from the leading edge. It also satisfies (45) for y small and has ψ0 (y) → 0 as y → ∞. This structure is similar to that seen previously in other boundary-layer solutions as they approach a singularity, see Merkin and Mahmood (1998); Merkin (1994) for example and the time-dependent solution near a stagnation point

123

410

J. H. Merkin, I. Pop

described above. This structure close to the singularity at x = x s should perhaps be expected as it appeared originally in the Goldstein–Stewartson singularity at boundary-layer separation Stewartson (1958). The form of the singularity is consistent with the results shown in Fig. 8a for N = 4. From (41) θ0 should become large, of O((xs − x)−1/5 ) as x → xs , whereas δ should remain finite, as can be seen in the figure. In Fig. 8b we plot θ0−5 against x. This should give straight line behaviour, at least close to the singularity, as can be seen in the figure where θ0−5 appears to decrease almost linearly for all x < xs . Although for N > 3 the numerical solution of (39, 40) breaks down at a finite value xs of x with θ0 becoming large, the values of xs are relatively small even for small values of β. Thus this breakdown in solution occurs well before the temperature-dependent boundary condition (1) has come fully into play. In fact from Fig. 8a, xs 0.108 giving, with β = 0.2, e−xs 0.979. Thus the boundary condition has changed by only a small amount from the prescribed heat flux condition at x = 0, for which there is a well-defined solution (Fig. 6a). This suggests that even small perturbations (in x) to this heat flux condition of the type given in (40) can cause the solution to break down when N > 3, with the nature of the solution near the singularity being given by (41, 42) though not necessarily satisfying boundary condition (43). No such problems are encountered when N < 3/2. When 3/2 ≤ N ≤ 3 we cannot use the above approach as there is now no solution at x = 0. To get some idea of what might happen in this case we put ψ = xh f (x, y) in Eqs. 5, 6, resulting in

 2 ∂ 3h ∂h ∂h ∂ 2 h ∂ 2h ∂h ∂ 2 h +h 2 − =x − (47) ∂y ∂ y ∂ x∂ y ∂ x ∂ y2 ∂ y3 ∂y and apply the boundary conditions   N  ∂h ∂ 2h −βx N −1 h = 0, =− e +x 2 ∂y ∂y

on y = 0,

∂h → 0 as y → ∞ ∂y (48)

From (47, 48) the solution starts at x = 0 with h = 1 − e−y with the effects of the temperature-dependent boundary condition (1) becoming more significant as x increases. Our numerical integration of (47, 48) again broke down at a finite value xs of x. We illustrate   this in Fig. 9 with plots of θ0 against x for representative values of N , here θ0 = ∂h x . In each case θ0 grows rapidly as xs is approached, with the value of xs increasing ∂ y y=0 as N is increased. For these plots we used a relatively large value for β of β = 10.0 so that the initial prescribed wall heat flux effect quickly reduced. We found it easier to compute a solution for the smaller values of N > 3/2 in this case with θ0 achieving higher values in the numerical solution. Behaviour similar to that shown in Fig. 9 was found for the other values of N and β tried. 4.3 Critical Values, N = 3/2, N = 3 We consider first the case N = 3/2. We have already identified that N = 3/2 is a critical case when deriving a dimensionless version for our model, see expression 2. In this case to gβ K k 2 make the model dimensionless we are forced to take a length scale  as  = and are αm νγ 2

123

Natural Convection Boundary-Layer Flow

411

0

25

1.6

20

15

10

1.75

5

2.0 0 0.0

0.1

0.2

0.3

0.4

x

0.5

Fig. 9 The wall temperature θ0 plotted against x for the values of N labelled on the figure, obtained from the numerical solution of (47, 48) with β = 10.0 for the flow on a vertical plate, showing the development of a singularity at finite value of x

then unable to specify the velocity and temperature scales U0 and Tref independently. They are related by U0 =

gβ K Tref ν

(49)

This means that, in effect, we could specify some scale for the wall temperature as Tref with the velocity scale given by (49). The result is that any transformation ˜ T = a 2 T˜ , ψ = a ψ,

y˜ = ay

(a constant, a > 0)

(50)

leaves the dimensionless problem invariant. Here, since we have a length scale forced by the model, we have to leave x unscaled. Thus, for N = 3/2, we require some further constraint to specify the problem uniquely. A consequence of this indeterminacy is that a similarity solution can be sought in the form ψ = eβx ,

ρ = eβx y

(51)

leading to    + β   − 22 = 0, (0) = 0,  (0) = − (0)3/2 ,  → 0 as ρ → ∞ (52) where primes now denote differentiation with respect to ρ. We note that the constant β cannot be scaled out of (52) and that, since 3β

∞ ( )2 dρ =  (0)3/2

(53)

0

we must have β > 0.

123

412

J. H. Merkin, I. Pop

However, the problem given by (52) is not unique as we can specify  (0) arbitrarily. Suppose that  (0) = θw (say) where we can assume that θw is given. Then, putting  = θw1/2 β −2/3 

ρ = θw1/2 β 1/3 ρ

(54)

transforms (52) into 







 +   − 22 = 0, (0) = 0,  (0) = −1,  → 0 as ρ → ∞

(55)

where primes now denote differentiation with respect to ρ. The numerical solution to (55)   gives  (0) = 0.84746 = c0 . Since, from (54),  (0) = θw β −1/3  (0) = θw β −1/3 c0 and 3  −1/3 having assumed that  (0) = θw , we have c0 β = 1 or β = c0 = 0.6086, with θw still undetermined. We now examine the case when N = 3. We solved (47, 48) numerically for N = 3 with plots of θ0 and δ = xh(x, ∞) against x being shown in Fig. 10a. Here we took β = 10.0 to allow the prescribed wall heat flux part of (48) to reduce rapidly as x increased. For this case we see that the numerical solution proceeds to much larger values of x than for the smaller values of N , compare Figs. 9 and 10a, before breaking down at a finite value of x, here xs 6.476 compared (say) to xs 0.346 for N = 1.75, with θ0 increasing very rapidly as x gets close to xs . The values of δ increase almost linearly with x and do not show any sign of developing a singularity. For this value of xs the temperature-dependent part of the boundary condition is dominant so that it is this effect that is clearly forcing the singularity. The nature of the singularity in this case is still given by (41–43) and, for N = 3 we can obtain an explicit solution

Ai  (Y ) + K 0 Bi  (Y ) Ai  (−d0 ) 2/3 G 0 (Y ) = −6 b0 , K0 = −  (56) Bi (−d0 ) Ai(Y ) + K 0 Bi(Y ) where G 0 (0)

= b0 ,

b0 Y = 1/3 6



1 Y− 2 b0

,

d0 =

1 b0 61/3

and Ai, Bi are Airy functions Jeffreys and Jeffreys (1962). We are unable to determine the constant b0 directly and this reflects the asymptotic nature of this solution Stewartson (1955). 1/2 + · · · for Y large and so from (41) ψ ∼ 62/3 b0 y 1/2 + · · · From (56), G 0 ∼ 62/3 b0 Y for y small in the outer region, in agreement with (44). Also, from (41), θ0 ∼ b0 ξ −1/3 for ξ = xs − x small, so that a plot of θ0−3 against x should show straight line behaviour, at least for values of x close to xs . We show this plot in Fig. 10b, confirming this prediction and giving an estimated value for b0 as b0 0.218.

5 Conclusions We have considered the free convection boundary-layer flow that can result within a fluidsaturated porous medium when the wall heat transfer is related to the wall temperature by the power-law relation (1). We have considered two cases, stagnation-point flow and flow on a vertical surface. In both cases we were able to reduce the steady problem to similarity form, using respectively transformations (10) and (30), and then make a further transformation (12) which reduced the problem to standard free convection problems that have been treated previously, namely (13) and (34). This enabled us to express the (dimensionless) wall temperature θ0 in terms of known solutions. The results are shown in Fig. 1 for stagnation-point

123

Natural Convection Boundary-Layer Flow

413

(a) 1.6 1.4 1.2 1.0 0.8 0

0.6 0.4 0.2

x

0.0 0

1

2

3

4

5

6

(b) 300

-3 0

250

200

150

100

50

x

0 4.0

4.5

5.0

5.5

6.0

6.5

Fig. 10 a The wall temperature θ0 and δ = xh(x, ∞) plotted against x for N = 3 obtained from the numerical solution of (47, 48) with β = 10.0 in (48) for the flow along a vertical surface, showing the solution approaching a singularity at a finite value of x. b A plot of θ0−3 against x for N = 3. Here • corresponds to values of x used in the numerical integration

flows and in Fig. 6b for flows on a vertical surface. It also identified values of N , N = 3/2 for stagnation-point flows and 3/2 ≤ N ≤ 3 for vertical surfaces, where there was no solution to our original problem. For the stagnation-point flow we also considered a time-dependent problem (7–9) which showed that, for N < 3/2, the time-dependent solution approached the corresponding similarity solution given by (11) for large times, see Fig. 2. However, for N ≥ 3/2 the numerical solution to the time-dependent problem developed a singularity at a finite time ts , with the wall temperature θ0 becoming large as t approached ts , see Figs. 3 and 5a. The value of ts was seen to decrease as the value of N increased. The nature of the singularity at

123

414

J. H. Merkin, I. Pop

t = ts was shown to be given by transformation (18), leading to a wall temperature in the form given by (24) and (27) for N ≥ 3/2, confirmed by the plots shown in Figs. 4 and 5b. We considered two types of boundary condition for the vertical surface problem. For N < 3/2 and N > 3 we applied transformation (30) and solved Eq. 39 subject to boundary condition (40). The numerical solution started at x = 0 with the appropriate prescribed heat flux solution (Fig. 6a) and proceeded to the corresponding similarity solution (in Fig. 6b) for N < 3/2, see Fig. 7. For N > 3 we saw that the numerical solution quickly broke down, with the wall temperature θ0 becoming large (Fig. 8a). The nature of this singularity appeared to be described by (41, 42), Fig. 8b. For 3/2 ≤ N ≤ 3 we had to use a different approach as there is no prescribed heat flux solution for N in this range. We now solved (47, 48), where there is a prescribed wall heat flux solution at x = 0, again finding that the numerical broke down at a finite value for x (Figs. 9, 10a) with the wall temperature θ0 becoming large. This gives the strong suggestion that the boundary-layer flow on a vertical surface cannot sustain even small perturbations to a prescribed heat flux boundary condition along the lines of that given in (1) when N ≥ 3/2.

References Chaudhary, M.A., Merkin, J.H.: Free convection stagnation point boundary layers driven by catalytic surface reactions: II times to ignition. J. Eng. Math. 30, 403–415 (1996) Chaudhary, M.A., Merkin, J.H., Liñan, A.: Free convection boundary layers driven by exothermic surface reactions: critical ambient temperatures. Eng. Ind. 5, 129–145 (1995) Copson, E.T.: Theory of Functions of a Complex Variable. Oxford University Press, Oxford (1935) Ingham, D.B., Pop, I. (eds.): Transport Phenomena in Porous Media III. Elsevier Science, Oxford (2005) Jeffreys, H., Jeffreys, B.S.: Methods of Mathematical Physics. Cambridge University Press, Cambridge (1962) Lesnic, D., Ingham, D.B., Pop, I.: Free convection boundary-layer flow along a vertical surface in a porous medium with Newtonian heating. Int. J. Heat Mass Transf. 42, 2621–2627 (1992) Magyari, E., Pop, I., Postelnicu, A.: Effect of the source term on steady free convection boundary layer flows over a vertical plate in a porous medium, Part I. Transport. Porous Media 67, 49–67 (2007) Mahmood, T., Merkin, J.H.: The convective boundary-layer flow on a reacting surface in a porous medium. Transport. Porous Media 32, 285–298 (1998) Mealey, L., Merkin, J.H.: Free convection boundary layers on a vertical surface in a heat generating porous medium. IMA J. Appl. Math. 73, 231–253 (2008) Merkin, J.H.: Natural-convection boundary-layer flow on a vertical surface with Newtonian heating. Int. J. Heat Fluid Flow 15, 392–398 (1994) Merkin, J.H.: Free convection boundary layers on a vertical surface in a heat generating porous medium: similarity solutions. Quart. J. Mech. Appl. Math. 61, 205–218 (2008) Merkin, J.H.: Natural convection boundary-layer flow in a heat generating porous medium with a prescribed heat flux. Z. Angew Math. Phys. 60, 543–564 (2009) Merkin, J.H., Chaudhary, M.A.: Free convection boundary layers with temperature dependent boundary conditions. J. Theo. Appl. Fluid Mech. 1, 85–100 (1996) Merkin, J.H., Mahmood, T.: Convective flows on a reactive surfaces in a porous media. Transport. Porous Media 33, 279–293 (1998) Merkin, J.H., Pop, I.: Free convection near a stagnation point in a porous medium resulting from an oscillatory wall temperature. Int. J. Heat Mass Transf. 43, 611–621 (2000) Merkin, J.H., Zhang, G.: The boundary-layer flow past a suddenly heated vertical surface in a saturated porous medium. Wärme-Und Stoffübertragung 27, 299–304 (1992) Nield, D.A., Bejan, A.: Convection in Porous Media (3rd edn). Springer, New York (2006) Slater, L.J.: Confluent Hypergeometric Functions. Cambridge University Press, Cambridge (1960) Stewartson, K.: On asymptotic expansions in the theory of boundary layers. J. Math. Phys. 13, 113–122 (1955) Stewartson, K.: On the Goldstein theory of laminar separation. Quart. J. Mech. Appl. Math. 11, 399–410 (1958) Vafai, K. (ed.): Handbook of Porous Media (2nd edn). Taylor and Francis, New York (2005) Wright, S.D., Ingham, D.B., Pop, I.: On natural convection from a vertical plate with a prescribed surface heat flux in porous media. Transport. Porous Media 22, 181–193 (1996)

123