Studies on flow instability of helical tube steam generator with Nyquist criterion

Nuclear Engineering and Design 266 (2014) 63–69

Contents lists available at ScienceDirect

Nuclear Engineering and Design journal homepage: www.elsevier.com/locate/nucengdes

Studies on flow instability of helical tube steam generator with Nyquist criterion Fenglei Niu a,∗ , Li Tian a , Yu Yu a , Rizhu Li b , Timothy L. Norman c a b c

State Key Laboratory of Alternate Electrical Power System with Renewable Energy Sources, North China Electric Power University, Beijing 102206, China Institute of Nuclear and New Energy Technology, Tsinghua University, Beijing 100084, China Westinghouse Electric Company, Madison, PA 15663, USA

h i g h l i g h t s • • • •

Density-wave oscillation in helical-tube steam generators was studied. The multi-variable frequency domain method was used for the modeling. The flow stability was evaluated by the Nyquist stability criterion. The calculated results are consistent with the experimental results.

a r t i c l e

i n f o

Article history: Received 9 February 2013 Received in revised form 15 October 2013 Accepted 24 October 2013

a b s t r a c t The steam generator of the 10 MW High Temperature Gas-Cooled Reactor (HTR-10) in China consists of a series of helical tubes where water/steam flows inside and helium flows outside. It operates under middle pressure, which tends to cause the flow instability. Density-wave oscillation is the most common type of two-phase flow instability in the steam generators. This paper presents the research on flow instability for the HTR-10 steam generator. The drift flux model was used for two-phase flow analysis. The transfer matrix was obtained by using linearized perturbation and Laplace transformation on the conservation equations. The flow stability was evaluated by the Nyquist stability criterion. The results obtained from frequency domain method were compared and discussed with the results from the time domain method and the experimental results. © 2013 Elsevier B.V. All rights reserved.

1. Introduction The helical tube steam generator as shown in Fig. 1 has the following advantages (Yunlong et al., 1992; Qincheng et al., 1998; Wen, 1982): (1) vortex flow takes place inside the tube because of the centrifugal force, which leads to a disturbance in the boundary layer, and a better heat transfer performance than the straight tubes; (2) the heat transfer for the compact layout is improved due to the increased surface area; (3) better thermal expansion adaptability which improves the reliability and safety of the equipment. For these reasons the helical tube steam generators were used in some commercial power plants (Olson, 1980; Elter and Franke, 1986; Dee and Macken, 1977; Henri Fenech, 1981) as well as China’s 10 MW High Temperature Gas-cooled Reactor (HTR-10) (Huaiming et al., 2004). The HTR-10 steam generator consists of 30 helical tubes. The secondary loop operates under

∗ Corresponding author. Tel.: +86 10 61773160; fax: +86 10 61773156. E-mail address: [email protected] (F. Niu). 0029-5493/$ – see front matter © 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.nucengdes.2013.10.026

middle parameters (Weihua, 1996) (e.g. the inlet water at 6.0 MPa and 104 ◦ C, and the outlet steam at 4.0 MPa and 440 ◦ C) at 100% power rate, so the flow instability that usually occurs in the negative slope portion of the hydrodynamic curve is a problem. Density-wave oscillation is one of the most representative instable flows in the steam generators. This instability mode is due to the feedback and interaction between the various pressure drop components and is caused by the lag introduced due to the finite speed of propagation of kinematic density waves (Lahey, 1977), which is often found in the boiling channels. Flow instability may cause a lot of problems, such as forced mechanical vibration of components, boiling crisis, wall temperature oscillation and system control problems, which may eventually lead to equipment failure and safety problems (Zhaoyi and Ruian, 1992; Jijun and Dounan, 1993). Previous and current research on flow instability mainly focused on vertical or horizontal tubes, whereas the research on helicalcoiled tubes has received much less attention. Several available reports mention a single variable frequency domain model (such as NUFREQ) and neglect the influences of other heat transfer channels and external loop, where the system is simplified into a single

64

F. Niu et al. / Nuclear Engineering and Design 266 (2014) 63–69

Nomenclature A CK c0 DK D(s) f G0 g h I j Ki p Pw q

