408

zany1.

Heat Transfer to Magnetohydrodynamic Flow in a Flat Duct By L.

E. ERICKSON, C. S. WANG, C. L. I-IWANG, a n d LIANG-TSENG FAN, Kansas State University, Manhattan, Kansas, USA

1. Introduction

Recently much interest and effort has been focused on the problem of heat transfer in an electrically conducting fluid flowing within a magnetic field. A rectangular duct, which is considered in this work, has applications in devices such as magnetohydrodynamic generators or pumps. SIEGEL [1] 1) investigated heat transfer to the region where the temperature distribution is fully developed and the heat flux at the wall is uniform. ALPRER E2~ investigated the same problem, but assumed that the duct wails were electrically conducting. A case that neglects JouLE's heating in the fluid was investigated b y GERSHUNI and ZKUKHOVITSKII[3] and by REGIRER [4]. YEN [5] investigated the effect of wall electrical conductance for the case with the constant heat flux. The case with the constant heat flux, and with viscous and electrical dissipation, was investigated by PERLS~UTTER and SIEGEL [6i. The case with the constant wall temperature and with viscous and electrical dissipation, which is considered in this paper was previously investigated by NIGAM and SINGt~ [7]. However, the JOULE'S heating term in their paper was wrongly represented [6J, rendering their results invalid. The purpose of this paper is to present the results of the exact solution of the same problem [71, with the use of a finite difference analysis. The fluid properties are assumed to be constant. A magnetic field is imposed perpendicularly to the duct walls, and there is a net electrical current flow parallel to the walls and perpendicular to the flow direction with a variable external resistance connecting the two end plates. The velocity profile is the fully developed H a r t m a n n profile [8]. The temperature profile is uniform at the beginning of the thermal entrance region which is under consideration. Viscous and electrical dissipation effects are not neglected. The solution is obtained by substituting the H a r t m a n n velocity into the energy equation. The finite difference analysis is then employed to obtain the temperature profiles on an IBM 1620 digital computer. The solution, for Graetz number ranging from 10 to 10000, for Prandtl number of 1.0, and for H a r t m a n n number of 0,4, and 10 is presented graphically. 2. Basic E q u a t i o n s

The geometry under consideration, illustrated in Figure 1, consists of two semiinfinite parallel plates extending in the x- and z-directions. The fluid flows in the x-direction; the magnetic field is imposed in the y-direction ; and the electric current flows in the z-direction. Furthermore, the following assumptions are made: 1) Numbers in brackets refer to References, page 418.

Vol. 15, 1964

H e a t T r a n s f e r to M a g n e t o h y d r o d y n a m i c F l o w in a F l a t D u c t

1. The flow is laminar, 2. All the fluid properties, 9, 3. The permeability, #e, and quantities, 4. R a p i d oscillations do nQt negligible, 5. The effect of gravitational

409

@, k, and v are constant, the electrical conductivity, a~, are constant scalar exist and, therefore, the displacement current is force is negligible.

ltlllfl

I

f

~X

f~=const

z

i

ffllftl

i

r

Figure I Parallel p l a t e c h a n n e l w i t h t r a n s v e r s e m a g n e t i c field.

Under these assumptions, the basic equations of m a g n e t o h y d r o d y n a m i c s in emu units m a y be written as follows [91 : c u r l H = 4 Jr J ,

OH

curlE=--#e divJ-

& ,

(1) (2)

0,

(3)

divH = 0.

(4)

O h m ' s law for a moving fluid is

J = ~ (E + W x ~ n ) .

(5)

Continuity equation is div W = 0 .

(6)

The modified Navier-stokes equation is c)W

0r + ( W . grad) W = - - - 1 grad p + v i72 W + 1 (J x / 4 H ) . Q 9

(7)

The fully developed velocity profile used in this work was, as stated previously, originally obtained b y HARTI~(ANN I81. The development is briefly reviewed here. As indicated in Figure 1, the two electrically non-conducting infinite plates are at rest at y = • a. Since it is assumed t h a t the velocity u is in the x direction and the uniform magnetic field H o is imposed parallel to the y axis, it is possible to assume t h a t

w = (~, o, o)

(8)

H = (Hx, Hv, 0)

(9)

and

410

L . E . ERICKSON, C. S. ~vVANG, C. L. H~,VANG, a n d LIANG-TSENG FAN

ZAMP

are functions of y only. Because the fluid near the center line where y -- 0 moves faster than t h a t near the walls where y -- ~ a, it tends to pull out the linear of force in its direction of motion. Thus the field acquires a component H , parallel to the motion, and also a uniform electric field E, which, in general, acts in the z-direction. Using Equations (8) and (9) in Equations (1), (2), (4), (5), and (7) and employing the b o u n d a r y conditions 1. u = O at y= Jza, (10)

Ou

2.

bT=o

one obtains I9] pM

U m

at

[coshM--

y=O,

(11)

cosh(M y/a)]

(12)

The average value of u between y = • a is given b y a

j u dy uo

-~

P

[M c o t h M - 1] .

(13)

Then u = U = M ~ooshM cosh (M yla)] uo [ -)~rcoshM ~ sinhM J '

