An exact solution for the moving boiling boundary problem

An exact solution for the moving boiling boundary problem

Nuclear Engineering and Design 203 (2001) 243– 248 www.elsevier.com/locate/nucengdes Technical note An exact solution for the moving boiling boundar...

101KB Sizes 0 Downloads 321 Views

Nuclear Engineering and Design 203 (2001) 243– 248 www.elsevier.com/locate/nucengdes

Technical note

An exact solution for the moving boiling boundary problem D. Hasan, Y. Nekhamkin, V. Rosenband, E. Elias *, A. Gany, E. Wacholder Faculty of Mechanical Engineering, Technion-Israel Institute of Technology, Haifa 32 000, Israel Received 30 May 2000; accepted 5 July 2000

Abstract An exact transient solution of the fluid velocity and temperature fields in a one-dimensional incompressible flow within a non-uniformly heated channel is presented. The first order partial differential equations for mass and energy conservation governing this problem are solved using Laplace transform technique. An analytical expression for the boiling boundary (BB) location as a function of time is derived from the temperature field solution when the saturation temperature is inserted. Results obtained reveal the interesting behavior of the temperature field and BB in space and time due to a step change in the fluid inlet mass flow rate. © 2001 Elsevier Science B.V. All rights reserved.

1. Introduction The problem of convective boiling in a heated channel is of practical importance for the design and operation of heat exchanger and for the analysis of reflood phenomena in a reactor core during accident conditions. This paper considers sub-cooled incompressible fluid undergoing convective boiling while flowing upward in a non-uniformly heated channel. A schematic description of the system and its boiling boundary (BB) is presented in Fig. 1. The fluid enters as single-phased liquid at the bottom of the flow channel. At a

* Corresponding author. Tel.: + 972-4-8293263; fax: + 9724-8324533. E-mail address: [email protected] (E. Elias).

certain level along the channel boiling starts and a two-phase (vapor –liquid) region is established. The position of the interface between the singlephase region and the two-phase region is defined as the position at which the fluid phase changes due to the commencement of boiling, and is therefore referred to as a ‘boiling boundary’. It is thermodynamically defined as the position where the enthalpy of the single-phase fluid is equal to the saturation enthalpy. Transient thermal-hydraulic evaluation of boiling systems requires the accurate prediction of the distributions of flow and enthalpy at each point in space and time. In practice, this is normally accomplished by a simultaneous numerical solution of the coupled two-phase conservation equations in conjunction with adequate two-phase thermophysical properties’ tables, flow pattern assump-

0029-5493/01/$ - see front matter © 2001 Elsevier Science B.V. All rights reserved. PII: S 0 0 2 9 - 5 4 9 3 ( 0 0 ) 0 0 3 4 2 - 3

244

D. Hasan et al. / Nuclear Engineering and Design 203 (2001) 243–248

tions, and constitutive relations. An important requisite of the analysis is the determination of the transient position of the BB. Across this boundary there is discontinuity in several different physical conditions and thermo-physical properties such as compressibility, velocity of sound, vapor quality, fluid velocity, constitutive relations and heat transfer regime transitions (Wacholder et al., 1985). To establish the transient position of the BB, one must study the dynamics of the single-phase portion of the heated channel as is demonstrated in the present work. The notable work of Gonzalez-Santalo and Lahey (1973) provides a solution by the method of characteristics for both the single-phase and the two-phase regions. However, their work is limited to a simplified case of one-dimensional flow, with constant (in space and time) heat flux. Solutions for a wide range spectrum of moving boundary problems involved in melting and freezing processes are comprehensively studied in Alexiades and Solomon (1993). In the present work, the single-phase region in a one-dimensional flow channel with arbitrary spatial distribution of heat flux has been solved using Laplace transform. The analytical solution of the fluid temperature field provides means for algebraic derivation of analytical expressions for the instantaneous position of the moving boundary and for the time-span to reach it.

Fig. 1. Schematic of a variably heated channel cooled by fluid flow.

