International Journal of Engineering Science 39 (2001) 1881±1897
www.elsevier.com/locate/ijengsci
Eect of surface conditions on ¯ow of a micropolar ¯uid driven by a porous stretching sheet N.A. Kelson a, A. Desseaux b,* a
Centre in Statistical Science and Industrial Mathematics, Queensland University of Technology, Australia b Institut Universitaire de Technologie, GMP, Le Mount Houy, 59313 Valenciennes Cedex 9, France Received 10 July 2000; received in revised form 4 October 2000; accepted 15 November 2000
Abstract Self-similar boundary layer ¯ow of a micropolar ¯uid driven by a stretching sheet with uniform suction or blowing through the surface is considered. A perturbation analysis is used to derive closed form solutions, and a number of numerical solutions are used to validate the analysis. In order to investigate the eects of dierent microrotation boundary conditions, results are obtained here which prescribe a ®xed ratio between the microrotation and the shear stress at the surface. Ó 2001 Elsevier Science Ltd. All rights reserved.
1. Introduction The micropolar model of ¯uid ¯ow has attracted considerable attention from researchers since the fomulation of the model by Eringen [8]. Eringen's micropolar model includes the classical Navier±Stokes equations as a special case, but can cover, both in theory and applications, many more phenomena than the classical model. Examples of industrially relevant ¯ows that can be studied using the theory include the ¯ow of low concentration suspensions, liquid crystals, blood, porous media, lubrication and turbulent shear ¯ows. Extensive reviews of the theory and applications can be found in the articles by Ariman et al. [1,2] and the recent book by Lukaszewicz [16]. The potential importance of micropolar boundary layer ¯ow in industrial applications has motivated a number of previous studies, of which those of Chiam [3], Hady [10], Heruska et al. [13] and Hassanien and Gorla [12] are most directly relevant to the present case where the ¯ow is driven by a porous stretching sheet. These studies obtained solutions for a limited range of values
*
Corresponding author. Tel.: +03-27-29-13-18; fax: +333-27-51-14-39. E-mail addresses:
[email protected] (N.A. Kelson),
[email protected] (A. Desseaux).
0020-7225/01/$ - see front matter Ó 2001 Elsevier Science Ltd. All rights reserved. PII: S 0 0 2 0 - 7 2 2 5 ( 0 1 ) 0 0 0 2 6 - X
1882
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
Nomenclature j micro-inertia density (m2 ) s microrotation boundary condition parameter u; v cartesian velocity components (m s 1 ) x; y cartesian coordinates parallel and normal to sheet (m) C1;2;3 =N1;2;3 dimensionless parameters dimensionless streamfunction and microrotation F1 ; F2 D; N stretching factor and microrotation (s 1 ) dimensionless mass transfer parameters V1 ; b g dimensionless normal distance from surface l; j dynamic viscosity and coupling coecient (kg m 1 s 1 ) m kinematic viscosity (m2 s 1 ) microrotation/spin-gradient viscosity (m kg s 1 ) ms q ¯uid density (kg m 3 ) for the physical parameters that govern the ¯ow, and, inevitably, the conclusions drawn about the ¯ow behaviour are necessarily of limited validity. In recent works by the present authors [6,15], similarity solutions for the present ¯ow were analysed using both numerical methods and a perturbation analysis. The latter was used to obtain solutions for an extended range of parameters not previously considered, and certain limiting solutions to the governing equations were also derived. In the studies of this ¯ow mentioned above, it was assumed that the microrotation vector is zero on the solid surface. Only Hassanien [11] considered alternative conditions in his study of ¯ow due to an impermeable stretching sheet. In order to investigate the eect of dierent surface conditions, the purpose of this work is to obtain closed form solutions (both exact and approximate) for micropolar ¯ow over a stretching sheet with a ®xed ratio of the microrotation component and the shear stress at the surface. It is desirable to compare the model predictions against experimental data, but we are unaware of any for the present ¯ow. Instead, we compare our results with the corresponding results for Newtonian ¯ow in order to highlight the dierences between the micropolar and classical models. To these ends, this paper is organised as follows. In Section 2 the governing equations and similarity transformations are presented, and appropriate values for the physical parameters are discussed. In Section 3, a perturbation analysis is used to derive closed form solutions to ®rst order, and the numerical methods used here are brie¯y described in Section 4. The results are discussed in Section 5, and an analysis of the near-wall behaviour for this complex ¯ow is given. Finally, the conclusions of this study are summarised in Section 6.
2. De®ning equations Consider ¯ow in the region y > 0 driven by a plane porous surface located at y 0 with a ®xed end at x 0. We suppose that the surface is stretched in the x-direction such that the x-com-
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
1883
ponent of the velocity varies linearly along it, i.e., u
x; 0 Dx where D is an arbitrary constant and D > 0. A uniform velocity V0 through and normal to the stretching surface is also considered. The simpli®ed two-dimensional equations governing the ¯ow in the boundary layer of an isothermal, steady, laminar, incompressible micropolar ¯uid in an otherwise quiescent medium are (see e.g., [18]) 2 ou ou j ou j oN ;
1 u v m 2 ox oy q oy q oy oN oN j ou ms o2 N u v 2N :
2 ox oy qj oy qj oy 2 Compared with classical Newtonian ¯uids, the governing equations include the microrotation or angular velocity N whose direction of rotation is in the x±y plane and the material constants j; j and ms . When these constants are zero, the governing equations reduce to the usual Navier±Stokes equations for boundary layer ¯ow. The appropriate physical boundary conditions are u
x; 0 Dx;
v
x; 0 V0 ;
N
x; 0
s
ou on y 0; oy
3a
and u ! 0;
N ! 0 as y ! 1;
3b
where V0 < 0 corresponds to suction, and V0 > 0 corresponds to injection. There is no concensus as to what boundary condition should be used for the microrotation, and a number of plausible alternatives are discussed in [16]. In order to investigate the eect of different surface conditions for the microrotation, a linear relationship between N and the surface shear stress is used here, where the boundary parameter s varies in the range 0 6 s 6 1. The governing equations (1) and (2) can be expressed in a simpler form by introducing the physical parameters N1 j=
qm; N2 Dms =
qm2 and N3 Dj=m and a stream function W such that u oW=oy and v oW=ox. Speci®cally, we introduce the transformations p g D=m
1 N1 y; p
4 W Dm
1 N1 xF1
g; p N D3 =m
1 N1 xF2
g; into (1) and (2) to obtain
F10 2
F1 F100 F1000
F10 F2
F1 F20
N1 F 0; 1 N1 2
N1 00 N2
F1 2F2 F 00 : N3 N3
1 N1 2
5
6
1884
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
For comparison with previous studies, a rescaling of parameters using C1
N1 ; 1 N1
C2
N2 N1
1 N1
and C3
N3 1 N1
7
is convenient. Eqs. (5) and (6) can be written as F1000
F1 F100
F10 2
C1 F20 ;
C3
F1 F20 C1
C2 F200 F100 2F2
8 F2 F10
and the appropriate transformed boundary conditions are r 1 C1 0 V0 V1 ; F2
0 sF100
0; F1
0 1; F1
0 Dm lim F10
g 0;
g!1
lim F2
g 0;
g!1
where V1 < 0 corresponds to suction, and V1 > 0 corresponds to injection. For micropolar boundary layer ¯ow, the wall shear stress is given by [17] ou : sw q
m j=q jN oy y0
9
10a
10b
11
Using U0 Dx as a characteristic velocity scale in (11), the dimensionless wall shear stress or skin friction coecient Cf is Cf
sw 1 sC1 p Rx 1=2 F100
0; 2 qU0 1 C1
12
where Rx Dx2 =m is the Reynolds number based on the ¯uid kinematic viscosity. For C1 0, the above reduces to the corresponding expression for a Newtonian ¯uid, as originally analysed by Crane [4]. 2.1. Parameter values The problem as de®ned above contains the ®ve parameters s; V1 ; C1 ; C2 and C3 . Of these, we note that C1 may only vary in the range 0 6 C1 < 1 while C2 and C3 must be non-negative for micropolar ¯uids [8]. The range for parameter C1 suggests that it can be used as a perturbing parameter, as described in Section 3. In general, the parameter values may depend on the material properties of the working ¯uid and the ¯ow conditions, and their determination is a dicult matter [14]. To study the micropolar model, and for consistency with many other studies, all parameters are taken as independent apart
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
1885
from the above restrictions. However, with this assumption it should be noted that certain limiting behaviour must also be assumed. For example, the limiting case C1 0 is of interest because it corresponds to a decoupling of Eqs. (8) and (9). For this case, the limits for C2 and C3 =C1 are assumed zero or ®nite. We note further that the limiting value of C1 1 is inadmissible, a point apparently not appreciated in earlier works. This condition, which represents ``in®nite coupling'' between the microconstituents and the bulk ¯ow, leads to critical conditions at the surface. Firstly, from Eq. (10a), C1 1 implies that V1 0 for any arbitrary (but ®nite) mass transfer through the surface. However, both Chiam [3] and Heruska et al. [13] computed numerical solutions using C1 1 and non-zero values for the surface mass transfer parameter. Secondly, from Eq. (12) with s 6 1, the condition C1 1 leads to either an in®nite wall shear stress or F100
0 0. If instead s 1 in (12), then this implies that either the shear stress is identically zero or F100
0 is in®nite. The underlying problem is this: the case of in®nite coupling corresponds to a ¯ow dominated by viscosity, and the assumptions leading to the simpli®ed governing equations for boundary layer ¯ow are unlikely to be valid for C1 1. With regard to the other parameters, we note that only very small values of C3 are of interest, and changes in the prescribed value of this parameter appeared to have only a slight in¯uence on the bulk ¯ow predictions of Chiam [3] and Heruska et al. [13]. Following Hassanien and Gorla [12] and Hady [10], who did not include the left-hand side terms of (2) or (5) at any stage in their analysis, we set C3 0 and consider the four parameter (V1 ; s; C1 ; C2 ) system F1000
2
F1 F100
F10
C1 F20 ;
C2 F200 F100 2F2
13
along with the boundary conditions given in (10a) and (10b). In Table 1, the range of parameter values used in previous studies are given. Clearly, a limited range of C2 values have been investigated, with values near C2 2 being most frequently chosen. It should be noted that Heruska et al. [13] concluded that little change in the ¯ow pattern is found Table 1 Range of parameter values used in other studies Study
s
C1
C2
C3
Chiam [3] Hersuka at al. [13] Hassanien and Gorla [12] Hady [10] Hassanien [11]b
0 0 0 0 0 1/4, 1/2 0 0
0:01 ! 1 0:001 ! 1 0.2 0.2 0:01 ! 1 0.1 0 ! 0:9 0 ! 0:1
1!4 1 ! 10 2 2 1 ! 10 1 2 0:125 ! 2
0 ! 0:05a 0 ! 0:05 0 0 0 ! 0:05 0 0 0
Kelson et al. [15] Desseaux and Kelson [6] a
V1
0 0 0 0
4!0 5!5 0:7 ! 0:7 0:7 ! 0
The range shown is for the variable K3 C3 =C1 used by Chiam. Hassanien's results were limited to the case of an impermeable sheet only due to his choice of similarity transforms. Also, Hassanien's de®nitions of C2 and C3 dier from those used here.
b
1886
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
as C2 varies. Similarly, Hassanien [11] concluded that variations in C2 cause only slight variations in the wall shear stress. However, in a recent paper we have shown, via analytical and numerical results, that distinctly dierent behaviour for the stream function and microrotation can be expected for the limiting cases of either small or large C2 [15]. Conceptually, the present work is mainly concerned with the case of C2 small. Our experience indicates [6] that this case, and the case of strong injection, represents the computationally more dicult situation. Here, a perturbation approach is employed, allowing departures from Newtonian ¯ow behaviour due to the coupled (C1 > 0) case to be clearly identi®ed and analysed. 3. Perturbation analysis Choosing C1 as the perturbing parameter, we expand the similarity functions using F1 f0 C1 f1 C12 f2 ;
F2 g0 C1 g1 C12 g2
14
and substitute them into (13). By collecting terms in powers of C1 , a hierarchy of dierential equations for the functions fn and gn are obtained. For fn , these are f0000 f0 f000
f00 2 0;
fn000 f0 fn00
2f00 fn0 f000 fn RHSn
for n > 0;
15b
g00 . The microrotation functions gn are obtained via
where e.g., RHS1 C2 gn00
15a
2gn fn00
for n P 0:
16
The appropriate boundary conditions for Eqs. (15a), (15b) and (16) are f0
0
f00
0 1;
V1 ;
fn0
0 0;
fn
0 0;
g0
0 gn
0
sf000
0;
sfn00
0;
lim f 0 g!1 0 lim f 0 g!1 n
lim g0 0 for n 0; g!1
lim gn 0 g!1
for n > 0:
17a
17b
An exact closed form solution to (15a) for f0 is [15] f0
V1
1
exp
bg b
1 b
x;
18a
where 1 b 2
q exp
bg and x x
g : V1 V12 4 b2
18b
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
1887
Eqs. (18a) and (18b) also describe Newtonian ¯uid ¯ow over a stretching sheet, as originally analysed by Crane [4], Gupta and Gupta [9] and others subsequently. From (4), (7) and (18a), (18b) the leading order solution indicates that the boundary layer has constant thickness yT of order r 1 m 1 yT O ;
19 or in dimensionless terms gT O b D
1 C1 b where the eect of suction (b > 1) or injection (b < 1) is to respectively, thin or thicken the boundary layer compared with the impermeable case. Also, the estimate for the boundary layer width given in (19) increases without bound as C1 ! 1, consistent with our earlier remarks. p Introducing the parameters c
1=b
2=C2 and t
s; c 1 2s
2s=c2 the solution of (16) for the zeroth-order microrotation g0 can be expressed as b3 g0
4s bgx when c 1; 4 b3 c2 g0
x t
s; cb2c 2 xc 2 c2 1
20a when c 6 1:
20b
In the limit as c ! 1, (20b) simpli®es to (20a), as expected. For C1 0 the solutions for f0 and g0 yield exact values for the boundary conditions F100
0 and 0 F2
0 at the surface. From (18a), (18b), and (20a) and (20b), these are F100
0
f000
0
b;
F20
0
g00
0
b2 c c 2 1c
2s
1 c
all c:
21
For n > 0, a closed form general solution to (15b) for the case n 1 is presented below, and is given by the sum f1 hc p1 , where hc is the complementary function (CF) and p1 is a particular integral (PI). The CF can be expressed in terms of three independent solutions as hc a1 hc1 a2 hc2 a3 hc3 ; hc3
hc1 x;
hc2 1 bgx; Z 1=b2 Ei
t 2 exp
x Ei
x 2xEi
x x dt; t x
22
where a1 ; a2 and a3 are constants and Ei
x is the exponential integral function [15]. In order to introduce algebraic simpli®cations, and without loss of generality, a modi®ed de®nition of R 1=b2 Ei
x x
e t =tdt is used here. The PIs for n 1 are found by substituting (20a) and (20b) into (15b) and solving. Following our earlier work [6], integer values of c were taken here in order to greatly simplify the resulting expressions and the PIs were found with the aid of a computer algebra system. The solutions for p1 for c 1; 2; 3; 4 are given in Table 2. For b 1 (an impermeable sheet), these values for c correspond, respectively, to C2 2; 1=2; 2=9 and 1/8. The range considered here can be extended
1888
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
Table 2 Particular integral for n 1 and c 1 ! 4 c
Particular integral p1
1
b p1 3 4
2
p1
2b h 3
x log
x t
s; 2
3
p1
9b h 16
x log
x t
s; 3
24b4 3b4 x2 12b4 log
x 24b4 x log
x
4
p1
8b h 15
x log
x t
s; 4
4s
bg 4b2
144b6
2b2 log
x
i 4b2 x log
x b2 xlog
x2
18b6 x2 b6 x3
72b6 log
x
i 6b4 xlog
x2
i 144b6 x log
x 36b6 xlog
x2
to larger integer values for c, but was deemed sucient to analyse the eect of changing this parameter on the nature of the solution. Once the CF and PI are known, the values of the constants in the general solution can be determined, and the required values are given in Table 3. Note that in the table, and subsequently, a0 1 exp
1=b2 and b0 log
1=b2 . The functions and constants given in (22) and Tables 2 and 3 completely specify the solution of the ®rst-order stream function f1 for the values of c considered here.
4. Numerical method Two methods were used for the numerical solution of (13) subject to the boundary conditions given in (10a) and (10b), namely, a shooting method using a fourth-order Runge±Kutta algorithm, and a quasilinearisation scheme. For the speci®cs of the quasilinearisation method, see [5,6]. For both methods, the convergence criteria and step sizes for the integration were taken small enough so that results with discretisation error of less than 0:5 10 3 % (i.e., more than ®ve signi®cant digits) could be reported. To initiate the shooting method, initial guesses for the surface conditions F100
0 and F20
0 are needed, and for the quasilinearisation scheme, an initial approximate solution is required. For the majority of results obtained here, we used the solutions given in (18a), (18b), (20a), (20b) and (21) to provide initial approximate values. In both schemes, the boundary conditions (10b) were implemented on a ®nite domain of length gc . For some cases, convergence was best achieved by initiating computations on a given domain length gc using solution estimates obtained on shorter domains. The domain length was then incrementally increased by Dgc P 0:5 until an asymptotically valid solution was obtained or nonconvergence was detected. Using the above strategies, the problem of rapidly increasing solutions leading to over¯ow [12] was almost never observed. Occasionally, however, apparently spurious numerical solutions were obtained, as discussed in our earlier paper [6].
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
1889
Table 3 Values for the constants ai for n 1 and c 1 ! 4 c1
c2
a1
b3 1 4
1 b2
4s a0
2 b2
a2
h b3 1 2 4
1 b
4s a0
a3
b 4
a1
h 2b b0 b2
1 b0 t
s; 2
b2
4b0 2 3
1 b
a2
h 2b 2 3
1 b
c4
1 t
s; 2
2b2
b20 b4
4 4a0 4b0
i 2 b0 2b4
a0 b0
a1
9b b0 b2
1 b0 t
s; 3
16
1 b2
3b2 6b4
a2
9b 16
1 b2
2 b0
a3
27b5 t
s; 3 4
a1
8b b0 b2
1 b0 t
s; 4
15
1 b2
1 t
s; 3
3b2
36b8
4 4a0 4b0
a2
i b20 2b6 a0
4b3 t
s; 2 3
a3
c3
i
8b 2 15
1 b
a3
12b4
1
4b0 b20
6b6
4 4a0 4b0
b20
12b8 a0
12b6
a0 b0
b2 15b4
36b6
1
4b0 b20
b20 72b10 a0
1 t
s; 4
2b2
18b4 72b6
2 b0 72b8
a0 b0
192b7 t
s; 4 5
5. Results and discussion 5.1. Results for C1 0 We begin by discussing the uncoupled case C1 0, as a number of interesting features emerge from a study of the exact solutions given in (18a), (18b) and (20a), (20b). Some of these will be found to be directly relevant to the coupled case C1 > 0. Before proceeding, we document the accuracy and convergence properties of the numerical schemes used here. In the following, we mainly present details for the more frequently used
1890
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
shooting method [3,10,12]. A detailed comparison of the two numerical methods used here is given elsewhere [7]. Initially, numerical results for C1 0 were compared with the exact solutions given in (18a), (18b) and (20a), (20b) for a range of parameter values for b; s and C2 . Computed results with less than 0:5 10 3 % error could be obtained in many cases. However, as C2 or b was decreased, a shorter domain length gc was usually needed to avoid non-convergence. The use of shorter domains can potentially introduce signi®cant cutto errors, and to determine error trends a number of computations using Fortran double precision were performed. In Table 4, the maximum domain size
gc max that could be used for dierent values of the parameters C2 ; b and s is reported. Note that results are tabulated in terms of b
gc max rather than
gc max because the former quantity is, from (19), proportional to the fraction of the boundary layer width traversed by the computational domain. For convenience, the cutto error is reported in terms of the percentage error in the computed value of F1
g at the cuto point gc , as the errors involved in computed quantities near the cutto were always observed to be greater than those in the vicinity of the surface. Note that the tabulated results were found to be the same for s 0; 1=2 and 1 (cf. the case C1 > 0 below). From Table 4, it can be seen that for ®xed b 1, reducing C2 from 2 ! 1=8 is accompanied by a reduction in b
gc max from 23 ! 5:5, and a corresponding rise in cuto error. At C2 1=8, the cuto error of 1.3% implies that computed results are accurate to 2±3 signi®cant digits only. For ®xed C2 2, decreasing b from 1 ! 1=4 results in a similar trend. However, for this sequence b
gc max decreases more rapidly and the loss of accuracy is more severe. For b 1=4, the recorded cuto error of 47% indicates that no signi®cant digits can be guaranteed in the solution arrays away from the surface. From the above discussion, it is clear that for C1 0 the exact solutions given in (18a), (18b) and (20a), (20b) provide a better description of the ¯ow, especially for the computationally more dicult cases of small C2 or b. Subsequently, Eqs. (18a), (18b) and (20a), (20b) were used to obtain the solution curves discussed below. In Fig. 1, the in¯uence of the parameters s; b and C2 on the stream function F1 and microrotation F2 is shown when C1 0. In Fig. 1(a), solution curves for the stream function F1 f0 for suction (b 2) and injection (b 1=2) are compared with the impermeable case. The eect of increasing the suction parameter b is to increase the wall value F1
0 while decreasing the boundary layer width, consistent with Table 4 Cuto error for dierent C2 and b (s 0; 1=2 and 1) c b2 C2
1 2
2 1/2
3 2/9
4 1/8
b1
C2 b
gc max Cuto error (%)
2 23 < 0:5E-03
1/2 12 0.004
2/9 7 0.4
1/8 5.5 1.3
C2 2
b b
gc max Cuto error (%)
1 23 < 0:5E-03
1/2 10.8 0.06
1/3 5.8 5.3
1/4 4.8 47
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
(a)
(b)
(c)
(d)
1891
Fig. 1. Pro®les of stream function F1 f0 and microrotation F2 g0 for C1 0 and dierent b, C2 and s: (a) f0 ; (b) g0
b 1; s 0; (c) g0
b 1; C2 1=8; (d) g0
C2 1=8; s 0. Dashed lines in (b)±(d) are the limiting solution curve F2 f000
g=2 (when C2 ! 0).
Eq. (19). Note that the leading order solution for F1 is independent of s and C2 , and is incapable of capturing the non-Newtonian behaviour often associated with micropolar ¯uids. In Fig. 1(b), the eect of decreasing C2 from 2 ! 1=8 on the solution for F2 g0 is shown for ®xed values of b 1 and s 0. Clearly, as C2 decreases the near-wall maximum increases, and its location moves closer to the surface. In Fig. 1(c), the eect of increasing s from 0 ! 1 on the solution for F2 is shown for the case C2 1=8. Changing this parameter signi®cantly aects the near-wall pro®le. The behaviour of the solution curves given in Figs. 1(b) and (c) is best understood by considering Eq. (16) from which they are obtained. Away from the surface, the term involving C2 in (16) becomes small in the limit C2 ! 0, and the solution behaves like F2 f000
g=2, independent of s and C2 . This limiting far-®eld or outer solution curve is shown as a dashed line in the ®gures. As C2 decreases, the outer solution provides an increasingly better approximation for F2 , and the region of validity extends closer to the surface. However, immediately adjacent to the surface a boundary layer develops due to the boundary condition g0
0 sb. Clearly, the shape of the F2 pro®les within this boundary layer will depend on the value of s, as shown in Fig. 1(c). Finally, in Fig. 1(d), solution curves for the microrotation F2 for suction (b 2) and injection (b 1=2) are compared with the impermeable case. The eect of increasing the suction parameter b is to again decrease the boundary layer width. For the case s 0 shown, a maximum occurs in the near-wall layer. As b increases, the maximum increases, and its location moves closer to the surface.
1892
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
The behaviour of the pro®les in Figs. 1(b) and 1(d) for the case s 0 suggest that the height of the near-wall maximum may approach ± but will not exceed ± the surface value f000
0=2 b=2 of the outer solution. Naturally, the observed boundary layer behaviour raises interesting questions about instability, an analysis of which is outside the scope of this paper. 5.2. Results for C1 6 0 We now consider values of C1 > 0. Compared with the case C1 0, computations using correspondingly shorter domains were generally needed to avoid non-convergence. The need for shorter domains was especially noticeable for non-zero values of the boundary parameter s, presumably due to the increased coupling in both the governing equations and the boundary conditions. In the absence of an exact solution, this would suggest that the cuto errors would be even greater than those reported earlier in Table 4 for the case C1 0. In Table 5, we report some representative results which indicate the above-mentioned trends. The tabulated data for s 0 and the case C2 1=8 reveals that correspondingly shorter domains are always needed when C1 > 0. However, for the larger values of C2 1=2 and 2 given in the table, the maximum domain length could be increased slightly as C1 ! 1. In contrast, the data for C2 2 with s 1=2 or s 1 reveal that drastically shorter domains are required as C1 ! 1. In terms of the expected cuto errors in the results, the tabulated data suggests that the worst scenario for the case s 0 is when C2 ! 0 and 0 < C1 < 0:3, whereas for s ! 1 the worst scenario is C1 ! 1. In the absence of an exact solution, a direct comparison of analytic and numerical results was used to gauge the accuracy of the perturbation approximation. Fig. 2 illustrates one such comparison, where the computed pro®les for the case s 0; V1 0; C1 0:2 and C2 0:5 are compared with the perturbation analysis using s 0; b 1 and c 2. From Fig. 2(a), the computed pro®le for F1 and the solution curve obtained from the ®rst-order perturbation approximation F1 f0 C1 f1 are almost coincident. Although not illustrated here, this behaviour was also noted for C2 2; 2=9 and 1/8 (i.e., c 1; 3 and 4). As a ®nal comment regarding Fig. 2(a), we also noted that corresponding pro®les with s 1=2 and 1 were also virtually coincident close to the wall. However, away from the surface small dierences were observed, which we again attribute to the increased coupling when s > 0. Incidentally, excellent near-wall agreement was also observed in the corresponding pro®les for the case s 0 using much larger values for the perturbing parameter (C1 ! 1, see [15]). For C1 0:9, for example, the perturbation approximation and the computed value of F100
0 diered by less Table 5 Eect of C1 on b
gc max (b 1, s 0 except where indicated)
b
gc max for
C1
0.0
0.1
0.2
0.3
0.5
0.7
0.9
C2 2 (s 0:5) (s 1:0) C2 1=2 C2 1=8
23 23 23 12 5.5
21.5 24 22 11 4.5
23 23.5 22 12 4
23.5 24 22 11.5 4.5
25 25 6 11.5 4.5
25.5 12 4 14 5.0
24.5 6 2.5 14 4.5
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
(a)
1893
(b)
Fig. 2. Solution curves for s 0, V1 0, C1 0:2 and C2 0:5: (a) comparison of computed pro®le of streamfunction F1 (cts) and ®rst order perturbation approximation f0 C1 f1 (dashed); (b) Comparison of computed pro®le of microrotation F2 (cts) and leading order perturbation approximation g0 (dashed).
than 0.02%. Here again, only small dierences in the far-®eld pro®les were noted. Finally, we mention that a reliable comparison could not be made where both C1 ! 1 and s ! 1, due to the likely cuto errors involved in the computations. The totality of our results just discussed suggest that the contribution from the second and higher-order functions (fn ; n > 1) to Eq. (14) is signi®cantly smaller than f1 . In Fig. 2(b), the computed pro®le of the microrotation F2 and the leading order perturbation approximation F2 g0 are compared for the same case as in Fig. 2(a). Immediately next to the wall the two pro®les are virtually coincident, with small dierences becoming noticeable away from the surface. This behaviour was also noted for the cases C2 2; 2=9 and 1/8. Clearly, for C1 0:2 the agreement between the analytic and computed pro®les is very good. In fact, even for C1 ! 1 or s ! 1 agreement remained good. However, somewhat larger dierences in the nearwall pro®les were observed in these cases. For example, for the case s 1; C1 0:2 and C2 2, the maximum dierence between the zeroth-order and computed pro®les for F2 was about 10% (occurring at the surface). Since the primary concern of this work is on the bulk ¯ow behaviour, a continuation of the perturbation analysis to obtain the ®rst-order function g1 did not seem warranted here. As a ®nal comparison, results are reported in Table 6 for the case C1 0:2 and C2 2 examined in the works of Hassanien and Gorla [12] and Hady [10]. Our computed and analytic results for the surface conditions F100
0 and F20
0 are both in excellent agreement with the numerical solution of Hassanien and Gorla [12]. By comparison, Hady's [10] approximate analytic solution is signi®cantly worse. Table 6 Values of F100
0 and F20
0 for b 1, s 0, C1 0:2 and C2 2
F100
0 0 F2
0
Computed (this work)
Analytic (this work)
Computed ([12])
Analytic ([10])
0.99078 0.25120
0.99080 0.25000
0.99081 0.25121
0.97500 0.25316
1894
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
5.3. Analysis of near-wall behaviour We now turn our attention to an investigation of the complex near-wall behaviour of the micropolar ¯ow considered here, using the perturbation solutions for the analysis. In Fig. 3, the pro®les for the ®rst order stream function f1 for the case s 0 and b 1 and dierent values for C2 are compared. Fig. 3(a) gives a clear illustration of the departure of the stream function from the Newtonian pro®le given earlier. Fig. 3(b) shows in more detail the nearwall behaviour of the pro®les. Referring to the ®gure, for the case s 0 and b 1, we observe a narrow near-wall region where f1 P 0 and f100
0 > 0 in each of the pro®les. Physically, this nearwall behaviour corresponds to a reduction in the skin friction, as observed in e.g., ¯uids with extremely small polymeric additives [8]. As C2 decreases, we note that the width of the region where f1 P 0 also decreases, and the location of the maximum moves closer to the surface. In each of the illustrated pro®les, the height of the near-wall maximum is very small and begins to decrease for C2 6 0:5. Although not entirely obvious from the detail given, the value of f100
0 (which is directly proportional to the reduction in skin friction) actually increases as C2 decreases. We shall discuss this behaviour in more detail shortly. In Fig. 4, the eect of changing either the suction parameter b or the boundary parameter s on the pro®le for f1 is shown.
(a)
(b)
Fig. 3. (a) Comparison of f1 pro®les for C2 2; 0:5; 2=9 and 1=8 (s 0, b 1); (b) near-wall behaviour.
(a)
(b)
Fig. 4. Comparison of f1 near-wall pro®les for C2 0:5 and (a) s 0 and b 1=2; 1; 2; (b) b 1 and s 0; 1=4; 1=2; 3=4; 1.
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
1895
In Fig. 4(a), the eect of changing b only on the near-wall pro®le for f1 is illustrated for the case C2 0:5 and s 0. Increasing b (i.e., increasing suction) corresponds to a narrowing of the nearwall region where f1 P 0 in the pro®les. In Fig. 4(b), the eect of increasing s from 0 ! 1 is shown for the the case C2 0:5 and b 1. Of the pro®les given, only the curve for s 0 shows a nearwall region where f1 P 0. As s increases, the value of f100
0 decreases, corresponding to an increase in the skin friction. The above results demonstrate that under some circumstances a stress reduction is observed. To complete our investigation of the near-wall ¯ow, we now turn our attention to a more detailed analysis of the in¯uence of the various parameters on the wall shear stress for micropolar ¯ow (i.e., C1 > 0). To this end, it proves convenient to de®ne a percentage reduction ratio Ps for the wall shear stress as follows: Ps 100
f100
0 f000
0
100
f100
0 : b
23
The quantity C1 Ps gives the percentage change in the wall shear stress compared with Newtonian ¯ow. Here, Ps < 0 corresponds to a reduction, whereas Ps > 0 corresponds to an increase. In Table 7, we report values of Ps for s 0 and various b and C2 values. For a ®xed value of b, the tabulated columns of data indicate that the percentage stress reduction is greater for smaller values of C2 . If instead the value of C2 is ®xed, the tabulated rows of data suggest a more complicated trend. For suction (b > 1), the data suggests that increasing b reduces the percentage shear stress reduction. However, for injection (b < 1), the percentage stress reduction initially increases, but then begins to decrease for stronger injection levels. For practical applications, the trends suggest that if maximum wall stress reduction is desired, then ¯ow con®gurations should be chosen which correspond to small values of C2 ; s 0 and low injection rates through the surface. As a ®nal illustration, we give the near-wall pro®les for f1 in the cases of strong injection or suction. In Fig. 5(a), the near-wall pro®les for f1 for the case s 0 and C2 1=8 are shown, where the suction parameter b is increased from 1 ! 4. As b increases, the near-wall region where f1 P 0 is reduced, and is almost completely extinguished for the strongest suction level, b 4. In Fig. 5(b) the near-wall pro®les for f1 for the case s 0 and C2 2 are shown where b 6 1, i.e., injection. As the injection rate increases (i.e., b decreases), the width of the near-wall region where Table 7 Values of Ps for various b and C2 (s 0) b Ps for C2 9=2 2 1/2 2/9 1/8
4
3
2
3/2
1
2/3
1/2
)4.63
)0.72
)1.21
)2.30 )4.18
)6.14
)4.60 )8.81 )11.7 )13.7
)13.2
1/3
1/4
)5.35 )8.30 )14.9
)8.24
1/6 )3.88
)7.45
1896
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
(a)
(b)
Fig. 5. Comparison of f1 near-wall pro®les for suction (b > 1) and injection (b < 1): (a) s 0, C2 1=8 and b 1; 4=3; 2; 4; (b) s 0, C2 2 and b 1; 1=2; 1=3; 1=4.
f1 P 0 increases. However, the near-wall maximum does not continue to increase, and we note that the value of f100
0 initially increases, then decreases, as b decreases. 6. Conclusion In this paper, a perturbation analysis has been used to study self-similar micropolar ¯ow over a porous stretching sheet. Conceptually, our approach diers from previous studies by focusing on small values of C2 , where boundary layer eects are more pronounced and numerical solutions are more dicult to obtain. Some numerical diculties involved in a solely computational approach have been documented here. Using the analytic solutions, departures from Newtonian ¯ow behavior were identi®ed and the complex near-wall behaviour of this ¯ow was investigated for various conditions at the surface, including dierent mass transfer rates and dierent microrotation boundary conditions. We also investigated the micropolar wall shear stress reduction and enhancement eects as a function of the various surface conditions and a number of trends were identi®ed and discussed. To conclude, we mention one such trend that seems of particular interest. If maximum wall stress reduction is desired in a practical application, then our results suggest that ¯ow con®gurations should be chosen which correspond to small values of C2 ; s 0 and low injection rates through the surface. References [1] T. Ariman, M.A. Turk, N.D. Sylvester, Microcontinuum ¯uid mechanics ± a review, Int. J. Eng. Sci. 11 (1973) 905± 930. [2] T. Ariman, M.A. Turk, N.D. Sylvester, Applications of microcontinuum ¯uid mechanics, Int. J. Eng. Sci. 12 (1974) 273±293. [3] T.C. Chiam, Micropolar ¯uid ¯ow over a stretching sheet, Z. Angew. Math. Mech. 62 (1982) 565±568. [4] L.J. Crane, Heat transfer on continuous solid surfaces, Z. Angew. Math. Phys. 21 (1970) 645±647. [5] A. Desseaux, Numerical contribution to a problem of convection in porous media with lateral mass ¯ux, in: Proc. Ninth Int. Symp. on Transport Penomena in Thermal-Fluids Eng., Singapore, 1996, (pp. 994±999).
N.A. Kelson, A. Desseaux / International Journal of Engineering Science 39 (2001) 1881±1897
1897
[6] A. Desseaux, N.A. Kelson, Flow of a micropolar ¯uid bounded by a stretching sheet, ANZIAM J. 42 (E) (2000) c536±c560. [7] A. Desseaux, N.A. Kelson, Numerical solutions of a micropolar ¯uid ¯ow along a stretching wall (submitted). [8] A.C. Eringen, Theory of micropolar ¯uids, J. Math. Mech. 16 (1966) 1±18. [9] P.S. Gupta, A.S. Gupta, Heat and mass transfer on a stretching sheet with suction and blowing, Can. J. Chem. Eng. 55 (1977) 744±746. [10] F.M. Hady, On the solution of heat transfer to micropolar ¯uid from a non±isothermal stretching sheet with injection, Int. J. Num. Meth. Heat Fluid Flow 6 (1996) 99±104. [11] I.A. Hassanien, Boundary layer ¯ow and heat transfer on a continuous accelerated sheet extruded in an ambient micropolar ¯uid, Int. Comm. Heat Mass Transfer 25 (1998) 571±583. [12] I.A. Hassanien, R.S.R. Gorla, Heat transfer to a micropolar ¯uid from a non±isothermal stretching sheet with suction and blowing, Acta Mech. 84 (1990) 191±199. [13] M.W. Heruska, L.T. Watson, K.K. Sankara, Micropolar ¯ow past a porous stretching sheet, Comput. Fluids 14 (1986) 117±129. [14] C.K. Kang, A.C. Eringen, The eect of microstructure on the rheological properties of blood, Bull. Math. Biol. 38 (1976) 135±158. [15] N.A. Kelson, A. Desseaux, D.L.S. McElwain, Limiting solutions and linearised analysis of micropolar ¯uid ¯ow driven by a porous stretching surface, in: Proc. Seventh Australasian Heat Mass Trans. Conf. (G.B. Brassington, J.C. Patterson eds.) pp. 181±186. Chalkface Press, Western Australia. [16] G. Lukaszewicz, Micropolar Fluids: Theory and Applications, Birkh auser Boston, 1999. [17] J. Peddieson, R.P. McNitt, Boundary-layer theory for a micropolar ¯uid, Recent Adv. Eng. Sci. 5 (1970) 405±426. [18] P.S. Ramachandran, M.N. Mathur, S.K. Ojha, Heat transfer in boundary layer ¯ow of a micropolar ¯uid past a curved surface with suction and injection, Int. J. Eng. Sci. 17 (1979) 625±639.