(14)

The general form of the energy equation is derived b y PAl I10]. For one-dimensional flow of an incompressible fluid with constant properties the energy equation can be simplified to Ot~ (O~t 0 % 10~l~+

J-~

(15)

For steady flow with negligible heat conduction in the fluid flow direction, E q u a t i o n (15) becomes Ot k O~t # J Ou \2 j'~ (16)

u

c; ok

c,

It can be shown [11] t h a t E q u a t i o n (5) simplifies to ]

(17)

GEEo+#
and E 0 is taken as --e u o Bo, where e is the electric field factor and varies between zero and one, corresponding to the external resistance varying from zero to infinity. Since ,ue H o = Bo, Equation (17) becomes f=Uo~eBo

-e+

~;0 '

With this value for J , the energy equation becomes Ot

" ox = Q

k

02t

#

+

/Ou\2

t oy ) +

ug%B~

(

_e

(19)

+

,,o

Vol. 15, 1964

411

Heat Transfer to M a g n e t o h y d r o d y n a m i e Flow in a Flat Duct

Introducing the dimensionless parameters 0-

t-q

p ~ _ #C~ P r a n d t l n u m b e r ,

t o -- t I '

Y-

Y

M = a B o]/-~

~'

Hartmann number

vT

[Z X

X-

k

o~

a2

fl =

'

__ x / a

Uo

Rea'

(2o)

u~

G (to - t,)-'

Equation (19) becomes U O0 OX

1 020 [OU~2+ Pr OY2 + fi \ O Y ] M 2 fi (e -

(21)

u) ~

The term,

Uo~ Cp (~0 -- 11) '

fl

is a criterion of negligibility of the viscous dissipation [121. The viscous dissipation can be neglected when the magnitude of this term is of an order which is smaller than 1.0. When the uniform fluid temperature at the entrance is higher than the wall temperature and the wall temperature is constant, the boundary conditions are 1.

0=1

at

X=0

00 2. 0 - ' 7 = 0

at

Y

3.

at

Y=I

0=0

and

0<

and

0
and

0
0

Y<

1,

(22)

3. Solution of the Energy Equation In order to solve the energy equation, the velocity profile is first determined from Equation (14) and the energy equation is, then, solved by employing a finite difference analysis. The approximating equations are (see the mesh network of Figure 2) Oi, k + I - OL k-1

O0 u=

g~,k,

(0i+1, k+l

0~0 OY a

O0

Oi+l,

ox

k--

Oi, k

~o-?- =

2A v

'

AN

-- 2 0 i + l , k +

0i+1, k - l )

q- (0j, k + l -- 2 0j, k -- 0j, k - l )

(23)

2 (A Y)~ ( u j=.l,k+l

OU

dY

-

q+1,~-1)

2AY

Substituting these expressions into Equation (21), considering the values of 0 with subscript j + 1 as unknowns, and taking k from i to n, one obtains n simultaneous equations with n unknowns ( Os + l, 1; 05 +1, ~ ; . . . ; 0~ +1, n). These linear simultaneous equations are solved by using THOMAS' method [131 and the IBM 1620 computer. The boundary conditions for the difference equations become: 0at

1

X>O,

at

X = 0, and

Y=O,

and

0 ~ Y % 1, 0j+1,~+1=0

O j + l , e = Oj+l,o

at

X>O,

V=l.J

/

(24)

412

L. E. ERICKSON, C. S. \~r

C. L. HWANG, and LIANG-TSENG-FAN

ZAMP

Y-l+n+ln~ / 17-1 R*7 - -

zl g-o ~ l

7

j-7 j j+7

r

X=O Mesh

Figure 2 network for difference representations.

4. Heat Transfer P a r a m e t e r s The bulk temperature (or mixing mean temperature) is evaluated after the temperature, 0 = O(X, Y), is determined. The defining equation 1

f u(Y) o(x, Y) d r

(25)

Oh(X) = 2

f u(v) dv

0

for the bulk temperature in finite difference form becomes at X = (J + 1) d X , /$

~Y'0j+l, k gj+l, ~ aY Oh(X) =

(26)

~=1

U/+I, k A Y k=l

Since

2

Uj+ s, e d Y = 1,

(27)

k=l

Equation (26) becomes o~(x) = ~

% ~ ~, ~ ~ +~,~ ~ Y .

(2s)

k=l

For the case of constant wall temperature the local Nusselt number, N u~, is of secondary importance, but the average Nusselt number, N u m, is of primary importance because it can be used to evaluate the temperature of the fluid leaving the heat exchanger. The average Nusselt number, N u,~, is defined as N u,, -- h~ De k

(29)

Vol. 15, 1964

Heat Transfer to Magnetohydrodynamie Flou Tin a Flat Duct

413

For a duct with spacing of 2 a the equivalent diameter, De, is equal to 4 a. The average heat-transfer coefficient for a duct with length of X and with unit width is, therefore, given b y Q A (~t)

hm -

2 ~(1)Q (~t) "

(30)

The log-mean temperature difference is defined as (Z~r

:

(lO-- tl) -- (tb, x - - tl) l o -- t 1

In

(31)

tb, x -- ll

where to, tl, and tb, x are the temperature at the entrance, temperature of the wall, and the bulk temperature of the fluid respectively. In dimensionless form Equation (31) becomes (32)

1 - - Ob, x 1

(AO)I -

ln--

Ob, x

The total heat flux Q from the entry to x is x

Q = -2/k

ot

(33)

0

Substituting Equations (33) and (30) into Equation (29), the mean Nusselt number is obtained as 4 N um

X(AO)

00 --

0u

dX

,

(34)

0

in which all the terms are dimensionless. The integration is accomplished by using the trapezoidal rule as follows: N U m --

4 X(AO)

Z

01' n+l) -{- (OJ+i' n -

2AY

Oy+l, ~+1) zJ

5. Stability and Convergence The stability and convergence of Equation (21) after Equation (23) is substituted into it, are treated in Reference [14]. It states that such substitution always produces stable conditions, but GRANOY [15] indicated that, as A X and A Y simultaneously and independentIy approach zero in the finite difference equations, the solutions of the difference equations do not necessarily approach to the solution of the differential equations uniformly. It has also been stated that, in m a n y cases, a functional relation must be maintained between A X and A Y as they both approach zero. Thus, it is important to achieve convergence to the true solution of the differential equation within the available computer memory capacity.

414

L . E . ERICKSON, C. S. Waxo, C. L. Hwa.xG, and LIANG-TSENC; ]FAN

ZAMP

The last two terms of E q u a t i o n (21), fl(OU/OY) ~ and M ~ fl (e - U) 2, are constant with respect to temperature. W i t h these two terms omitted, E q u a t i o n (21) becomes U 00 = 1 OX

0~0

(41)

Pr OY ~ "

The stability of E q u a t i o n (41) is practically unaffected b y the addition of the constant terms [1@ The truncation error, resulting from substitution of E q u a t i o n (23) into Equation (41) can be written as [14]

4 0 .... ,) = 0 ( ~ x ) + 0(~Y~).

(4e)

However, if the values of A X and A Y are chosen so t h a t the value of Pr U(AY)2/ 12 (AX) is of a smaller order than 1/2, the truncation error becomes (43)

e[O.... ,] = O ( d X z) + O ( z I Y 4 ) .

In order to obtain the truncation errors of the higher order smallness, the value of P r U(AY)~/12 (AX) is kept less than 0.05. T h o u g h the velocity U is in the range of 0 < U < 1.5, it is taken as 1.0 in calculating the values of P r U(AY)Z/12 (AX). The Prandtl number, Pr, is 1-0. The mesh sizes employed in the calculation are shown in the Table. 6. R e s u l t s a n d D i s c u s s i o n

Figures 3, 4, and 5 show the bulk fluid temperature and wall temperature as functions of the dimensionless distance, X. The wall temperature in dimensionless form is zero. In these figures the parameters are the viscous dissipation factor, fl, and H a r t m a n n number, M. The efficiency of a M H D generator m a y be defined as the ratio of the electrical power to the flow power, which is identical to the value of the electric field factor, e [16]. The value of e for the m a x i m u m power is 0-5 L161. In general a reasonable 74 13 12 ll 79 9

;2

~=7,M~7o

z

,~,.o.s,M:lo

5 4

b'-l, tt~

3 B=~.5,M=4 ] ~ 6