v vgj w ˛  c (s) o (s)  ı g ˝ ˚f0 ˚(x)

area (m2 ) kinematic-wave velocity void distribution parameter derived density return difference matrix friction coefficient mass flux acceleration of gravity specific enthalpy unit matrix volume velocity of two phase mixture coefficient pressure drop wetted perimeter heat flux density velocity drift velocity rate of flow void fraction density closed-loop characteristic polynomial opened-loop characteristic polynomial transfer function disturbance quantity vapor produce rate frequency of phase change two phase friction multiplier two phase local resistance multiplier

channel model. Some reports analyze the two-phase section using uniform flow model or fixed uniform flow model which simplifies the calculation process but leads to greater error. This paper analyzes density-wave oscillation in HTR-10 steam generator using thermodynamic equilibrium conservation equations, drift flux model linearized perturbation theory, Laplace transformation and multi-variable frequency domain method. In the drift flux model, the relative motion between the two phases and the non-uniform distribution are taken into account. In the instability analysis process, the frequency-domain method and the time-domain method were often used. The time-domain method uses the nonlinear conservation equations, thus, are theoretically able to predict limit cycle phenomena as well as being of use for general transient analysis. The real motivation for linear frequency domain stability analysis is that exact solutions are possible, the computer costs are far less, and the threshold of instability can be predicted at least as accurately as with a nonlinear time-domain analysis (Lahey, 1977). The linear frequency-domain method can also eliminate the mathematical instability. The multi-input to multi-output transfer matrix is then obtained. The coupling factor is introduced to reflect the interrelationship between the parallel channels and the external loop. The stability of the system can be evaluated by the Nyquist criterion. These analysis results will be compared with the time-domain results and the experimental results. 2. Fundamentals of the model

(1) In the secondary loop, the heat transfer and flow parameters of every sub-segment control unit are uniform, which are represented by the mean values. (2) The single phase flow is incompressible and the two phase flow is under thermodynamic equilibrium. (3) The semi-empirical formulae of flow and heat transfer under steady state condition are applicable for this case. (4) The fluid is assumed to be ideal fluid. (5) The potential energy and the kinetic energy in the energy conservation equation are neglectable. (6) The heat conduction in the axial is neglected. 2.2. Mathematical models 2.2.1. Single phase section The fluid in the subcooled section and superheated section are single phase. The conservation equations in the single phase section of the secondary loop are: Mass equation: ∂v =0 ∂x

(1)

Energy equation:

 

∂h ∂h +v ∂t ∂x



=

qPw Ax

(2)

Momentum equation:



−dp = 

∂v ∂v +v ∂t ∂x



dx +

  v2 fv2 dx + gdz + kn ı(x − xn ) 2 2d n

(3) where v represents the velocity of the fluid, h is the specific enthalpy, q is the heat flux for the inside wall of the helical tube, Pw is the wetted perimeter, Ax is the cross sectional area of the flow channel, f is the friction coefficient, kn is the local friction coefficient, and ı(x − xn ) is a unit step function. 2.2.2. Two-phase section Mass equation: Liquid phase: ∂ ∂ [(1 − ˛)f ] + [(1 − ˛)f vf ] = −g ∂t ∂x

(4)

Gas phase: ∂ ∂ (˛g ) + (˛g vg ) = g ∂t ∂x

(5)

Energy equation: ∂ ∂ Pw q (6) [(1 − ˛)f hf + ˛g hg ] + [(1 − ˛)f hf vf + ˛g hg vg ] = Ax ∂t ∂x Momentum equation:



−dp =

+ +



∂ v2 ∂ (v) + (v2 ) dx + f˚2f dx + gdz 0 2 d ∂t ∂x f ∂ ∂x



f −   − g

 (v2 ) kn

n

2f



f g 



[vgj + (c0 − 1)]2

dx

˚(x)ı(x − xn )

(7)

2.1. Assumptions

where g is the vapor generation rate, which can be calculated by the thermodynamic equilibrium assumption g = AqPhw . ˛ is the

The following assumptions are made in the model derivation process (Rakopoulos et al., 1980; Fenglei, 1997).