2. Problem formulation The thermo-physical properties and flow model assumptions made in the subsequent analysis are: “ Constant density (incompressible flow) and constant thermodynamic properties in the liquid phase. “ Constant system pressure. “ One-dimensional flow. “ Non-uniform stationary heat rate per unit length of channel wall: q%= q%0 sin (yz/L). “ Negligible viscous dissipation in the fluid flow. “ Negligible heat capacity of the channel wall. “ The flow conditions at the channel bottom inlet (i.e. mass flux and temperature) are specified as a function of time. With these assumptions, the mass and energy conservation equation for the liquid flow can be written as: (V =0, (z

(Continuity equation)

 

(q (q yz +V =i sin , (t (z L

(1)

(Energy equation) (2)

where, q(z, t) = (T− T0)/(Ts − T0) and i=q%0/ zACf(Ts − T0). The parameters z and Cf are the liquid density and specific heat, respectively. The flow within the heated channel is assumed to be initially at steady state. This is achieved by maintaining an initial liquid inlet mass flow rate, m; 0, which balances the amount of liquid evaporated from the single-phase/two-phase interface and thus maintains a constant BB at a level v0L (v z/L). The inception of the system transient is due to a step change in the incoming mass flowrate from m; 0 to um; 0, where u is a prescribed factor which may be smaller or larger than 1. Since the fluid in the channel is assumed incompressible, any velocity or mass flow rate disturbance in the channel inlet propagates instantaneously through the liquid region (i.e. infinite wave speed of propagation). The corresponding initial conditions of the problem at t= -0, 05z/L5 v0 are therefore: V(z, -0)=V0 = m; 0/Az

(3)

q(z, 0)=



 

D. Hasan et al. / Nuclear Engineering and Design 203 (2001) 243–248

yz 1 1 −cos L €

(4)

245

qs(0, s)= 0, at z= 0

(10)

The general solution of Eq. (9), subject to the boundary condition (10) is:

where, 1 iL = € yV0

   

qs(s, z)= a+b sin

The boundary conditions (at z =0, t ] 0) are: V(0, t)= Vi =V0u

(5)

q(0, t)= 0

(6)

where, Vi is the fluid inlet velocity at the heated channel bottom at time t ]0. m; 0 and V0 are the initial fluid inlet mass flow rate and velocity, respectively, at time t B 0. The BB is defined, in temperature terms, as the fluid height coordinate v1L, at which the fluid reaches saturation temperature,Ts. The BB relative height, v0, in the steady state flow before the transient inception, is obtained from a singlephase region heat-balance: 1 v0(tB 0)= cos − 1(1 − €) y

 

yz yz + c cos L L

+ d exp −

sz V

(11)

where, the coefficients a, b, c, d are: Vy −i €L , b= − y 2 2 s + V L iL 1 1 iL s c= − − − Vy s € Vy 2 y s + V L

    

1 a= , s€



d=



iL 1 − Vy €

<

(7)

1 − s



s

s2+ V

y L

2

,

= 2

Substituting the inverse transform expressions of all the functions of s in Eq. (11) yields the following solution for the non-dimensional temperature field:

3. Method of solution

<

3.1. Temperature distributions

q(z, t)

The velocity, V(z, t), is obtained from the continuity Eq. (1) and boundary condition (5) as:

yz Á cos Ã1 L u− 1 z = Í 1− − cos y − € u u L à Ä

V(z, t) =Vi =m; 0u /Az =const.

(8)

Substituting this value into Eq. (2) yields an uncoupled first order partial differential equation for the solution of the temperature field. The Laplace transform operator is then applied to each term in Eqs. (2), (4) and (6), in order to obtain an ordinary differential equation for the solution of qs(s, z) — the transformed temperature function of q(t, z), i.e.

 



 

=

(9)

 



 

z z 1 yz t5 = trs 1− cos V L u€ L

z z t] = trs V L (12)

where, trs is the asymptotic residence time defined as trs =

i yz 1 yz dqs s + qs = sin + 1 −cos dz V L L Vs V€

Here s is the Laplace transform parameter. The corresponding boundary condition is:

t trs

 

zAL L = . um; 0 V

The first branch of the temperature solution (Eq. (12), for t5(z/V)= (z/L)trs), describes the temperature distribution, at each point in space and time, of all those liquid particles which al-

D. Hasan et al. / Nuclear Engineering and Design 203 (2001) 243–248

246

ready were in the liquid region at transient inception, i.e. at t= 0. The second branch of the solution (i.e. for t ] (z/V) =(z/L)trs), describes the new steady state temperatures at different z coordinates at times larger than z/V. It refers to all those liquid particles which entered the channel after inception of the transient, i.e. at t]0. It also describes the new temperature spatial distribution in the liquid region, established after transient termination. The range of definition of z in both temperature solution branches is: 05 z5 v1(t)L, where v1(t) is the instantaneous relative height of the BB.

3.2. BB tracking Variations in time of the BB relative height, v1(t), is tracked down by the solution of Eq. (12), for the case T =Ts. Here, Ts is the saturation temperature of the fluid at the channel pressure. Substitution of T = Ts into the first branch of the temperature solution of Eq. (12), yields the relative height of the BB, v1(t), as a function of time t, in terms of the mass flow rate step-factor, u, the asymptotic residence time, trs, and€: 1 F2(t) − F 22(t) − F1(t)F3(t) v1(t)= cos − 1 y F1(t) where,

  n 

F1(t)= 1+(u− 1)2 +2(u −1) cos



(13)

yt trs

F2(t)= u(1−€) 1 +(u −1) cos

yt trs

F3(t)= u 2(1− €)2 −(u− 1)2 sin2

yt trs

The second branch of the temperature solution yields the relative height of the BB for the new steady-state, for every mass flow rate step-factor u and €: 1 v1 = cos − 1(1 − u€) y

(14)

Thus, it is summarized that the value of v1 is changing from its initial value given in Eq. (7), through the transient according to Eq. (13), until reaching its final value given in Eq. (14), i.e.:

1 1 cos − 1(1−€)“ v1(t)“ cos − 1(1−u€) (15) y y The time-span for the transient BB relative height, Eq. (13), to reach a specified value, v1(t), is given by





1 u(1− €)− cos(yv1) t= trs v1 − cos − 1 y u− 1

n

(16) This time is bounded by the maximum transient time, (trs/y)cos − 1(1−u€), which is the time required for a fluid particle entering the channel bottom inlet at velocity Vi = V0u to reach the new steady-state BB height. After this time, the BB remains at a constant position, as expressed by Eq. (14).

4. Results and discussion Quantitative results of the present analysis are illustrated by the example of a partially uncovered AP600 type nuclear reactor core channel, 720 s after reactor scram (shutdown). The assumed initial condition refers to single-phase liquid entering at the bottom of the channel. The inlet flow maintains a constant single-phase level of v0L, as shown in Fig. 1. The initial fluid inlet mass-flow rate, m; 0 = 1.884×10 − 4 kg/s, is determined to balance the amount evaporated, when the top of the two-phase region is at the relative height of v2 = 0.55 (Hasan, et al., 1999). The respective fluid inlet velocity is V0 = 1.84 cm/sec. The pressure within the channel is atmospheric corresponding to the reflooding phase of a hypothetical accident. The liquid inlet and saturation temperatures are T0 = 300 K and Ts = 373 K, respectively. The coolant heat of vaporization at the channel pressure is r=2.44×106 J/kg. The maximum input heat rate per unit length, transferred to the coolant, due to the reactor decay heat, is q%0 =358 W/m, with a total active core height of L=3.66 m, and flow channel cross sectional area of A= 0.1026×10 − 4 m2. Coolant (water) density and heat capacity are z= 1000 kg/m3 and Cf = 4.187 J/kg K. The initial BB relative height under these conditions is at v1(t=0)= 0.169. The transient starts when at t= 0 the inlet flow is instanta-

D. Hasan et al. / Nuclear Engineering and Design 203 (2001) 243–248

Fig. 2. Propagation of boiling boundary (BB) in space and time for different values of u.

247

larger. Once the single-phase region is filled with new fluid particles, which entered the channel after transient inception, a new steady state is established. The greater the u, the shorter the ensuing transient time and the higher the new steady BB position. The dimensionless form of the single-phase fluid temperature q =(T−T0)/(Ts − T0) is analogous to the concept of the thermodynamic steam-quality that is valid for sub-cooled liquids as well as for the two-phase region. This dimensionless temperature’s spatial distributions are depicted in Fig. 3, for the case of u= 1.4 and various times (i.e. t=0, 5, 10 and 28.65 s). The bounding time for the new steady state, is in this case t=28.65 s, as obtained from Eq. (16). After this time the temperature spatial distribution remains unchanged, as could also be obtained from the first branch of the temperature solution, Eq. (11) for t=28.65. A different view of the dimensionless temperature time-histories at various relative heights (i.e. z/L= 0.1, 0.12, 0.14, 0.16, 0.169, 0.177, 0.185, 0.193 and 0.2012), is depicted in Fig. 4 for the case of u= 1.4. The curves describing the temperature –time variations at heights larger than 0.169 refer to positions which were within the two-phase region at the transient inception. The temperature at these points preserve their saturation temperature values (i.e. dimensionless temperatures equal

Fig. 3. Non-dimensional temperature spatial distributions at different time instances, u= 1.4.

neously multiplied by a factor of u. Consequent temperature distribution and BB position histories are now derived from Eqs. (11) and (12), respectively. Fig. 2 depicts BB position histories for different step size in the inlet mass flow rate, u. The BB level is stabilized at new heights after residence times (trs/y) cos − 1(1 − u€), corresponding to the upper transient time bounds of Eq. (16). The locus of these points, depicted as a broken curve in Fig. 2, shows that the new steady state is achieved earlier with increasing u. This is explained by the fact that the single-phase region is filled with new liquid particles faster when u is

Fig. 4. Non-dimensional temperature variation in time at different heights, u = 1.4.

248

D. Hasan et al. / Nuclear Engineering and Design 203 (2001) 243–248

to one) until the moving BB reaches their position. After that, the temperature decreases in time until it eventually converges to its ultimate steadystate value. The locus of the time points where new steady states were established, depicted here as an ascending curve, describes the main particle path characteristic z =Vt.

References Alexiades, V., Solomon, A.D., 1993. Mathematical Modeling

.

of Melting and Freezing Processes. Hemisphere, Washington, USA. Gonzalez-Santalo, J.M., Lahey, R.T. Jr., 1973. An exact solution for flow transients in two-phase systems by the method of characteristics. ASME J. Heat Transfer (95)4, 470– 476. Hasan, D., Nekhamkin, Y., Rosenband, V., Elias, E., Wacholder, E., Gany, A., 1999. Nuclear Fuel Cores’ Thermal Hazards’ Criteria, Proc. INS-20, The Dead Sea, Israel, December, 1999. Wacholder, E., Kaizerman, S., Toomarian, N., Cacuci, D.G., 1985. An exact sensitivity analysis of a simplified transient two-phase flow problem. Nucl. Sci. Eng. 89 (1), 1 – 35.