An analytical solution to fuel-and-cladding model of the rewetting of a nuclear fuel rod

Nuclear Engineering and Design 61 (1980) 101-112 O North-Holland Publishing Company

AN A N A L Y T I C A L S O L U T I O N T O FUEL-AND-CLADDING M O D E L OF T H E R E W E T T I N G OF A N U C L E A R F U E L ROD Hsu-Chieh Y E H Westinghouse Electric Corporation, Nuclear Energy Systems, P.O. Box 355, Pittsburgh, PA 15230, USA

Received 21 April 1980 An exact solution of the quasi-steady two-dimensional conduction equation for the rewetting of a nuclear fuel rod in water reactor emergency core cooling is obtained for a fuel-and-cladding model. A method of solving non-separable differential equations is presented, which is used in the present analysis. The recently developed theorem of orthogonality of piecewise continuous eigenfunctions is also used to handle the composite rod in the present model. The present analysis reveals that the wet front velocity increases with the increase of the gap resistance between the fuel and the cladding, and approaches a limiting value, which is equal to the wet front velocity of the tube of cladding alone, as the gap resistance becomes infinite. For convenience in practical application, the results of the present analysis are correlated in simple expressions.

indrical rod model. T h e exact solution of the two-dimensional cylindrical rod model was obtained by Yeh [13] using his method of solving heat conduction for complicated geometries [14, 15]. This is a powerful mathematical tool, as illustrated by the solution of a n u m b e r of difficult problems [13-16]. A survey of the literature on rewetting was presented by McAssey and Bonilla [17]. A m o r e realistic model for the rewetting of nuclear fuel would consist of inner fuel and outer cladding, with gap resistance between the fuel and the cladding. Unfortunately the analysis of the fuel-and-cladding model is very difficult and only an approximate solution has been obtained [18]. This p a p e r extends the work of ref. [13] to analyze the two-dimensional cylindrical fueland-cladding model by using a recently developed theorem [19, 20] and the method of solving non-separable differential equations described in this paper. O n e of the complications in solving the rewetting problem of the fuel-and-cladding model is due to the fact that the fuel rod is a composite rod having gap resistance between the fuel and the cladding. In the general boun-

1. Introduction

During a postulated loss-of-coolant-accident ( L O C A ) of a nuclear reactor, the overheated fuel rods are cooled by the emergency cooling water which enters the core by either top spray or b o t t o m flooding. The rod surface cannot wet until the rod t e m p e r a t u r e is brought down to the Leidenfrost temperature. The prediction of the velocity of advance of the wet front is important in evaluating the performance of the emergency core cooling system in a postulated LOCA. Most studies on the rewetting problem are for the one-dimensional conduction model [1--6]. Although the results of these analyses predict the wet front velocity quite well for low coolant flow rates, it is felt that the two-dimensional conduction model should be considered for high coolant flow rates. Because of mathematical difficulty most two-dimensional analyses are either numerical or approximate ones. These include T h o m p s o n ' s numerical solution [6, 7] and asymptotic expression [8], approximate analyses of Duffey and Porthouse [9], C o n e y [10] and Tien and Yao [11] for a slab model, and Blair's approximate solution [12] for a cyl101

102

H.C. Yeh / Fuel-and-cladding model of rewetting of a nuclear fuel rod

dary value problems with composite media, the coefficient of the differential equation and/or the derivative of the dependent variable are discontinuous. Therefore, the Sturm-Liouvelle theorem does not apply and separation of variables will not yield the desired orthogonal eigenfunctions. The theorems of references 19 and 20 can be used to obtain the orthogonal eigenfunctions for the boundary value problems with composite media in two ways: (a) By transforming the dependent variable with an undetermined constant and determining the undetermined constant by the theorems of references 19 and 20 such that the solutions of the new dependent variable form a complete orthogonal set of eigenfunctions. (b) By using these theorems [19, 20] to obtain a weighting function such that the eigensolutions are orthogonal with respect to the weighting function. Since method (a) has been illustrated in references 19 and 20, this paper demonstrates the application of the weighted-orthogonality method (b). For convenience in practical application, the analytical results of the uniform cylindrical rod model [13] and the present fuel-and-cladding model are correlated in simple expressions in the present paper.

2. Formulation of the problem

Consider a cylindrical fuel rod consisting of the inner uranium dioxide pellets and the outer cladding. For simplicity the gap between two pellets is ignored and the pellets are considered as an uniform cylindrical rod, which will be called the 'fuel' (fig. 1). For mathematical simplification, the inside diameter of the cladding and the outside diameter of the fuel are assumed to be equal, but there is a finite gap conductance, hgap (or gap resistance 1/hgap), between cladding and fuel. It is assumed that the quasi-steady state exists and there are two distinct heat transfer regimes on the outer surface of the cladding: a uniform high heat transfer