void fraction, g , f ,  are the steam density, liquid density, and the average density of the two phase flow, respectively. hg , hf are

x fg

F. Niu et al. / Nuclear Engineering and Design 266 (2014) 63–69

65

Fig. 1. The steam generator and its helical-tube assembly of HTR-10.

the specific enthalpy of the steam phase and liquid phase, ˚f0 is the two-phase friction multiplier, ˚ is two-phase local resistance multiplier. vgj is the drift velocity given by

 

vgj = vg − c0 j where c0 is the void distribution parameter. Since the drift flux model is used, the Zuber type void fraction transfer equation must be included as follows (Prosperetti, 1990):

∂ + ∂x



f f − g

 − c0 ˛

+ ˝

(8)

With  = (1 − ˛)f + ˛g , it can be deduced that: ∂ ∂ + CK = −˝DK ∂t ∂x

(9)

where CK is the kinematic-wave velocity and given by

 CK = c0 j + vgj + ˛

dvgj d˛

+j

dc0 d˛



˝ is the frequency of phase change and given by ˝=

∂ (v0 ıG + G0 ıv) dx ∂x

G0 G0 2 G0 2 ıG + ˚2f ,0 ıf + f0 ı˚2f 0  d 0 0 2 d 2 f f fd

f0 ˚2f

+

∂˛ ∂˛ = + Ck ∂x ∂t



−dıp = sıG +



f g

02 + f g − 2f 0 02 ( − g )2

 ı(x − xn ) kn

n

2f

2

 dx + gdzı

[vgj + (c0 − 1)] ı

 dx

(G0 2 ı˚ + 2G0 ˚ıG)

(10)

Applying the above equations on the number i control unit, the following expression can be derived:

ıpi = ıpia + ıpif + ıpig + ıpid + ıpil = K1i ıDK,i−1 + K2i ıCK,i−1 + K3i

ı˝tp,i ˝tp,i−1

+ K4i

ıGi−1 ıv + K5i i−1 + K6i ıhi−1 Gi−1,0 vi−1,0

(11)

where Ki is the coefficient and ˝tp = c0 ˝.

g (f —g ) f g

and DK = c0  + (1 − c0 )f . With the equations above, the following expression is obtained with the Laplace transformation and linearized perturbation on the momentum equation:

3. Transfer matrix The SG flow channels are divided into a number of control units. In the control units, the relationship between pressure drop

66

F. Niu et al. / Nuclear Engineering and Design 266 (2014) 63–69

disturbance ıp and flow rate perturbation ıw can be expressed as: ıpl,j (s) = l,j (s)ıwin,j (s) ıptp,j (s) = tp,j (s)ıwin,j (s)

(12)

ıpv,j (s) = v,j (s)ıwin,j (s) The relationship between the disturbance of external loop pressure drop ıpt (s) and the flow rate perturbation ıwt (s) can be expressed as: ıpt (s) = t (s)ıwt (s)

ıwt (s) =

n 

(13)

ıwin,j (s)

(14)

j

Fig. 2. Schematic of a multi-channel feedback system.

results in the change of the flow and pressure drop in the external loop. On the other hand, the change of the pressure drop of the external loop will lead to the change of the inlet flow of each channel. These effects make a multichannel feedback system, which is shown in terms of the block diagram given in Fig. 2.

Using the linearized perturbation theory and the Laplace transformation, the system transfer matrix can be expressed as:

4. Stability criterion

ıWin (s) = −1 (s) · ıPl (s) l

According to the multichannel feedback system model derived above, the return difference matrix (Rohrs, 1993) of the system can be expressed as:

(15)

where ıPl (s) = ıPloop (s) − (s) · ıWin (s) (s) = tp (s) + v (s) + t (s)

(s) D(s) = I + (s) · −1 l and the closed-loop characteristic polynomial is given by

T

ıWin (s) = [ ıwin,1 (s)

ıwin,2 (s)· · ·ıwin,n (s) ]

c (s) = o (s) det D(s) ıPl (s) = [ ıpl,1 (s)

ıpl,2 (s)· · ·ıpl,n (s) ]