SU/Tfa, Cg ' '

'

B:O.I,i'i~ ' '

'

/J=O.7, H - 4 ' "

'

0.01 L~Oa aDZ 0.04 a05 O.OS 0.97 aa8 ao9 ala o.11 o.12 ~13 0.]4 0#5 x ,X/a ~

Figure 3 Bulk fluid temperature and surface temperature for constant wall temperature,

t),r

=

1.0, d = 0.5.

H e a t T r a n s f e r to M a g n e t o h y d r o d y n a m i c F l o w in a F l a t D u c t

Vol. 15, 1964

415

Z6 Z4 ZB

J

J

J

J

/]=0.5,kidO

I 7o

B=I.0,~/=~ ~5, M=4

I

B =

B=O.ZM=IO B:O.Z,H~4 q~

0.6 0.4 112 surface

0

o.'oJ o.e ~o;-o.o4 eo~ ~;e o.o~ o.;8 a;g ~1o o.;~ o.72x Figure 4

B u l k fluid t e m p e r a t u r e a n d s u r f a c e t e m p e r a t u r e for c o n s t a n t wall t e m p e r a t u r e , P r ~ 1-0, t; = 0.8.

Z0

0.8

/3=7,M-4

O.6

22 ,E

0.4

B=tZS,M=4