REGION A (WET) : REGION B ('DRY) :

Z <.0 z > 0

I a,

AI

,

u

J

_ _ .

I o

, [

GLADDING U

A

I

B

P I

_.J ta]

(b)

(=)

Fig. 1. Physical model: (a) bottom flooding, (b) spray cooling and (c) rod cross-section.

coefficient over the wet surface and a negligible heat transfer coefficient over the dry surface. The origin of the coordinate system, which is moving with the wet front, is on the rod axis at the wet front elevation. As in ref. [13], the whole rod (fuel plus cladding) is divided into two regions, A and B, with the cross section at the wet front as the dividing surface (fig. 1). Region A is on the wet side and region B is on the dry side of the wet front. The quasi-steady conduction equation in the coordinate system moving along with the wet front at the constant velocity u, and with the conventional assumption of negligible internal heat generation, is given by

(-3T)+O2T PiC~uOT=o, rWr -b-U 4 ae i = 1,2,

(1)

where the subscript i = 1 denotes the quantities for the fuel and i = 2, for the cladding. Eq. (1) is valid for both regions A and B.

103

H.C. Yeh / Fuel-and-cladding model of rewetting of a nuclear [uel rod

The external boundary conditions for region A are

In terms of the new variables, conditions (2) through (6) become

TA(L --oo) = finite,

OA(r, --~) = finite,

TA(0, ~) = finite,

(2a,b)

kE(OTA/Of)¢=, 2 = -h[Ta(~2, 2 ) - %1

(2c)

TB(L 0o) = T=,

TB(0, 2) = finite,

(3c)

The internal boundary conditions are (4a)

k1(OTIO0,=,,- = k2(OTIO0,=,,+,

(4b)

which apply to both regions A and B and the subscripts A and B on T are understood and are omitted. On the dividing surface between regions A and B, the temperature and the heat flux must be continuous, i.e.,

(OTA'~

\0~/~=0

Os(O, z) = finite,

O(r l , z ) = O(r~f , z kl

hgap[T(~{, ;~)- T(~i~, ;01 = -kz(OT/OF),=,~+,

TA(L O)= T.(?, 0),

(10c)

(3a,b)

(OTa/O?)e=,2 = O.

= (OTa'} \O;~/~=0"

(10a,b)

(aOA/ Or)r=r 2 = --BOA(r2, z ),

On(r, ~) = 1,

and those for region B are

0A(0, Z) ----finite,

) --

{00s] = O. \ Or It=r2 (lla,b,c)

/ 00 x Rgap~-~ ) r=rl+,

(12a)

O0

(12b)

\O-r],= r,+ 0A(I, 0) = 0s(r, 0),

(00A~ = (00B~ , \3z/~=o \Oz/~=o

(13a,b)

0B(r2, 0) -----0c

0A(r2, 0) = 0c,

(14a,b)

or

where B is the Biot number, 0¢ is the dimensionless Leidenfrost temperature and Rgap is the dimensionless gap resistance defined as B = hrE/k2,

Oc =

( T ¢ - Tf) (T~ - Tt)'

k2 Rgap = (hgap~2).

(15)

(5a,b) Besides the above conditions the temperature on the cladding surface at the wet front must be equal to the Leidenfrost temperature Tc, TB(~z, 0) = Tc

or

Tg(f2, 0) = Tc.

(6a,b)

The above expressions can be put in dimensionless form by introducing the following variables O(r, z ) - ( T - Tf)/(T~ - Tr), r -----~/~2,

(7a)

z - z/r2,

(7b,c)

where the subscripts A and B on 0 and T are understood and are omitted. With definitions (7), eq. (1) becomes 1 0 (r 00~ 020 , ~ , \ 00 r Or\ Or/+~ +u (r)~z=O,

(8)

in which the dimensionless wet front velocity u*(r) is a step function of r in either region A or B as defined by u (r), _

{* ul

= plclug2/kl, u~ = pzc2uFz/k2,

0 <- r < rl, rl < r <- r2.

(9)

Eqs. (12a,b) are valid for both regions A and B, and the subscripts A and B on 0 are understood and omitted. Because u*(r) in eq. (8) is a function of r, a direct separation of variables is impossible. In the following a method is presented to resolve this difficulty.

3. Method of solving non-separable differential equations The linear differential equations, which is non-separable due to the variable coefficient such as u*(r) in eq. (8), can be solved by the following method. Step 1. Obtain eigenfunctions R, (r) by solving the equation which can be obtained by separating the variables of the original differential equation, eq. (8), with u*(r) being temporarily assumed to be constant. Step 2. Express the solution in the series of eigenfunctions R, (r). For the present case this

104

H.C. Yeh / Fuel-and-cladding model of rewetting of a nuclear fuel rod

means

kt

k-~R ~(r?) = R ~,(ri~),

0 = ~ R,(r)Z,(z). n=l

(16)

Step 3. Substitute eq. (16) in the original differential equation, eq. (8), with u*(r) being a function of r. The resultant equation will contain terms of the product of the variable coefficient and the eigenfunction, u*(r)R,(r). Step 4. Expand the product u*(r)R,(r) in the equation obtained by step 3 in the series of eigenfunctions. The resultant equation can be reduced to a form which has zero on the one side of the equal sign, and on the other side has a series of eigenfunctions with coefficients being functions of Z~(z) and its derivatives. Since the eigenfunctions are functions of r, while the coefficients are functions of z, the equation can be satisfied only if each coefficient is zero, which leads to a system of simultaneous equations for

Z.(z). Step 5. Solve the system of simultaneous equations obtained from step 4 for Z.(z) with appropriate boundary conditions.

(19b)

where the prime denotes the derivative with respect to the argument of the function. The solution of eq. (17) can be expressed by = ~CA1J0(Ar)+ DA1 Y0(Ar), 0 <--r < rl,

RA

[CA2Jo(Ar)+D~.Yo(Ar),

To satisfy condition (18a) DA~ must be zero. The other constants can be determined from conditions (18b) and (19a,b) which yield a t 2 C a 2 4" al3DA2 = 0,

(21)

a21CA1 + a 2 2 C A 2 a31CA1 +

+

a23Dm = 0,

(22)

a32Cg2+

aa3DA2 = 0,

(23)

where a:2 = )tJ~(Ar2) + BJo(Ar2),

(24a)

a~3 = AY~(Ar2) + BYo(Ar2),

(24b)

azl = J0(Arl),

a22 = -Jo(Ar0 + RgapAJ~(ArO (25a,b)

a23 = - Yo(Arl) + RgapAY/,(Arl),

a3~ = (kdk2)JD(Ar~), 4. Solution

d (r dRA~ + A2rR A = 0,

-~r\

dr]

(17)

which is obtained by separating the variables of eq. (8) with u*(r) being temporarily assumed to be constant. The boundary conditions required for solving eq. (17) can be obtained from conditions (10b,c) and (12a,b) by using eq. (16), giving RA(0) = finite,

a32 = -J~)(Arl),

a33 = - Y~(Ar~).

The problem formulated above can be solved by the method just described, the theorem of ref. [20], and the method of refs. [14] and [15]. Consider region A. According to step 1 of the method described above, the eigenfunction Ra(r) can be obtained by solving the following equation

R'A(r2)= -BRA(r2),

RA(ri-) = RA(rT) -- RgapR k(ri~),

(18a,b) (19a)

(20)

rl
(25c) (26a,b) (26c)

Eqs. (21}-{23) are homogeneous simultaneous equations for CAI, CA2, and DA2. Non-trivial solutions exist if the determinant of the coefficients is zero, i.e., A --- a21a32a13 + a31al2a23 - a13a22a31 - a33a12a~l = 0.

(27)

Eq. (27) can be solved for the eigenvalue A. For each value of A, A,, which satisfies eq. (27), only two of the three equations (21)--(23) are independent, and two of the unknowns, say CA1 and DA2 can be solved from these two equations in terms of the other, say CA2. The results are

CAI = fA,CA2,

DA2 = gA,CA2,

(28)

where fA, = (a12az3 -- a,3az2)/a13a21,

(29)

gA, = -a12/aj3.

(30)

H.C. Yeh / Fuel-and-cladding model of rewetting of a nuclear fuel rod

The function RA., aside from the constant multiplier CA~, is

105

B,iB2, - B3,B4i = P(X~f ) b(xF)

p(x?) b(x~)' (35)

i=1,2 ..... N-1.

[fAjo(A,r),

(31)

RA,(r) = [J0(A,r) + gA, Y0(A,r), where A, is the positive roots of eq. (27). The functions {RA,(r)} do not form an orthogonal set, because the function RA,(r) and its first derivative are discontinuous at r = r~ as indicated by eqs. (19a) and (19b) and, therefore, the Sturm-Liuovelle theorem of orthogonality does not apply. The functions, however, can be made orthogonal to each other with respect to a proper choice of weighting function which can be found by using the theorem of ref. [20]. This theorem is as follows.

The proof of the theorem can be found in ref. [20]. To apply the theorem, re-write eq. (17) as

~-r (FAir~--r ) + AFAirR = 0 , rH
i = 1,2,

(36)

where r0 -= 0 and FAi is constant within the interval rH < r < r~ yet unknown. It is desired to determine the unknown FAi from the above theorem such that the solution RA,(r) we have obtained be orthogonal with respect to the weighting function FAir. That is

Theorem. Given the differential equation

~---~[p(x)dd-~x]+[q(x)+ Aw(x)]y = 0

(32)

i=l

f Rg.(r)RA,.(r)FAirdr=O,

if n e r o ,

ri-I

where p(x), q(x) and w(x) may be discontinuous at x = xl, x2. . . . . xN-~ in the closed interval xo
(33a)

C l y ( x N ) - C2y'(xN) = 0,

(33b)

where the prime denotes the derivative and A l, A2, C1 and C2 are arbitrary constants, and these solutions possess the following discontinuities at ~

Xl,

X2,

. . . ~ Xi,

• • • ~ XN-I:

BN = 1,

B21 = k2/kb

p = FAir,

+

(34a)

[b(x)yl'=x; = Bzi[b(x)y]'=;~ + B4iy(xT), i = 1, 2 . . . . . N - 1,

(34b)

in which Bli and B2i are the constants and b(x) and b'(x) may be discontinuous at x = xl, x2. . . . . XN-I, and if Yl, Y2, Y3, -.. are the solutions corresponding to these values of A, then the functions {y,(x)} form a system orthogonai with respect to the weight function w(x) over the interval (Xo, xN) if

(37)

b = 1,

(38a,b,c) (38d,e)

w = FAir,

for ri-1 < r < ri,

i = 1, 2.

(39a,b)

Substituting eqs. (38)-(39) in the condition of orthogonality (35) yields FAI/FA2 =

t

y (x ,) = B liy (x ~) + B3, [b (x)y ] ~=~,,

ifn=m.

To determine FAi, compare eq. (36) with eq. (32) and eqs. (19a,b) with eqs. (34a,b), we have the following relations

B31 = -Rgap, B41 = 0,

and

X

rs0,

kl/k2.

(40)

Eq. (40) implies that the ratio of FA1 to FA2 must be equal to kdk2. The proportional constant is immaterial as it does not affect the orthogonal relation (37). Fgi will be taken as the dimensionless quantity: FA1 = kl/k2,

FA2 = 1.

(42a,b)

Now, express 0A in the series of the eigenfunction R A n :

Og = f~, RA,(r)ZA,(Z) n=l

(43)

106

H.C. Yeh / Fuel-and-claddingmodel o[ rewettingof a nuclearfuel rod

Substitute eq. (43) in eq. (8), giving

£ 1-d(rdRg"'~ ,=lrdr\

dr ] Z A " + ~nR= lA "

+ .=1 ~, u*(r)Rg, dZA, dz

where A , and s are constants to be determined. Substituting eq. (50) in eq. (49) gives

d2Zg, dz 2

N

s2A. + ~ UAi,SAi - AZ,A, = 0, i-1

= O.

(44)

According to step 4 of the above method, the product u*(r)RA,(r) is expressed in the series of eigenfunctions, giving

u*(r)RA,(r) = ~_~ UA,iRAi(r),

(45)

i-1

where UA,~ are the constants defined as

n = 1,2 . . . . . N.

(51)

Eqs. (51) are simultaneous linear equations for the unknowns A.. Since they are homogeneous equations, non-trivial solutions exist only if the determinant of the coefficients is zero. That is

(s2 + UA.S -- * ~) U.2~s UA,~s (s 2 + UAa2S - ~ ~) U~3s U~23s

.oo

.o.

=0.

..°

UA.i = ~',..,

(52)

U*(r)RA.(r)RAi(r)Fa~r dr

]=1 rjt

/~_l f R2i(r)Fajrdr.

(46)

rjl

Substituting eq. (45) in eq. (44) and utilizing eq. (17), one obtains

--

A nRAnZAn + n-I

RAnZ;n n=l

(47)

"l- £ 2 ZrAnUAniRAi =0. n = l i=1

Changing the dummy indices n and i, n ~ i and i-~ n, in the last term, eq. (47) becomes

A, ZA. + n=l

)

UA,.Z'A, -- A2.ZA. = 0. i=1

(48)

Since, in general, RA. is not zero, we must have

Z';.. + ~ , U.~.Zk~ - A , Z A . = O , '='

n=1,2 ..... (49)

Eqs. (49) represents infinite many differential equations which can be solved for infinite many variables ZA.. AS common practice in all series solutions, numerical calculation is carried out by taking a finite number of terms only. Similarly, in solving eqs. (49), only finite number of equations and ZA,, say n = 1, 2 . . . . . N, are considered. Since eqs. (49) are linear, the solutions can be written as

ZAn (Z) = A. e ~,

(50)

In general eq. (52) can be solved to obtain N positive roots s,. For each root, say s = si, only N - 1 equations in eq. (51) are independent, which can be solved for the N - 1 unknown A, in terms of the remaining one, say Aj, giving A , = a,y4j,

n = 1, 2 . . . . . N,

(53)

where a# = 1. Now, ZA,(Z) is the sum of all solutions, that is N

ZAn(Z) = Z a.;Ai e'& j=l

(54)

For a finite number of terms N, solution (43) can be written 0A =

RA.(r)

a,e% e s~ •

(55)

n=l

The solution for region B can be obtained in the similar manner. The resultant expression for 0s is

OB= 1 + ~. RB.(r)

b.jBj e-'~ ,

(56)

where

f,fsjo(/3,r),

R s . ( r ) = tJ003,r)+ ga, Yo(I3.r).

(57)

The eigenvalues/3, are the positive roots of an equation having the same expression as eq. (27) and fs. and gs. are given by eqs. (29) and (30),

H.C. Yeh / Fuel-and-cladding model of rewetting of a nuclear fuel rod

107

respectively, in which the a's are given by eqs. (24) through (26) with A replaced by/3 and B = 0. The orthogonal relation similar to eq. (37) also holds for region B. The weighting function for both regions A and B are the same and subscripts A and B on F factor will be omitted. That is

Multiplying eq. (61) by the eigenfunction of region A, RA.(r), and weighting function Fir and integrating it over the interval (0, rE), gives

Fj={FI=FAI=Fm=-~z,

0--
NAm ~ a~iAi = Pm +

F2 = FA2 = FB2 = 1.

rl < r --< r>

sianiAi

n=l

tibniBi RB.(r).

RAn(r)=-n=l

(62)

)

bniBi Qm., n=l

(58)

m = 1, 2 . . . . . N,

The constants tj in eq. (56) are the positive roots of an equation similar to eq. (52) with s, UA~n and A replaced by t, -UB~n and/3, respectively, where UB~, are given by

(63)

where q

NAm ~ ~l f FjrR2Am(r)dr'

(64a)

rj-I

rj 2

Uain=j~_1 f u*(r)Ran(r)RBi(r)Firdr Prn

rii

=--j~_~ f FjrRA~(r)dr, =l

(64b)

rj-i rj

rj

rRA (r)Ran(r)dr. ri-I

(64c)

rj-1

The constants b.j are obtained in the same manner as for a~. That is for each root of t, say t = ti, solve any N - 1 equations of the following N equations

Similarly, multiplying eq. (62) by ~rRBm(r) and integrating it over the interval (0, r2) gives

sianiAi O~,. =--NBm

n=l

"=

tibmiBi,

N

t2Bn - ~ Uai.tBi -/32~B.

= O,

m = 1, 2 . . . . . N,

(65)

where n = 1, 2 . . . . . N,

(59)

2f

for N - 1 unknown B, in terms of the remaining one, say B~, giving

Na~ - ~.

B, = b,~Bj, n = 1, 2 . . . . , N,

Eqs. (63) and (65) represent 2N equations which can be solved simultaneously for 2N unknown coefficients Ai and Bi. Condition (14) can be written as

(60)

where b# = 1. The coefficients A, a n d / 3 , in eqs. (55) and (56), respectively, can be determined by the matching conditions (13a,b) at the dividing surface between regions A and B using the orthogonal properties of eigenfunctions [14, 15]. Substitute eqs. (55) and (56) into eqs. (13a,b) giving

n=l

=

an.,Ai RAn(r) = 1 + .-. n=l

)

b.~Bi

Ran(r), (61)

jTM =l

FirR 2~(r) dr.

(66)

rjl

1 + ~'=t = b.,B,

RB.(r2) = 0¢,

(67a)

or N n=l

The computational procedure is as follows. First, a value of u is selected, and A. and B. are

108

H.C. Yeh / Fuel-and-cladding model of rewetting of a nuclear fuel rod

determined from eq. (27). Then eqs. (63) and (65) are solved simultaneously for Ai and Bi. Condition (67a) or (67b) serves as a check of the correctness of the assumed value of u. If condition (67a) or (67b) is not satisfied, a new value of u is selected and the calculation is repeated. It is noted that the above solution becomes exact as the number of terms N in the series expansion approaches infinity.

5. Correlation of the analytical results Since the calculation of the above analysis is quite complicated, for convenience in practical application, the results of the above analysis and the previous work [13] for an uniform cylindrical rod model are correlated as follows: u**

=-

( u ~ - u2..~)/(u2.tube - u2.~,)

= ug* + (1 - u,~*) x [1 -

e -1"464b - - 0 . 5 ( U ~ * / U f f I * )

1"322

× e -2°7b sin(3.429b)],

(68)

b = RgapBF(O~),

(69)

F(0c) = 0.16710c - 2.54702 + 14.460~ - 27.5880~ + 18.720~, u~* = u81"

1

(70)

r /klx m

, +0.15

(71) u'~2 = u '~/u~ = plctk2/(p2CEkO,

(72a)

u~l* = -Uz.~t/(u2.tube- Uz.cyO.

(72b)

The parameter UExy~is the dimensionless wet front velocity of an uniform cylindrical rod of ref. [13] having the same material properties as the cladding of the present fuel-and-cladding model. The results of ref. [13] (fig. 1 of ref. [13]) can be correlated, in terms of the cladding material properties of the present model, as follows:

where E = 0c[2/(1 - 0c)]1/2,

(74)

and B and 0¢ have the same definitions as above. Eq. (73) reduces to a one-dimensional solution [13] when B <~ 1: * = U 2,cyl

O~[2B/(1

-

Oc)]1/2.

(75)

The parameter u 2,tube*is the dimensionless wet front velocity of a tube of cladding alone with the absence of the inner fuel. The solution for u 2,tubecan be obtained by following the analysis of ref. [13] with eigenfunction J0(A.r) replaced by J0(A,r) + g A n Y0(A,r) for wet region A and Jo([3,,r) replaced by Jo([3,,r)+ gB,,Yo(!3,,r) for dry region B. In ref. [13] it is noted that for small Biot number B the two-dimensional solution for the wet front velocity of an uniform cylindrical rod coincides with Yamanouchi's [1] onedimensional solution for a tube if the thickness of tube wall 72- ?~ is replaced by half of the radius of the cylindrical rod ~2/2 (in terms of the present notations). For large B the two-dimensional solution deviates considerably from the one-dimensional solution. Thus one expects that the two-dimensional tube solution for small Biot number B should coincide with Yamanouchi's one-dimensional tube solution. In other words, one expects that fig. 1 of ref. [13] with the radius r2 (= ~0 of ref. [13]) in the definition of u* and B replaced by two times the tube wall thickness 2(~2- ?l) is valid for a tube for small Biot number B. This is indeed true. Now, one may ask: "Is fig. 1 of ref. [13] for large Biot number B also valid for a tube if ~2 in the definitions of u* and B is replaced by 2(72f 0 ? " Fortunately, this is also true. Therefore, eq. (73), which is a mathematical representation of fig. 1 of ref. [13], can also be used to evaluate U2,tube*with B replaced by 2(F2- r t ) B / r 2 and U2,cy I* by 2(f2- ?l)U~,t,be/~2. That is: , U2,tube =

r2

0.6003 2(~2--- r l ) x exp[ l

~ 1.4443 l~ 0.75 a_, *"t

{log(O.1299BtE°4443)}2 + O.05] 1/2,

U 2,cyl* ---- 0.6003EL~XB °75

(76) [ 1 - {log(0.1299BE04~3)} 2 + 0.05]l ' a , x e x p t16

(73)

where

H . C Yeh / Fuel-and-cladding model of rewetting of a nuclear fuel rod

properties of the fuel being the same as those of * is the special case of the the cladding, and u2,t,~ present model with Rgap = ~, in which the cladding is thermally isolated from the fuel and the physics of rewetting is identical to the case of a tube. These special cases of the present model have been checked with good agreement with the solution for a tube and the results of ref. [13] for an uniform cylindrical rod. Fig. 2 plots u** against Rgap for 0c = 0.3, B = 1, u?2 = 1, kl/k2 = 0.3, 1, 3, and 10. The special case of a uniform cylindrical rod of ref. [13] corresponds to u** = 0, and the special case of a hollow tube of cladding alone corresponds to u** = 1. It is seen from fig. 2 that all curves approach the asymptotic value of u** = 1 as Rgapbecomes large. This is expected as mentioned above, since for this limiting case the cladding is thermally isolated from the inner fuel. Fig. 3 shows the dimensionless temperature distribution for the dimensionless gap resistance of Rgap = 1.103, U~2 = 4.938, 0c = 0.5, rl = 0.885, B -- 10, and k d k : = 0.2784, and fig. 4 shows the dimensionless temperature distribution for the

I .6 .6 \

.4 .2 0 -.2

--.4 -.6 I

-.8

I

I

I

J

I

J

I

I

I00

0

J

200

I

I

I

I

J

I

I

;300

[

I 4.00

Rgap Fig. 2. D i m e n s i o n l e s s w e t f r o n t v e l o c i t y u** v e r s u s d i m e n s i o n l e s s g a p r e s i s t a n c e R o p for 0c = 0.3, B = 0.1, u 1"2= 1 a n d rl = 0.885.

Bt

2(~2- ~1) B - 2(~z - ~a)h ¢2

--

109

(77)

k2

It is noted that u ~.cy~is the special case of the present model with Rgap = 0 and the material

!i ! Z

Fig. 3. D i m e n s i o n l e s s t e m p e r a t u r e 0 as a f u n c t i o n of d i m e n s i o n l e s s a x i a l d i s t a n c e z for R ~ p = 1.103, ul*2 = 4.938, 0c = 0.5, rl = 0.885, B = 10 a n d kdk2 = 0.2784.

H.C. Yeh / Fuel-and-cladding model of rewetting of a nuclear fuel rod

1 I0

I

I

1

I

I

.9 .8 .7 .6

e .5 .4

.3

~CLADDING

.2

&r _ Fz _

q

.I 5

I

0 ~ -5

-4

-3

I

I

-2

-I

J

I

0

2

Z

F i g . 4. D i m e n s i o n l e s s B = 10, a n d

kl/k2

temperature

0 as a function of dimensionless

axial distance z for

Rgap =

0, u 12 = 4 . 9 3 8 , 0c = 0.5, r~ = 0 . 8 8 5 ,

= 0.2784.

same parameters as in fig. 3 except that Rgap = 0. Twelve terms have been used for these calculations. The fact that the temperature at z = 0 joins quite smoothly indicates that the series is convergent. In general, the number of terms required increases as Blot number B increases. In both figs. 3 and 4 the temperature domain of the fuel and that of the cladding are shaded for clarity. It is seen from fig. 3 that the gap resistance causes the temperature domains of the fuel and the cladding to split apart on the wet side and the fuel temperature to be higher than that for the case with zero gap resistance of fig. 4. In the limiting case of Rgap = ~, the fuel temperature would remain at T~ (i.e., 0 = 1), because the stored energy of the fuel cannot transfer out through the gap between the fuel and the cladding. From fig. 3 it is seen that the surface temperature of the cladding drops abruptly from practically equal to T= to practically equal to tf in a short distance from the wet front. Thus, although T= is mathematically the temperature at infinity, it is practically the temperature of the rod right before the wet front arrives. T~ is

commonly called the 'quench temperature'. In the analysis of a loss-of-coolant-accident, the fuel rod temperature on the dry side and hence T~ are computed with the internal heat generation and heat transfer coefficient. Therefore, although the heat transfer coefficient on the dry side and the internal heat generation are neglected in the present analysis, they have an indirect effect on the wet front velocity through the computation of T~. Physically it means that the surface heat flux and the internal heat generation have the effect of building up the stored energy in the fuel rod on the dry side, while the rewetting velocity depends on how fast the stored energy can be removed. This in turn depends on how much the stored energy has been built up by the action of the internal heat generation and the surface heat flux.

6. Concluding remarks The exact solution of the rewetting problem for the fuel-and-cladding model is made possible by the method of refs. [14] and [15] of solving

H.C Yeh / Fuel-and-cladding model of rewetting of a nuclear fuel rod

potential field problems for complicated geometries, the theorem of refs. [19] and [20] of the orthogonality of piecewise continuous eigenfunctions, and the method of solving nonseparable equations described in the present paper. The present solution is valid over the whole range of Biot number. The effect of the thermal resistance of the gap between fuel and cladding is to increase the wet front velocity. The wet front velocity approaches a limiting value, which is equal to the wet front velocity for a hollow tube of cladding alone, as the gap resistance increases indefinitely. For convenience in practical application, the results of the present analysis and that of ref. [13] for a uniform cylindrical rod are correlated in simple expressions. In the present paper fuel plus cladding are considered as one single region rather than two separate regions, as commonly considered in treating problems with composite media. Thus the whole rod consists of only two regions, A and B, rather than four regions. The eigenfunctions in each region are piecewise continuous. This problem of piecewise continuous eigenfunctions is treated elegantly and systematically by the theorem of orthogonality of piecewise continuous eigenfunctions [20]. It is interesting to note that the method of weighted-orthogonality and the method of refs. [19] and [20] of transformation of the dependent variable are equivalent. This can be seen from the orthogonality condition (37). The orthogonality condition (37) can be written as

The method of refs. [20] and [19] consists of changing the dependent variable by X, = G j R A . , and determining Gj by using the orthogonality theorem of piecewise continuous eigenfunctions. It can be shown that Gj thus obtained is exactly equal to F)/2 and the orthogonality condition is given by eq. (80). Thus the two methods are equivalent. However, with the method of weighted-orthogonality the solution is expanded in the series of eigenfunctions {RA,}, while with the method of transformation of the dependent variable, the solution is expanded in the series of eigenfunctions {X,}, which results in a slightly different expression.

Nomenclature

A, a12, a13, . . . , etc.

B. B

CA,, C~2 cl DA1, DA2

FA,, F.i, F/

fAn gA. h

hgap L ki

NAm, NBm

,~ i (F~ERA"XF~nRAm)rdr=O'if n#m, ,j-i

~ 0,

if n = m. (78)

1,.

By the transformation of the dependent variable Xn = F}/2RA.

(79)

the orthogonality condition (78) becomes

~ 'J-I

X.Xmrdr = 0,

if n ~ m,

#0,

if n = m .

r

(80)

111

to, r l , r2

coefficient of series expansion for 0A defined in eqs. (24)-(26) coefficient of series expansion for 0B Biot number, coefficients in eq. (20) specific heat coefficients in eq. (20) weighting function, eqs. (37), and (58) coefficient in eq. (31) coefficient in eq. (31) heat transfer coefficient gap conductance Bessel function of first kind of order n thermal conductivity defined in eqs. (64a) and (66), respectively defined in eq. (64b) defined in eq. (64c) radial coordinate 0 outside radius of fuel or inside radius of cladding outside radius of cladding dimensionless radial coordinate, ~/~2 0, ~dF2, and 1, respectively

h~2/k2

112

RA.(r), RB.(r)

Rgap

H.C. Yeh / Fuel-and-cladding model of rewetting of a nuclear fuel rod

eigenfunctions in regions A and B, respectively dimensionless gap resistance,

B 1

2

k2/(hgapr2)

T

Tf T~

Tc u

u*

temperature fluid temperature temperature of rod at far downstream of the dry side Leidenfrost temperature wet front velocity

piciuf2/ki

u

( u ~ - u 2 . ~ , ) / ( u 2 . , . ~ - u2.~,)

U 2,cyl

p2c2u~2/k2 for a uniform cyl-

U 2,tube

u~z Y.

z

An

~n 0

indrical rod pzc2u~2/k2 for a hollow tube of inner radius ~ and outer radius ~2 u ~/u~ (= plClk2/(.a2c2k,)) Bessel function of second kind of order n axial coordinate dimensionless axial coordinate, ~/~2 eigenvalue eigenvalue dimensionless temperature

( T - T,)/(T~- T,) Oc pi

dimensionless Leidenfrost temperature (To- Tf)/(T®- Tf) density

Subscripts A

quantities for region A

quantities for region B physical properties for fuel physical properties for cladding

References [1] A. Yamanouchi, J. Nucl. Sci. Tech. 5 (1968) 547. [2] B.D.G. Piggott and D.T.C. Porthouse, Berkeley Nuclear Laboratories, RD/B/N2692. [3] K.H. Sun, G.E. Dix and C.L. Tien, ASME J. Heat Transfer 96 (1974) 126. [4] K.H. Sun, G.E. Dix and C.L. Tien, ASME J. Heat Transfer 97 (1975) 360. [5] M.H. Chun and W.Y. Chon, ASME paper no. 75WA/HT-32, Winter Meeting, Houston, Texas, Nov. 1975. [6] T.S. Thompson, Nucl. Engrg. Des. 22 (1972) 212. [7] T.S. Thompson, Nucl. Engrg. Des. 31 (1974) 234. [8] T.S. Thompson, Proc. 5th Int. Heat Transfer Conference, Tokyo, Vol. IV (1974) 139. [9] R.B. Duffey and D.T.C. Porthouse, Nucl. Engrg. Des 25 (1973) 379. [t0] M.W.E. Coney, Nucl. Engrg. Des. 31 (1974) 246. [11] C.L. Tien and L.S. Yao, ASME J. Heat Transfer 97 (1975) 161. [12] J.M. Blair, Nucl. Engrg. Des. 32 (1975) 159. [13] H.C. Yeh, Nucl. Engrg. Des. 34 (1975) 317. [14] H.C. Yeh, J. Appl. Phys. 46 (1975) 4431. [15] H.C. Yeh, J. Appl. Phys. 47 (1976) 2923. [16] H.C. Yeh, Nucl. Engrg. Des. 32 (1975) 85. [17] E. McAssey and C. Bonilla, ASME paper 77-HT-93 (August 1977). [18] K.G. Pearson, B.D.G. Piggott and R.B. Duffey, Nucl. Engrg. Des. 41 (1977) 165. [19] H.C. Yeh, Nucl. Engrg. Des. 36 (1976) 139. [20] H.C. Yeh, J. Appl. Phys. 48 (1977) 4423.