ıPloop (s) = [ ıploop,1 (s)

⎡ −1 (s) l

⎢ ⎢ ⎢ =⎢ ⎢ ⎣

T

ıploop,2 (s)· · ·ıploop,n (s) ]

T

⎤

−1 l,1 (s)

⎥ ⎥ ⎥ ⎥ ⎥ ⎦

−1 l,2 (s)

..

. −1 l,n (s)

⎡

⎤

tp,1 (s)

⎢ ⎢ tp (s) = ⎢ ⎢ ⎣

⎥ ⎥ ⎥ ⎥ ⎦

tp,2 (s) ..

.

where o (s) is the open-loop characteristic polynomial. Linear system instability occurs when the system becomes selfexcited, so self-excited should be avoid (Lahey, 1977). According to the multi-variable frequency domain control theory, the Nyquist stability criterion is that the closed-loop system is stable if and only if the net number of clockwise encirclements of the unity point by the Nyquist diagram of det D(s) plus the number of zeros of o (s) in the right half-plane of the complex plane is zero. That is encD(s) = −n0 Here, according to Eqs. (12)–(14), n0 = 0. From the Eqs. (12)–(15), we get det D(s) =

tp,n (s)

⎡ ⎢ ⎢ v (s) = ⎢ ⎢ ⎣

⎤

v,1 (s)

⎥ ⎥ ⎥ ⎥ ⎦

v,2 (s) ..

. v,n (s)

⎡

t (s)

t (s)

···

t (s)

⎤

+

t (s)

···

t (s)

n×n

When perturbation of the inlet flow occurs, flow disturbance causes enthalpy disturbance, leading to the changes in the boiling boundary, superheating boundary, and density distribution inside the channel. The changes of flow, density and boundary positions cause the changes of pressure drop in the liquid area, two-phase area, and superheated steam area. The flow change of each channel

1+

⎧ ⎪ ⎪ ⎪ ⎪ n ⎪ ⎨  j=1

tp,j (s) + v,j (s)

j=1



l,j (s)

⎫ ⎪ ⎪ ⎪ ⎪ ⎪   n ⎬  tp,k (s) + v,k (s) t (s) 1+ l,j (s) l,k (s) ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ k=1 ⎪ ⎪ ⎪ ⎪ ⎩ ⎭ k= / j

With the assumption that: f1 (s) =

⎢ ⎥ ⎢ t (s) t (s) · · · t (s) ⎥ ⎢ ⎥ t (s) = ⎢ ⎥ ⎣ t (s) t (s) · · · t (s) ⎦ t (s)

n  

n  

1+

j=1

=

1 n

tp,j (s) + v,j (s) l,j (s)

⎧ ⎪ ⎪ ⎪ ⎪ ⎪ n ⎨ j=1

f2 (s)

⎫ ⎪ ⎪ ⎪ ⎪ ⎪   n ⎬  tp,k (s) + v,k (s) pt /wt 1+ kc l,k (s) l,j (s) ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ k=1 ⎪ ⎪ ⎪ ⎪ ⎩ ⎭ k= / j

=−



nt (s)wt pt

(16)

F. Niu et al. / Nuclear Engineering and Design 266 (2014) 63–69

Fig. 3. Nyquist curves of the system with the same coupling factor but different operating pressure.

Eq. (16) can be expressed as: det D(s) = f1 (s) − kc f2 (s)

(17)

where kc is named by the coupling factor. 5. Results and discussion Fig. 3 shows that system pressure affects density-wave oscillation. The Nyquist curves gradually get away from the origin as the system pressure increase, which means the system stability increases with the operating pressure, because the higher pressure reduces to the smaller density difference between the vapor and the water. When the system reaches critical pressure, the density of water is the same as the vapor, so flow instability will not happen. In addition, as the system pressure increases, the saturation temperature will rise, and therefore the bubble content in the channels and the pressure drop in two-phase area will decrease, which eventually improves the system stability. The Nyquist curves gradually get away from the origin as the coupling factor decreases, as shown on Fig. 4, which means the smaller coupling factor can improve the stability of the system. For the multichannel system, when the inlet oscillation occurs in channel A, the inlet flow of the channel B will also oscillate due to the change of the overall external flow. If the coupling factor decreases to zero, the multichannel system becomes a single channel system and, therefore, the influence of the external loop will disappear. Flow drift is prone to occur in the negative slope area of pressure drop-flow curve in heating channels, as shown in Fig. 5, leading to the two-phase flow instability. The flow rate has two distinct