~~

as

~-o.s,M.zo

o2 zll

o

--

surface\

B=~ p,~7.M,~o

o.oz aoe ao3 ao4 a)s aoe aOz ao8 aae o.70 a)z aze azs o.}4 0.75 x Figure 5

B u l k fluid t e m p e r a t u r e a n d s u r f a c e t e m p e r a t u r e for c o n s t a n t wall t e m p e r a t u r e , P r

1-0, e = 1.0.

compromise must be made between the conflicting requirement for the maximum efficiency and the maximum power; e = 0.8 is a generally accepted value ~161. Figure 3 shows the case of e ~ 0.5. It indicates that the bulk temperature depends strongly on the viscous dissipation factor,/3, and Hartmann number, M, and increases with the increase of both/~ and M. Figure 4 shows the case of e = 0.8. Comparison of Figure 3 with Figure 4 reveals that, as the electric field factor, e, increases (i. e., corresponding to the decrease of the net current flow), the bulk temperature decreases sharply. Figure 5 shows the case of e = 1, that is, for the case with no net electrical current flow in the channel. In this case, the bulk temperature is only the function of viscous dissipation factor,/~. Hartmann number, M , has little effect on the bulk temperature. The bulk temperature is always less than 1-0, i. e., it is always less than the entrance temperature. Figures 6, 7, and 8 give the average Nusselt number, N~, as a function of the Graetz number with/~ and M as parameters. The Nusselt number is based on the log-mean temperature difference defined by Equation (31). They indicate that, the