Fig. 4. Nyquist curves of the system with the same operating pressure but different coupling factors.

67

Fig. 5. The hydrodynamic characteristics in the heated water channels.

effects on pressure drop in the boiling section. The pressure drop caused by friction decreases with reduction of the flow rate. As the flow decreases, however, the boiling takes place and the twophase pressure drop significantly increases due to the contribution of the acceleration pressure drop. Thus, whether the pressure drop increases or not depends on which effect plays the major role. Similarly, increase of the system pressure can reduce the density difference between the water and the vapor, which results in the decrease of the acceleration pressure drop. However, higher system pressure can result in shorter two-phase and superheated steam length, which will reduce the effect of the negative slope area on Fig. 5, and therefore increase the stability of the system. Fig. 6 shows that flow rate affects the density-wave oscillation. The Nyquist curves gradually get away from the origin as the flow rate increase, which means the system stability increases with flow rate. Fig. 7 shows the Nyquist curves gradually get away from the origin as the power rate increases, which means the system stability increases with the power rate. That is because the inlet mass flow increases accordingly with the power rate. As discussed in Fig. 6, the system stability increases with flow rate. Fig. 8 shows the Nyquist curves gradually get away from the origin with the increase of outlet mass quality, which means the system stability increases with the outlet quality. Under the same operation pressure, flow rate, and subcooling, the outlet quality increases with the heating power. When the outlet quality increases, the two-phase and the superheated section will occupy more proportion in the whole heating channel, as the result, the instability will occur more likely.

Fig. 6. Nyquist curves of the system under the same power but different flow rate.

68

F. Niu et al. / Nuclear Engineering and Design 266 (2014) 63–69

Fig. 7. Nyquist curves of the system under the different power rate. Fig. 10. Single pipe flow pulsation amplitude vs. operating power.

Fig. 11. Single pipe flow pulsation amplitude vs. feed water pressure.

Fig. 8. Nyquist curves of the system for the different outlet quality.

Fig. 9 shows that inlet resistance coefficient affects the densitywave oscillation. The Nyquist curves gradually get away from the origin as the inlet resistance coefficient increases, which means the system stability increases with the inlet resistance coefficient. The single-phase pressure drop increases with the inlet resistance, so higher inlet loss coefficient can give longer positive slope section on Fig. 9. Thus, increasing the inlet resistance can serve as a means of suppressing the density-wave instability, but higher inlet resistance will lead to higher pump load. When the heat flux and inlet subcooling are fixed, the increase of flow rate will shorten the length of the two-phase section, which will reduce the possibility that the flow is in the negative slope area

Fig. 9. Nyquist curves of the system under the same power but different inlet resistance.

on Fig. 5. Thus, the increase of flow rate will enhance the stability of the system. The flow instability tests were performed with a full-scale steam generator in the Institute of Nuclear Energy Technology at Tsinghua University (Rizhu and Huaiming, 2002). At 20% power or below, the flow will be unstable and the pressure drop will oscillate between parallel pipes. At 30% power, the experimental results show that the flow is stable for outlet pressure from 2.5 to 4.0 MPa and is stable for inlet temperature from 75 ◦ C to 180 ◦ C, but the system will be unstable if the inlet orifice resistance is less than 40 kPa. The oscillation amplitude is more than 5%. Therefore, three throttling orifices are necessary at the inlet to maintain system stability. These experimental results are consistent with the calculation results above. For further verification of the frequency domain method used on this paper, a 2-D model based on the time domain method was developed. The experimental and the calculation results are shown in Figs. 10–13, which confirm that the flow stability will increase

Fig. 12. Single pipe flow pulsation amplitude vs. subcooling.

F. Niu et al. / Nuclear Engineering and Design 266 (2014) 63–69

69

results got by the frequency domain models are consistent with the experimental results and calculation results from the time domain method. Acknowledgments This work was supported by the National Natural Science Foundation of China (91326108, 51206042). References

Fig. 13. Single pipe flow pulsation amplitude vs. inlet throttling.

with the flow rate, operating pressure, inlet subcooling, and the inlet throttling degree. The inlet throttling degree defined as



=

Pi + Ppre Pboi + Psup + Pst

(18)

where Pi , Ppre , Pboi , Psup , and Pst are the pressure drops in the inlet, the subcooling section, the two-phase section, the superheated section, and the thermal isolation section, respectively. The time domain model has a good agreement with the frequency domain model in predicting that the flow rate, operating pressure, inlet subcooling, and inlet resistance are the important factors influencing the flow instability occurred in the helical tube steam generator, but the computation with the time domain method is more time consuming, and, moreover, the numerical method for solving partial differential equations may itself introduce instability problems to the system. 6. Conclusions The multi-variable frequency domain method was used to study the flow instability problems in the helical tube steam generators, and the transfer matrix of the multichannel system was obtained. The flow stability in the HTR-10 steam generator is strongly dependent on the operating pressure, the coupling factor, the inlet resistance, the flow rate, the power rate, and the outlet quality. The

Dee, J.B., Macken, T., 1977. Development of components for the gas-cooled fast breeder reactor program. Nuclear Engineering and Design, 13–26. Elter, C., Franke, R., 1986. Gas pressure test of THTR steam generators. International Journal of Pressure Vessels and Piping, 37–48. Fenglei, N.I.U., 1997. Flow Instability Analysis of 10 MW High Temperature Gascooled Test Reactor Steam Generator. INET, Tsinghua University, Beijing (Master degree dissertation). Henri Fenech, 1981. Heat Transfer and Fluid Flow in Nuclear System. Pergamon Press, pp. 510–520. Huaiming, J.U., et al., 2004. Experiment and verification test of the once-through steam generator of the 10 MW high-temperature gas-cooled reactor – flow stability of the once-through steam generator. Nuclear Science and Technology 41, 524–528. Jijun, X.U., Dounan, J.I.A., 1993. Boiling Heat Transfer and Steam Liquid Two-phase Flow. Atomic Energy Press, Beijing. Lahey, R.T., Moody, F.J., 1977. The Thermal Hydraulics of a Boiling Water Nuclear Reactor. American Nuclear Society, Illinois, pp. 331. Olson, H.G., et al., 1980. The Fort St. Vrain high temperature gas-cooled reactor: V. Steam system performance tests. Nuclear Engineering and Design, 303–314. Prosperetti, A., Satrape, J.V., 1990. Stability of Two-phase Flow Models. Two Phase Flows and Waves. Springer-Verlag New York Inc., New York, pp. 98–117. Qincheng, B.I., et al., 1998. Study on wall temperature profile and heat transfer deterioration in helical coils of steam generator of HTGR. Nuclear Power Engineering 10, 408–412. Rakopoulos, C.D., et al., 1980. An experimental-stochastic and theoretical analysis of the density wave instability in a helical monotube vapour generator. Nuclear Engineering and Design, 369–384. Rizhu, L.I., Huaiming, J.U., 2002. Structural design and two-phase flow stability test for the steam generator. Nuclear Engineering and Design, 179–187. Rohrs, C.E., et al., 1993. Linear Control Systems. McGraw-Hill Inc., New York, pp. 307–312. Weihua, L.I., 1996. Thermal Hydraulic Calculation Report of 10 MW High Temperature Gas-cooled Test Reactor Steam Generator. INET, Tsinghua University. Wen, G.A.O., 1982. High Temperature Gas Cooled Reactor. Atomic Energy Press, Beijing. Yunlong, Z.H.O.U., et al., 1992. A design procedure of steam generator with multistart helical tubes. Nuclear Power Engineering 6, 1–8. Zhaoyi, Z.H.A.O., Ruian, Z.H.U., 1992. Thermal Fluid Mechanics of Reactor. Tsinghua University Press, Beijing.