416

L E. J2~RICKSON, C. S. WANG, C. L. I['vVANG, and LIANG-TSENG FAN

ZAMP

Nusselt number increases appreciably with the increases of/3 and M, while it increases only slightly with the increase of e. Table

Mesh Sizes for Finite Di//erence Solution o~ the Energy Equation X 0 0.001 0.01 0.1 0.5 2.5

Pr U(A Y)~"

AX

AY

N

} }

0.0005 0.001

0.00625 0.0125

160 80

0-0065 0-013

} } }

0.005 0-01 0.05

0.025 0.05 0.1

40 20 10

0,01 0.02 0.017

12 A X

The Nusselt number for the case with M = 0, fl = 0 is also plotted on Figures 6, 7, and 8. Such a case was treated b y NORRIS and STREID [17] analytically, while it is solved by the finite difference method in this work. The results from both methods show excellent agreement with the m a x i m u m deviation of 3% in the range of Graetz nmnber between 10 and 10000. The results are not compared with those of NIGAN and SINGH [7] because of the reason stated in the introductory section. They should be valid, because the method employed has been proved to yield the correct results known for the case with M = 0, /3 = 0 [17]. I~WANG [18] used the same method to study the forced convection heat transfer in the entrance region where both the velocity and temperature profiles are developing simultaneously under the condition of negligible viscous dissipation. He obtained results which are in excellent agreement with others. 7O

Nu 50

,~

40 l

~

30 ?0 7O

G7 RedPr x/De

Figure 6 Variation of average Nusselt number Nu, with Graetz number for constant wall temperature, Pr = 1.0 and e = 0.5.

Vol. 15, 1964

H e a t T r a n s f e r to M a g n e t o h y d r o d y n a m i c Flow in a Flat D u c t

417

7O

50 ,.

~),

5 "

,~5, # ~ xg

40

80

40

7a ZOO

2'0g ' 4~00 700Z0'00 2000 '4J00' 700~10a00 s~ Figure 7

V a r i a t i o n of a v e r a g e Nusselt n u m b e r N u , w i t h G r a e t z n u m b e r for c o n s t a n t w a l l t e m p e r a t u r e , Pr = 1.0 and e ~ 0.8,

GO

30 ~0 10

7o

2o

' 4~'

'Jo'~oo 2~o ~oo 1oo~ooo

zeee

40ca ze007eeee

Gz Figure 8

Nu, Pr =

V a r i a t i o n of a v e r a g e Nusseit n u m b e r

w i t h G r a e t z n u m b e r for c o n s t a n t wall t e m p e r a t u r e , 1.0 a n d e ~ 1.0.

Acknowledgement The study was supported by the Air Force Office of Scientific Research, Grant AF-AFOSR-463-64. ZAMP 15/27

418

L.E. ERICKSON,C. S. WANG,C. L. HWANG,and LIANG-TSENGF&tr

ZAMP

REFERENCES [I] R, SIXGXL, EHects o~ fffagnetic Field on Forced Convection Heat Trans/er in a Parallel Plate Channel, Trans. ASME, J. Appl. Mech., Sept., 415-416 (1958). [2] R. A. ALPHER, Heat Transfer in ~4agnetohydrodynamic Flow between, Parallel Plates, Int. J. Heat Mass Transfer, Vol. 3, 109 113 (1961). [3] G. E. GERSHUNI and E. M. ZKUKHOVITSKII,Stationary Convective Flow o~ an Electrically Conducting Liquid between Parallel Plates in a Magnetic Field, Soviet Physics, J E T P Vol. 3d (7), No. 3 (September 1958). [4] S. A. REGIRER, On Convective Motion o / a Conducting Fluid between Parallel Vertical Plates in a Magnetic Field, Soviet Physics, JETP, Vol. 37 (10), No. 1 (January 1960). I51 J. T. YEn, Effect o/ Wall Electrical Conductance on Magnetohydrodynamic Heat Trans/er in a Channel, Trans ASME, J. of Heat Transfer, Nov., 371-377 (1963). [6] M. PERLMUTTER and R. SIEGEL, Heat Trans/er to an Electrically Conducting Fluid Flowing in a Channel with a Transverse Magnetic Field, NASA Technical Note D-875 (August 1961). [71 S.D. NI6A~Tand S. M. SINGH, Heat Trans/er by Laminar Flow between Parallel Plates under the Action o~ a 7~ansverse Magnetic Field, Quarterly J. Mech. Appl. Math., Vol. 8, Part 1, 85-96 (1960). [8] J. HAaT~ANN and F. LAZARMS, Kgl. Danske Videnskab. Selskab, Mat-fys. Medd. t5, Vol. 6 (1937). 1_-9] T. C. COWLING, Magnetohydrodynamics, New York: Interscience Publications, Inc., 2-15 (1957). [i01 S. I. PAl, Energy Equation o/ Magneto-Gas Dynamics, Physical Review, Vol. 105, No. 5, 1424-1426 (March 1957). [11] M. ROIDT and R. D. CEss, A n Approximate Analysis o/Laminar Magnetohydrodynamie Flow in the Entrance Region o/ a Flat Duct, Trans. ASME, J. Appl. Mech., March, 171--176 (1962). [121 E . R . G . ECKERT and R. M. DRAKt~, Jr., Heat and Mass Trans/er, New York: McGraw-Hill Book Company, Inc., 2nd ed., 171-173 (1959). I131 L. LAPlDUS, Digital Computation/or Chemical Engineering, New York: McGraw-Hill Book Company, New York, 254-255 (1962). I141 R. D. RICHTMigYER,Di//erence Methods/or Initial Value Problems, New York: Interscience Publishers, 91-101 (1957). [151 R. A. GRAND'Z, The Aeronautronic Hop Program/or Fluid Flow, AFWC-TN-61-29, Part 1, 88-95, AFWC Second Hydrodynamic Conference Numerical Methods of Fluid Flow Problems (May 1961). [161 W. C. ]V[OFFATT, Magnetohydrodynamics Power Generation - Its Principles and Problems, Tech. 2~eport MIT-29-P, Project Squid, Dept. of Aerospace Engr., University of Virginia (Jan. 1963). [171 R. H. NORRIS and P. P. STREIO, Laminar-Flow Heat-Trans/er Coe//icients/or Duets, Trans. ASME, Vol. 62, No. 6, 525-533 (August 1940). [18] C.L. HWANG, A Finite Di//erence Analysis o/Magnetohydrodynamic Flow with Forced Convection Heat Trans/er in the Entrance Region o/ a Fiat Rectangular Duct, P h . D . Thesis, Kansas State University (1962). Zusammen/assung ]Die Untersuchung befasst sich m i t d e r Wgrmeiibertragung yon den Wgnden eines flachwandigen Kanals auf eine elektriseh leitende Fliissigkeit bei erzwungener Laminarstr6mung und in Gegenwart eines quergerichteten Magnetfeldes. Betrachtet wird der Fall konstanter Wandtemperatur mit variierender innerer Wgrmeentwicklung dutch viskose nnd elektrische Energiedissipation. Die massgebende Differentialgleichung wird durch eine Differenzengleichung ersetzt und mit der elektronischen Rechenmaschine gel6st. Als Resultat wird die Nusseltzahl angegeben, fiir die Prandtlzahl 1, die Hartmannzahlen 0, 4, 10 und die Graetzzahlen von 10 bis 10 000, wobei die Kennzahlen fiir ZS.higkeit und elektrische Feldstgrke als Parameter auftreten. (Received : December "20, 1963.)