Steady-state and transient analysis of two-phase natural circulation systems

Dr3 ρ2f β T gQh H Ar μ3 Cp

Dr3 ρ2f β h gQh H Ar μ3

Ns P i51




fLs Ks 1 DA2 A2

 i

1 φ2LO

P

NLh NS




K tp fLB 1 DA2 A2

 i

1 φ2LO

Nh l

P

NLh




Ktp fLhl 1 DA2 A2

h

Gr m

 i

1 φ2LO

!   Ntp f Ltp c K tp P 1 DA2 A2 Nh l

i

i

Steady-state and transient analysis Chapter | 5

Coolant

Coolant

Condenser

Condenser

SD

2036

Steam drum 885

2210

1515

Preheater

2445 575

203

850

Heater

1180

Heater 1180

655

2020

3400

(A)

(B)

Ress

FIGURE 5.3 Schematic of the different loops considered: (A) 7, 9.1, 15.74, and 19.86 mm diameter loops; (B) 49.3 mm diameter loop.

10

6

10

5

10

4

1/2.75

Ress = 1.96 (Grm /NG )

10

TINFLO-S prediction RELAP5/MOD3.2 homogeneous model

3

10

8

10

9

10

10

10

11

10

12

10

13

10

14

Grm /NG FIGURE 5.4 Comparison of flow rate predictions.

of Eq. (5.43) for open loops could result in an error of B15%. Therefore one should use Eq. (5.55) for open loops for better accuracy. TINFLO-S code predictions are in reasonable agreement with the dimensionless flow equation.

5.2.4

Comparison of the homogeneous and the two-fluid models

A comparison of the predictions of the two-fluid and the homogeneous models using the RELAP5/MoD3.2 code was made by Gartia et al. (2006).

204

Single-phase, Two-phase and Supercritical Natural Circulation Systems

10

6

Ress

1/2" FPTIL-RELAP5 (two fluid) 3/4" FPTIL-RELAP5 (two fluid) 1" FPTIL-RELAP5 (two fluid) 1/2" FPTIL-RELAP5 (homogeneous) 3/4" FPTIL-RELAP5 (homogeneous) 1" FPTIL-RELAP5 (homogeneous)

10

5

4

10 12 10

10

13

10

14

Grm /NG FIGURE 5.5 Two-fluid and homogeneous model predictions.

The predictions were made for the two-phase NCL schematically represented by Fig. 5.3A with an inside diameter of 9.1, 15.74, and 19.86 and the results are shown in Fig. 5.5. The predictions show that the error introduced by the use of the homogeneous model is only marginal when compared to the twofluid model for the steady-state two-phase NC. This is attributed to the fact that there is practically no thermal nonequilibrium at the low flow rates typical of two-phase NC and the velocity difference between the vapor and liquid is not significant enough to cause a large difference in prediction. However, for the 9.1 mm loop, the difference is perceptible due to the larger difference in velocities existing in smaller diameter loops. It may be noted further that the homogeneous model is found to reasonably predict other two-phase NC phenomena like stability (Furutera, 1986; Lee and Lee, 1991; Wang et al., 1994a,b; Nayak et al., 1998; Van Bragt and Van Der Hagen, 1998). Thus there are indications that the use of the homogeneous model for the two-phase NC phenomena may not result in serious error. However, this may not hold true for extreme transients like LOCA with cold water injection where the nonequilibrium effects are dominant. Fig. 5.6 shows a comparison of the predicted dimensionless flow rate with data from uniform diameter and nonuniform diameter loops. While plot  ting the data β h 5 ρin 2 ρe =ðρm Δiss Þ was used. In addition, the McAdams et al. (1942) model was used to evaluate the two-phase friction factor. The steady-state data used included four uniform diameter loops from Dubey et al. (2004) and Naveen et al. (2000) and two nonuniform diameter loops

Ress

Steady-state and transient analysis Chapter | 5

10

6

10

5

10

4

FIGURE 5.6 Comparison with test data.

0%

0%

+4

10

205

–4

3

Ress = 1.96(Grm /NG )

0.364

UDL data NDL data

10

2

10

7

10

9

10

11

10

13

10

15

Grm / NG (Mendler et al., 1961) and the integral test loop (Rao et al., 2002). The data fall in the following parameter range: Loop diameter: 9.1
99.2 mm; quality: 0.4%
24%, Circulation length: 8.7
96.6 m; inlet subcooling: 0.1
29 K, SD level: 90% 6 5%; pressure: 0.1
13.8 MPa, and power: 0.3
300 kW. The measured dimensionless flow rate as a function of Grm/NG is plotted in Fig. 5.6 for various loops. Two lines representing 1 40% and 240% deviations from the dimensionless flow equation are also shown in the figure. As can be seen, the agreement between the data and the dimensionless flow equation is within 6 40% which is reasonable.

5.3 Parametric effects on two-phase natural circulation systems One of the drawbacks of a dimensionless equation is that it often masks the parametric trends. In fact the steady-state performance of two-phase loops is significantly influenced by the operating conditions (i.e., pressure, power and inlet subcooling), geometric conditions (loop diameter and elevation difference) and the model used to calculate the two-phase friction multiplier. The influences of orificing and acceleration term are also important for twophase NCSs. The above-mentioned influences can be studied with the help of the simple steady-state flow rate Eq. (5.25) given above with the resistance evaluated by Eq. (5.26) for open loops and that given in Table 5.1 for closed loops. The thermal expansion coefficient was evaluated using Eq. (5.61) 2 with φ2LO and φ LO evaluated by Eq. (5.59). The calculation of steady-state flow rate is iterative. The calculations start with an assumed flow rate. With the assumed flow rate, calculate the single-phase heated length (Ls), boiling

206

Single-phase, Two-phase and Supercritical Natural Circulation Systems

length (LB), exit quality, and enthalpy distribution using the steady-state solution of the energy equation. Subsequently calculate the enthalpy integral appearing in the momentum equation. Knowing the exit quality, the density 2 and φ2LO and φ LO are estimated. Subsequently, the total hydraulic resistance is estimated and then the flow rate using Eq. (5.25). If the estimated flow rate is different from the assumed flow rate, the calculations are repeated until a specified convergence criterion is satisfied.

5.3.1

Two-phase natural circulation flow regimes

An interesting aspect of two-phase NC behavior is illustrated by Fig. 5.7. It shows that the two-phase NC flow rate can increase (Fig. 5.7A), decrease (Fig. 5.7B), or remain practically invariant with power (Fig. 5.7C). Based on the nature of the variation of the steady-state flow with power, three different NC flow regimes can be observed for two-phase loops. These flow regimes are designated as gravity-dominant (or buoyancy-dominant), frictiondominant, and the compensating regimes. It may be noted that the gravitydominant and friction-dominant flow regimes in two-phase NCLs were reported first by Fukuda and Kobori (1979). In a NCL, the gravitational pressure drop (being the driving pressure differential) is always the largest component of pressure drop and all other pressure drops (friction, acceleration, and local) must balance the gravity (or buoyancy) pressure differential at steady state. However, the NC flow regimes are differentiated based on the change of the pressure drop components with quality (or power).

5.3.1.1 Gravity-dominant regime The gravity-dominant regime is usually observed at low qualities. In this regime, for a small change in quality there is a large change in the void fraction (see Fig. 5.8) and hence the density and buoyancy force. The increased buoyancy driving force is to be balanced by a corresponding increase in the retarding frictional force that is possible only at a higher flow rate. As a result, the gravity-dominant regime is characterized by an increase in the flow rate with power (compare the gravity-dominant regime in Fig. 5.9 with Fig. 5.7A). 5.3.1.2 Friction-dominant regime The friction-dominant regime is observed at low to moderate pressures when quality is high. At higher qualities and moderate pressures, the increase in void fraction with quality is marginal (see the curve for pressures in the range 1
15 bar in Fig. 5.8) leading to an almost constant or marginal increase in buoyancy force. However, the continued conversion of highdensity water to low-density steam at high qualities (due to the increase in power) requires that the mixture volumetric flow rate and hence the mixture

Steady-state and transient analysis Chapter | 5

FIGURE 5.7 Effect of power on the two-phase NC flow rate. (A) Flow increases with power; (B) flow decreases with power; (C) flow nearly constant with power.

(A)

Mass flow rate (kg/s)

2.1

49.3 mm i.d. loop Pressure: 9.6 bar Subcooling: 1 K

1.8

1.5 Equation (5.25) Experimental data RELAP5/MOD3.2 TINFLO-S

1.2

10

20

30

40

Power (kW) (B)

Mass flow rate (kg/s)

0.03

0.02

Experimental data TINFLO-S Equation (5.25)

0.01 9.1 mm diameter tube 12.6–14.6 bar (abs) 0.00

1500

2000

2500

3000

3500

4000

Power (W)

(C)

Mass flow rate (kg/s)

0.04

0.03

9.1 mm i.d. loop Δ Tsub : 8–16K

0.02

Pressure: 57 + 1 bar 0.01

0.00 2000

Experiment RELAP5/Mod 3.2 TINFLO-S Equation (5.25) 3000

207

4000

Power (W)

5000

208

Single-phase, Two-phase and Supercritical Natural Circulation Systems

1 b ar

ba

r

1.0

70

b ar

ba 0

0.6

17

Void fraction

r

15

0.8

22

0

r ba

1 22

.2

ba

r

0.4

0.2

0.0 0.0

0.2

0.4

0.6

0.8

1.0

Quality FIGURE 5.8 Effect of quality on void fraction.

FIGURE 5.9 Effect of pressure on flow regime.

Flow normalized to 100% FP

6 GDR

Friction dominant region

5 15 bar

Predicted by Eq. (5.26) Loop diameter = 9.1 mm

4 3 2

Compensating regime 70 bar

1

170 bar

0

0

20

40

GDR

60

80

100

Power (%) velocity must increase, resulting in an increase in the frictional force and hence a decrease in the flow rate. Thus the friction-dominant regime is characterized by a decrease in flow rate with an increase in power (compare the curve for the FDR for 15 bar given in Figs. 5.7B and 5.9).

5.3.1.3 Compensating regime Between the gravity-dominant and friction-dominant regimes, there exists a compensating regime, where the flow rate remains practically unaffected by

2.0 GDR

1.6

Flow normalized to 100% power

Flow normalized to 100% power

Steady-state and transient analysis Chapter | 5

Friction dominant region 15 bar

Compensating regime

1.2 a 50 b

r

0.8 17

GDR Predicted by Eq. (5.26) Loop diameter: 19.86 mm GDR: Gravity dominant regime

0.4

0.0

r 0 ba

0

20

40

60

80

100

Power (%) (A)

209

2.0 1.6

GDR ar 2b

1.2

Compensating regime

r 8 ba

0.8

bar 170

0.4 0.0

Friction dominant regime

0

20

GDR

Predicted by Eq. (5.26) Loop diameter = 49.3 mm GDR: Gravity dominant regime

40

60

80

100

Power (%) (B)

FIGURE 5.10 Effect of pressure on flow regime for different loop diameters.

an increase in power. In this regime, the increase in buoyancy force is compensated by a corresponding increase in the frictional force, leaving the flow rate practically unaffected (compare Fig. 5.7C and the curve for 7 MPa in Fig. 5.9) in spite of an increase in quality.

5.3.1.4 Effect of pressure and loop diameter on the flow regimes The NC flow regimes depend strongly on the system pressure (see Figs. 5.9 and 5.10). In fact, at high pressures only the gravity-dominant regime (as in single-phase NC) may be observed if the power is low. The frictiondominant regime shifts to low pressures with an increase in loop diameter (see Fig. 5.10B). Knowledge of the flow regimes is important to understand the stability behavior of two-phase loops. 5.3.2

Effect of operating conditions

The operating conditions for two-phase NCLs depend on the type of loop considered. For example, in the case of open loops of the type shown in Fig. 5.2B, the operating conditions include power, pressure, inlet subcooling, and level. In the case of open loops of the type shown in Fig. 5.2A and closed loops (Fig. 5.1), the inlet subcooling is a dependent parameter. Therefore the operating parameters for such open loops include the feed water flow rate and its enthalpy (or temperature), whereas for the closed loops the heat sink conditions (pressure if it is a steam generator, coolant flow rate, and its inlet temperature) also become additional operating conditions.

5.3.2.1 Effect of power The effect of power on the two-phase NC flow rate is shown in Fig. 5.11A and B. The results given in Fig. 5.11A were generated for the open loop

210

Single-phase, Two-phase and Supercritical Natural Circulation Systems 0.15

30 bar 70 bar 100 bar

Flow rate Exit quality

0.6

Mass flow rate (kg/s)

Flow rate (kg/s), Exit quality (fraction)

0.8

0.4

0.2

0.0 0.0

6.0x10

4

1.2x10

5

1.8x10

Power (W)

(A)

5

2.4x10

5

3.0x10

5

0.10

0.05 o

ID = 14 mm, 15 C subcooling 0.00

0

10

20

30

40

50

60

70

Power (kW) (B)

FIGURE 5.11 Effect of power on flow rate and exit quality. (A) Flow rate and exit quality variation; (B) effect of pressure on flow rate.

shown in Fig. 5.3A, whereas the results given in Fig. 5.11B were generated for a uniform diameter closed loop of the type shown in Fig. 5.1. As expected, the flow rate increases with power in the buoyancy-dominant regime and it decreases with power in the friction-dominant regime. The exit quality is found to increase monotonously with power.

5.3.2.2 Effect of pressure It may be noted from Fig. 5.11B that the pressure has contrasting effects in the buoyancy-dominant and friction-dominant regimes. In the buoyancydominant regime increasing the pressure tends to decrease the flow rate while in the friction-dominant regime increasing the pressure tends to increase the flow rate. In the buoyancy-dominant region, increasing pressure tends to decrease the void fraction and hence the buoyancy force, leading to a decrease in the flow rate. In the friction-dominant regime, increasing the pressure decreases the void fraction increasing the two-phase density thereby decreasing the friction pressure drop leading to an increase in the flow rate. The results given in Fig. 5.12A were generated for the open loop shown in Fig. 5.3A with 19.9 mm loop diameter. The results given are similar to what is given in Fig. 5.11B for the closed loop. The effect of pressure on two-phase NC flow was studied experimentally and theoretically in the 49.3 mm inside diameter loop (see Fig. 5.3B) for a constant power and the results are shown in Fig. 5.12A. The trend of the data is in agreement with the predictions of Eq. (5.25). Initially increasing the pressure increases the two-phase NC flow rate. But at high pressures there is a significant decrease in the void fraction and hence the buoyancy force, leading to a decreasing flow rate with increasing pressure.

Steady-state and transient analysis Chapter | 5 2.0

Mass flow rate (kg/s)

0.45

Flow rate (kg/s)

0.36

0.27 50 bar 100 bar 150 bar

0.18

0.09

0.00 0.0

211

1.2 0.8

4

1.2x10

5

1.8x10

5

2.4x10

5

3.0x10

Eq. (5.31) Experimental data

0.4 0.0

6.0x10

Power: 30 kW

1.6

5

Power (W)

(A)

0

20

40

60

80

Pressure (bar)

(B)

FIGURE 5.12 Effect of pressure on two-phase natural circulation flow rate. (A) Flow rate versus power at different pressures; (B) effect of pressure at constant power.

5.3.2.3 Effect of feed water enthalpy Inlet subcooling in closed-loop two-phase NCSs is not an independent operating parameter. The independent parameter which can be controlled by the operator is the feed water enthalpy which will influence the inlet subcooling. The effect of feed water enthalpy on exit quality and flow rate for two-phase NCSs is shown in Fig. 5.13. These results were generated for the 19.86 mm inside diameter loop with the length scales shown in Fig. 5.3A. Increasing the feed water enthalpy increases the exit quality consistently. However, the effect of feed water enthalpy is found to be different at low and high powers. Increasing the feed water enthalpy in the buoyancy-dominant region (Fig. 5.13C) increases the flow rate, while in the friction-dominant region the flow rate decreases with an increase in feed water enthalpy (see Fig. 5.13D).

5.3.3

Effect of friction multiplier model

Apart from the homogeneous model, a large number of two-phase friction multiplier models exist. These models have a significant effect on the predicted flow rate. Fig. 5.14 gives the predicted two-phase NC flow rate as a function of power for the loop shown in Fig. 5.3A with inside diameter of 9.1 mm. The homogeneous model and the two-phase friction multiplier-based models of Martinelli-Nelson (1948), Lockhart-Martinelli (1949), Baroczy (1966), and Sekoguchi et al. (1970) were selected for this comparison. It is seen from Fig. 5.14 that for low quality (less than 0.9%) and low power, all two-phase friction factor models give practically the same flow rate. Beyond 0.9% quality, widely different trends in flow rate are predicted by the different friction factor models. In general, the homogeneous model predicts the trends correctly although there could be some deviations in the predicted magnitude.

212

Single-phase, Two-phase and Supercritical Natural Circulation Systems 0.32 10 kW 30 kW 50 kW

Exit quality (fraction)

Exit quality (fraction)

0.08

0.06

0.04

100 kW 125 kW 150 kW

0.24

0.16

0.08

0.02

3.50x10

5

7.00x10

5

1.05x10

6

1.40x10

6

3.50x10

5

7.00x10

5

1.05x10

6

1.40x10

6

Feed water enthalpy (J/kg)

Feed water enthalpy (J/kg)

(A)

(B) 0.400

0.40

Flow rate (kg/s)

Flow rate (kg/s)

0.375 0.35

50 kW 30 kW 10 kW

0.30

0.350

100 kW 125 kW 150 kW

0.325

0.25 3.50x10

5

7.00x10

5

1.05x10

6

1.40x10

0.300

6

3.50x10

5

7.00x10

5

1.05x10

6

1.40x10

6

Feed water enthalpy (J/kg)

Feed water enthalpy (J/kg)

(C)

(D)

FIGURE 5.13 Effect of feed water enthalpy on two-phase natural circulation flow rate. (A) Variation of exit quality at low power; (B) variation of exit quality at high power; (C) variation of flow rate at low power; (D) variation of flow rate at high power.

0.05

Homogeneous Baroczy Sekoguchi

Martinelli-Nelson Lockhart-Martinelli

0.03 0.02

Quality = 0.9%

Flow rate (kg/s)

0.04

FIGURE 5.14 Effect of twophase friction factor model.

0.01 0.00

0

2

Loop ID = 9.1 mm Pressure = 60 bar

4

6

8

10

Power (kW)

5.3.4

Effect of the acceleration term

Predictions for the two-phase NC flow rate were made with and without the acceleration term to study its influence. These predictions were obtained for

Steady-state and transient analysis Chapter | 5

213

0.5

Flow rate (kg/s)

0.4

0.3

0.2

With acceleration term Without acceleration term

0.1

0.0 0.0

6.0x10

4

1.2x10

5

1.8x10

5

2.4x10

5

3.0x10

5

Power (W) FIGURE 5.15 Effect of acceleration pressure drop on the predicted natural circulation flow rate.

the 19.9 mm loop with the length scales as shown in Fig. 5.3A for a system pressure of 70 bar. The effect of an acceleration pressure drop on the steadystate mass flow rate is shown in Fig. 5.15. The effect of the acceleration pressure drop is found to be significant at high powers in contrast to singlephase loops where its effect is negligible.

5.3.5

Effect of loop geometry

A large number of geometric parameters affect the two-phase NC flow rate. Important among these are the orientation of the source and sink, riser height (elevation), loop diameter, the inlet and exit orificing, etc. The effect of loop diameter and orificing only are presented here.

5.3.5.1 Effect of loop diameter The effect of loop diameter is predicted for the loop shown in Fig. 5.3A and the results are shown in Fig. 5.16. As seen from Fig. 5.16A, the effect of loop diameter looks more or less similar to that of single-phase loops at low powers. As seen from Fig. 5.16B, the gravity-dominant and frictiondominant regimes appear for high powers. 5.3.5.2 Effect of orificing Orificing has a significant location-specific influence on the predicted flow rate. The effect of orificing has been studied for the 19.9 mm loop with the length scales as shown in Fig. 5.3A. As expected, the inlet orificing (i.e., loss coefficient Ki in the single-phase region at the heater inlet) has

214

Single-phase, Two-phase and Supercritical Natural Circulation Systems FIGURE 5.16 Effect of loop diameter at medium pressure for a two-phase natural circulation loop. (A) Effect of loop diameter at low powers; (B) effect of loop diameter at high powers.

1.2

Mass flow rate (kg/s)

Pressure = 70 bar h feed = 572.1327 kJ/kg Elevation = 2.445 m

0.9

0.6

15.7 mm loop 19.9 mm loop 49.3 mm loop

0.3

0.0

0

1250

2500

3750

5000

Power (W) (A)

Mass flow rate (kg/s)

0.700

0.525

14 mm ID 20.7 mm ID 28 mm ID

0.350

0.175

0.000

o

Predictions at 70 bar and 15 C subcooling

0

75

150

225

300

Power (kW) (B)

much less influence on the flow rate compared to exit orificing (i.e., loss coefficient Ke in the two-phase region at the heater exit) as shown in Fig. 5.17. These calculations were carried out with a constant orifice loss coefficient of 30. Also, the predicted difference in flow rates widens with increase in power (Fig. 5.17).

5.3.6

Closure

It may be noted that all major parameters which influence the two-phase NC flow are studied here. However, the influences of some of the important

Steady-state and transient analysis Chapter | 5

215

0.40

Flow rate (kg/s)

0.32

0.24

Ki = 0 & Ke = 0 Ki = 30 & K e = 0 Ki = 0 & Ke = 30

0.16

0.08

0.00 0.0

6.0x10

4

1.2x10

5

1.8x10

5

2.4x10

5

3.0x10

5

Power (W)

FIGURE 5.17 Effect of inlet and exit orificing.

parameters such as the effect of elevation difference and the orientation of the source and sink are excluded as their effect is similar to that of single-phase NC.

5.3.7

Effect of parallel channels

Like a single-phase parallel channel system, two-phase parallel channel systems also exhibit multiple steady-state solutions with differing flow directions. Experiments in a parallel channel system with 10 channels showed that flow reversal can be induced by switching off the power in any channel. Linzer and Walter (2003) showed that the flow reversal in a two-phase parallel channel NCS can occur without reducing the power to zero. Based on their study, a criterion for flow reversal has been given by them for a threechannel system with two unequally heated parallel channels. The methodology for the analysis of parallel channel systems is similar to that discussed in Chapter 4 except that two-phase conditions are encountered. Hence, the effect of parallel channels is not discussed again here.

5.4 Transient performance of two-phase natural circulation systems For the transient case, the time-dependent equations are to be solved numerically as in single-phase loops. However, in two-phase NCSs the fluid density in two-phase regions is a function of phase fraction and the saturation values for the vapor and liquid phases. In single-phase regions, it is a function of pressure and fluid temperature. Similarly, the fluid enthalpy in two-phase

216

Single-phase, Two-phase and Supercritical Natural Circulation Systems

regions is also a function of phase fraction and saturation values for the vapor and liquid phases, while in single-phase regions it is a function of pressure and fluid temperature. Moreover, the estimation of two-phase pressure loss also involves the specification of a two-phase friction multiplier model. Earlier in this chapter, it was pointed out that the homogeneous equilibrium model which assumes equal velocity and equal temperature for both phases predicts the NCS behavior reasonably well. As already stated, although several advanced models like four-equation, five-equation, and six-equation models are available for predicting the behavior of two-phase systems, a detailed discussion of these advanced models is beyond the scope of this book and our discussion is restricted to the homogeneous equilibrium model. Naveen et al. (2017a) developed a system code NUTAN-Th for thermal hydraulic analysis of multichannel NCSs based on the homogeneous equilibrium model. The code has previously been used for analysis of single-phase NCLs (Naveen et al., 2014). For the homogeneous model, the transient mass and momentum balance equations used are the same as those described by Eqs. (5.1) and (5.2). The energy balance equation is somewhat modified as   @ 1 @ðwiÞ @P qξ 1 @ Ak @T @s ðρiÞ 1 2 5 1 ð5:64Þ @t A @s @t A A @s In the energy balance equation, the last term represents the heat transfer due to axial conduction in the fluid. Fluid axial conduction may not play a large role in two-phase flow, however, accounting this term does not introduce any numerical problems. Also, it has some minor effects during start-up from stagnant single-phase initial conditions especially for liquid metals. Since all the NCLs need to be started from single-phase conditions, this term has been included in the energy balance equation. In many closed NCLs, the heater consists of an electrically heated pipe wherein the heat is generated in the pipe wall itself. At the inner wall, heat is exchanged between the pipe and the fluid flowing through the heater. This heat transfer is governed by Newton’s law of cooling given by q 5 hi ðTw 2 T Þ Substituting for q from Eq. (5.65) into Eq. (5.64) yields   @ 1 @ðwiÞ @P hi ξ ðTw 2 T Þ 1 @ Ak @T @s ðρiÞ 1 2 5 1 @t A @s @t A A @s

ð5:65Þ

ð5:66Þ

In addition the equation of state is also needed, which can be written as ρ 5 ρðP; iÞ

ð5:67Þ

Eqs. (5.1), (5.2), and (5.66) form a nonlinear hyperbolic system of partial differential equations.

Steady-state and transient analysis Chapter | 5

FIGURE 5.18 Mesh cell labeling convention (Liles and Reed, 1978).

j+2 lk+2

k+2 k+1 k k-1 zk

zk-1

217

j+1 j

lk+1 lk j-1 j-2

k-2 Reference (i.e., z = 0)

5.4.1

Discretization of governing equations

For the numerical solution of the governing equations, a staggered mesh approach is used. The flow field is divided into nonoverlapping rigid control volumes called cells. The field variables, i, P, and ρ are defined at the cell centers, while the fluid mass flow rate, w, is defined at cell interfaces. The mesh cell configuration and labeling conventions for cell edges and cell centers are depicted in Fig. 5.18. The mass and energy equations are discretized over the mesh cells indicated by solid lines in Fig. 5.18 and the momentum balance equation is integrated over the dashed mesh cells. This forms a staggered spatial difference scheme. Discretization is done making use of a semi-implicit numerical technique. Eq. (5.1) is discretized as

n11 n11 w 2 w j j21 Δt ρn11 5 ρnk 2 ð5:68Þ k Ak lk The momentum conservation equation is discretized as

n n11 n 2Pn11 1 Pn11 ðj 5 1 to m 2 1Þ k Δt 1 kj 1 Cj wj k11 Δt 5 Oj ð5:69Þ where

1 2 0 !n !n 3 2 2 2 φ2LO 1 2 ϕ ϕ φ j flφ flφ j LO n LO LO A1 5ðjwjÞn Δt Cj 5 4Kj@ 1 1 j 4Dh A2 ρf 4Dh A2 ρf 2ρf A2k 2ρf A2k11 k

k11

ð5:70Þ

218

Single-phase, Two-phase and Supercritical Natural Circulation Systems

lk lk11 1 and ð5:71Þ 2Ak 2Ak11 2 3
2 n


2 n
2 n  w w w 1 1 n n 4 5 Oj 5 Ij wj 2 Δt 1 2 2 2 Δt ρA2 k11 ρA2 k 2ρ j A2k11 Ak      2 ρnk11 zk11 2 zj 1 ρnk zj 2 zk gΔt kj 5

ð5:72Þ   n For local losses the term jwj=2ρA2 j is evaluated at the smaller of the cross-section areas on the two sides of the junction. Hence,   1; if Ak =Ak11 , 1 ð5:73Þ ϕj 5 0; if Ak =Ak11 $ 1 The energy conservation equation is discretized as





n11 inj wn11 2 inj21 wn11 hni;k ξ k Tw;k 2 Tkn Δt   j j21 Δt n11 ρ n11 5 ρnk ink 2 1 Pn11 2 Pnk 1 k ik k Ak lk Ak 2 0 1 13 0 1

6 B n C C7 B n nC n C7 B B 2Δt 6 6Aj BTk11 2 Tk C 2 Aj21 BTk 2 Tk21 C7 C C7 B B lk Ak lk 6 4 @ lk11 1 lk A @ 1 lk21 A5 kk11 kk kk kk21 ð5:74Þ

Eqs. (5.68), (5.69), and (5.74) have four unknowns, namely, cell enthalpy, cell pressure, cell density, and junction mass flow rate. A solution of the governing equations in the time domain requires simultaneous solution of Eqs. (5.68), (5.69), and (5.74) and the equation of state (i.e., Eq. 5.67). In the present analysis, this is achieved by eliminating enthalpy from Eq. (5.74) using the equation of state. From the equation of state, the rate of change of density is written as

 @ρ @ρ @i @ρ @P 5 1 ð5:75Þ @t @i P @t @P i @t Using Eq. (5.75), Eq. (5.1) is written as

 @ρ @i @ρ @P 1 @w 1 1 50 @i P @t @P i @t A @s

ð5:76Þ

Integrating Eq. (5.76) over the control volume “k,” ð t1Δt ð sj
 ð t1Δt ð sj ð t1Δt ð sj
 @ρ @i @ρ @P 1 @w dsdt1 dsdt1 dsdt 5 0 P @i @t @P @t A @s sj21 sj21 sj21 t t t i ð5:77Þ

Steady-state and transient analysis Chapter | 5

219

 

 @ρ n  n11 n  @ρ n  n11 1 n11 l i 2 i 1 2 Pnk lk 1 wj 2 wn11 P;k k k k i;k Pk j21 Δt 5 0 ð5:78Þ   @i @P Ak Rearranging gives in11 2 ink k 52 Δt

1 n @ρ  @i 

wn11 j

2 wn11 j21

Ak lk

P;k

2

n @ρ  @P i;k n @ρ  @i P;k



Pn11 2 Pnk k Δt

 ð5:79Þ

The energy balance equation, that is, Eq. (5.66) is written in the form as follows: i


 @ρ @i 1 @ðwiÞ @P hi ξðTw 2 T Þ 1 @ @T 1ρ 1 2 5 1 Ak @t @t A @s @t A A @s @s

ð5:80Þ

Substituting for @ρ=@t from Eq. (5.1) into Eq. (5.80),
 i @w @i 1 @ðwiÞ @P hi ξ ðTw 2 T Þ 1 @ @T 1ρ 1 2 5 1 Ak 2 A @s @t A @s @t A A @s @s

ð5:81Þ

Integrating Eq. (5.81) over control volume “k,” gives Ð t1Δt Ð sj

@i dsdt 1 @t

ð t1Δt ð sj

1 @ðwiÞ dsdt t sj21 A @s ð t1Δt ð sj Ð t1Δt Ð sj i @w @P 2 t dsdt 5 sj21 A @s dsdt 2 sj21 @t t
 ð t1Δt ð sj Ð t1Δt Ð sj hi ξ ðTw 2 T Þ 1 @ @T dsdt 1 Ak dsdt sj21 t A @s t sj21 A @s t

sj21 ρ





  n11 n n11 n n n11 n11 n w i 2 w i i w 2 w n in11 j j j21 j21 k j j21 2 ik 1 ρk k 2 A k lk Ak lk Δt

 n11  n n11 n n h ξ T 2 T k i;k w;k k P 2 Pk 5 2 k 1 Ak Δt 2 0 1 0 13 B n B T n 2 T n C7 n C 1 6 6 B Tk11 2 Tk C B k21 C7 6Aj B C 2 Aj21 B k C7 lk A lk21 A5 @ lk Ak lk 4 @ lk11 1 1 2kk11 2kk 2kk 2kk21

ð5:82Þ

ð5:83Þ

220

Single-phase, Two-phase and Supercritical Natural Circulation Systems





   n11  n n wn11 inj 2 ink 2 wn11 inj21 2 ink n in11 j j21 2 i P 2 P k k k k ρk 2 1 Ak lk Δt Δt 2 0 1 0 13

n n11 n B n B T n 2 T n C7 n C hi;k ξk Tw;k 2 Tk 1 6 6 B Tk11 2 Tk C B k21 C7 1 5 6 Aj B C 2 Aj21 B k C7 lk A lk21 A5 @ lk Ak lk 4 @ lk11 Ak 1 1 2kk11 2kk 2kk 2kk21 ð5:84Þ Substituting for in11 from Eq. (5.79) into Eq. (5.84) gives, k "

n  @ρ  n n ij 2 ik 2 ρnk wn11 j  @i P;k "  # n n 
 n

 @ρ  n n11 Ak lk @ρ n @ρ  n n   Pn11  5 ij21 2 ik 2 ρ k wj21 1 1ρ k    @i P;k @P Δt @i P;k k i;k n  !

A l @ρ n n @ρ n @ρ  n k k n11 n    1 h ξ lk Tw;k 2 Tk 2 1ρ k Pn @i P;k i;k k @i P;k @P i;k k Δt 2 0 1 0 13 n B n B T n 2 T n C7 nC @ρ  6 6 BT 2 Tk C B k21 C7 1 2  6Aj B k11 C 2 Aj21 B k C7 lk A lk21 A5 @ lk @i P;k 4 @ lk11 1 1 kk11 kk kk kk21 ð5:85Þ

For convenience, Eq. (5.85) is written in the following form n n11 n11 Γnj;k wn11 5 Γnj21;k wn11 1 Γnw;k Tw;k 1 Θnk j j21 1 Ωk Pk

where

n  @ρ  n n 5 ij 2 ik 2 ρnk  @i P;k n

 @ρ  n ij21 2 ink 2 ρnk Γnj21;k 5  @i P;k Γnj;k

ð5:86Þ

ð5:87Þ ð5:88Þ

Steady-state and transient analysis Chapter | 5

Ωnk

Ak l k 5 Δt

 ! n n @ρ n @ρ   1ρ k  @i P;k @P i;k

@ρ n n h ξ lk @i P;k i;k i;k n n n !     @ρ A l @ρ @ρ n k k  hn ξ lk T n 2  1 ρ  Pnk and Θnk 5 2 k @P  @i P;k i;k i;k k @i P;k Δt i;k 2 0 1 0 13  B n B T n 2 T n C7 nC @ρ n 6 6 BT 2 Tk C B k21 C7 1 2  6Aj B k11 C 2 Aj21 B k C7 lk A lk21 A5 @ lk @i P;k 4 @ lk11 1 1 kk11 kk kk kk21 Γnw;k 5

221

ð5:89Þ ð5:90Þ

ð5:91Þ

The advantage of presenting the energy conservation equation in the form given by Eq. (5.86) is that this equation has enthalpy terms of old time step which are known. This allows for the solution of Eqs. (5.68), (5.69), and (5.86) using the tridiagonal matrix algorithm (TDMA). The solution of discretized energy balance equations, that is, Eqs. (5.74) and (5.86), needs fluid enthalpies at cell junctions. However, in the staggered mesh approach adopted here, the fluid enthalpy is not defined at the cell interfaces or junctions and, hence, needs to be evaluated in terms of known quantities. This can be achieved using the upwind scheme. As per this scheme, junction enthalpy is given by

 1 1 ηj n 1 2 ηj n n ij 5 ð5:92Þ ik 1 ik11 2 2 where

  ηj 5 wnj =wnj

ð5:93Þ

5.5 Solution procedure The solution of governing equations in time domain involves the following steps: 1. Eqs. (5.70)
(5.72) and (5.87)
(5.91) are evaluated first. These equations give the coefficients in Eqs. (5.69) and (5.86), respectively. The system of equations given by Eqs. (5.69) and (5.86) (in the case of flow through a channel, pressures at channel inlet and exit are given) constitute “2m 1 1” equations in “2m 1 1” unknowns. The unknowns are “m 1 1” junction flow rates and “m” cell pressures. The solution for this system of equations is obtained using TDMA.

222

Single-phase, Two-phase and Supercritical Natural Circulation Systems

2. Having obtained the junction mass flow rate for all the junctions at the new time step, Eq. (5.68) is solved for density in each cell at the new time step. Finally, the enthalpy in each cell at the new time step is obtained from the solution of Eq. (5.74) for each cell.

5.6 Typical transient solutions To illustrate the use of the above-described solution methodology, a few typical transients relevant to two-phase NC are solved here.

5.6.1

Start-up from rest

Unlike forced circulation systems where flow can be established before the heater is switched on, a NCS needs to be started from the state of stagnant initial conditions. A two-phase NCL passes through the state of single-phase NC before it reaches the state of a two-phase NCL. Naveen et al. (2010) studied the start-up of NCSs and studied the effect of loop diameter and operating conditions on their start-up. Fig. 5.19 shows the start-up predictions made using NUTAN-Th code for the loop shown in Fig. 5.3B from stagnant single-phase initial conditions at a pressure of 5 bar. At time t 5 0, the heater power is suddenly raised to 20 kW. To start with the single-phase NC sets in the loop till the fluid temperature reaches the saturation temperature. With the initiation of boiling, the buoyancy force increases significantly, which leads to a sharp rise in the loop mass flow rate. With the increased mass flow rate the heater exit enthalpy becomes less than the saturation value and as a result the system returns back to single-phase NC, resulting in a sharp reduction in flow. The reduced flow rate increases the exit enthalpy beyond the saturation value, resuming steam production and resulting in a sharp increase in the flow rate and the system returns back to single-phase again. After a few oscillations, steady two-phase NC is achieved. The inset in Fig. 5.19A shows the mass flow rate oscillations during the transition from single-phase natural circulation (SPNC) to two-phase natural circulation (TPNC). It may be noted that, although the riser is adiabatic, its steam quality is slightly larger (see Fig. 5.19B) because of the flashing induced by the static pressure decrease as the fluid rises up. It may also be noted that the quality appears first in the riser exit and then in the heater exit, which is again attributed to flashing. These oscillations could induce significant vibrations in the case of loops with very tall risers. Fig. 5.20 shows predictions for the start-up and power raising in the 15.7 mm inside diameter loop shown in Fig. 5.3A from stagnant single-phase initial conditions at 25 bar. At time t 5 0, the heater power is raised to 2.5 kW and is gradually increased in steps to 6 kW. To start with, the singlephase NC sets in the loop till the fluid temperature reaches the saturation value. With the initiation of boiling, the buoyancy force increases

Steady-state and transient analysis Chapter | 5

Pressure: 5 bar Power: 20 kW 2.5

2.0

2.0 Mass flow rate (kg/s)

Flow rate (kg/s)

2.5

1.5

223

TPNC

1.5

1.0

0.5

1.0 0.0 1400

1450

1500

1550

Time (s)

0.5 SPNC

0.0

1000

2000

Time (s) (A) 0.008

Steam quality (fraction)

Heater exit Riser exit

0.006

0.004

0.002

0.000 1400

1450

1500

1550

Time (s) (B) FIGURE 5.19 Start-up transient at 5 bar in the loop shown in Fig. 5.3B due to a sudden increase of power. (A) Flow variation; (B) quality variation.

significantly, which leads to a sharp rise in loop mass flow rate. In the gravity-dominant region, as described earlier, flow increases with an increase in power. Note that the flow rate oscillations are significantly reduced in this case due to the higher pressure.

5.7 Parametric influences In a NCS, system thermal-hydraulic behavior is highly dependent on system geometry and operating conditions. System operating parameters which

Single-phase, Two-phase and Supercritical Natural Circulation Systems

0.20

7000

Flow Mass flow rate (kg/s)

FIGURE 5.20 Predicted start-up behavior for 15.7 mm ID loop at 25 bar pressure.

8000 Pressure = 25 bar

0.15

6000 5000

0.10

4000 3000

Power 0.05

Power (W)

224

2000 1000

0.00 0

1000

2000

3000

4000

5000

0 6000

Time (s)

influence the transient behavior include power, pressure, inlet subcooling, feed water flow rate, and inlet temperature. The effects of some of these parameters are studied in this section.

5.7.1

Effect of initial pressure

The flow rate in a two-phase NCS is dependent on the operating conditions such as system pressure and heater power. This can easily be understood from the variation of steam void fraction with steam quality given in Fig. 5.8. A careful analysis of Fig. 5.8 shows that for the same value of steam quality, the void fraction value varies significantly with system pressure. At low pressure, a steam quality of 1% can result in almost 95% void fraction whereas the same quality near critical pressure leads to an almost equal value of void fraction. It was explained earlier that the driving force in a NCS depends upon the density difference and height difference between the sink and source. Therefore the larger the void fraction in the hot leg, the larger is the driving force. This effect becomes much more pronounced for systems having tall chimneys/risers. Fig. 5.21A and B show the start-up simulation for the 15.7 mm inside diameter loop with the length scales as given in Fig. 5.3A at pressures of 3 and 5 bar, respectively. It may be noted that the amplitude of oscillations during the transition from single-phase to two-phase NC decrease with increase in pressure whereas the steady-state flow at 5 bar is more than that at 3 bar indicating that both transients are in the friction dominant region. This becomes clear if we compare Fig. 5.21 with Fig. 5.12.

5.7.2

Effect of loop diameter

Fig. 5.22A and B shows the simulations for 25 bar pressure for 9.1 and 15.7 mm inside diameter loops, respectively. The length scales of the loops

Steady-state and transient analysis Chapter | 5 (A)

FIGURE 5.21 Model predictions for 15.7 mm ID loop at 250 W and 3 and 5 bar. (A) Start-up at 3 bar; (B) start-up at 5 bar.

0.20 Power = 250 W 0.15

Flow rate (kg/s)

225

0.10

0.05

0.00 15,000 15,500 16,000

16,500 17,000 17,500 Time (s)

18,000

(B) 0.20 Power = 250 W

Flow rate (kg/s)

0.15

0.10

0.05

0.00 12,000

12,500

13,000

13,500

14,000

14,500

15,000

Time (s)

correspond to that shown in Fig. 5.3A. A comparison of Fig. 5.22A and B shows that for the same operating pressure, the loops having lower diameter are more oscillatory. Increasing the loop diameter significantly reduces the amplitude of two-phase flow oscillations. Also as expected, the larger diameter loop shows significantly larger flows. Also, in the case of the 9.1 mm diameter loop, the peak flow rate is obtained at about 1.8 kW and subsequently the friction-dominant region sets in, whereas for the 15.7 mm loop the gravity-dominant region is observed even at 6 kW.

5.7.3

Effect of parallel channels on the transient behavior

Parallel heated channels are present in almost all nuclear as well as the fossil-fueled power plants. Hence, the study of parallel channel effects is important for transients in these plants. Parallel channels can show interesting behavior during the transients, including flow reversals. For simulating

Single-phase, Two-phase and Supercritical Natural Circulation Systems

Mass flow rate (kg/s)

(A)

9.1 mm ID

Pressure = 25 bar

0.07

8000

0.06

7000

0.05

6000

0.04

5000

0.03

4000

0.02

3000

0.01

2000

0.00

1000

–0.01 0

Power (W)

226

FIGURE 5.22 Model predictions for 9.1 mm ID loop at 25 bar pressure. (A) Start-up and power raising in 9.1 mm loop; (B) start-up and power raising in 15.7 mm loop.

0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000 5500 6000 Time (s)

(B) 0.20

8000 Pressure = 25 bar 7000

Flow

6000 5000

0.10

4000

Power (W)

Mass flow rate (kg/s)

0.15

3000

Power 0.05

2000 1000

0.00

0

1000

2000

3000

4000

5000

0 6000

Time (s)

parallel channel behavior in the time domain the governing equations remain the same as described earlier. However, most of the commercial computer codes are capable of handling only a few parallel channels during transient analysis, whereas NC-based boiling water reactors and fossil fuel plants have hundreds of parallel channels. For example, a typical configuration of the NC-based vertical pressure tube-type boiling water reactor, such as AHWR, has 452 parallel channels. Naveen et al. (2017a) developed a computer code, NUTAN-Th, to simulate the behavior of multichannel systems. As described earlier, the code is based on the homogeneous equilibrium model. In the code, the governing equations are discretized for each channel as described in Section 5.4 and pressure is assumed to be uniform in the inlet and outlet plena (headers). This implies that we need to solve the discretized equations for all the channels simultaneously. Naveen et al. (2017b) carried out cold start-up simulation of AHWR using computer code, NUTAN-Th. A typical configuration of the main heat transport system of AHWR considered for the analysis is shown in Fig. 5.23A and B.

227

Steady-state and transient analysis Chapter | 5 Steam

(A)

(B) Steam to turbine

Steam drum

Feed water

Steam drum

Feed water in

Downcomer Tail pipe Down comer

Distribution ring header Calandria housing active core

Tail pipes

Calandria

Ring header

Feeders

Feeder

FIGURE 5.23 Schematic of the pressure tube-type reactor. (A) Schematic of the MHTS; (B) schematic of one loop.

It has four loops (only two are shown in Fig. 5.23A), each with a steam drum from where four downcomers take off (only one is shown in Fig. 5.23A and B) and join a ring header which is common to all four loops. From the header, the feeders take off and join to an equal number of coolant channels where the subcooled water absorbs heat and becomes a two-phase mixture. The two-phase mixture rises through the tail pipes to the steam drum where it gets separated and the separated steam goes out. The separated water mixes with the subcooled feed water and enters the downcomers, thus completing the circuit. The density difference between the subcooled water in the downcomer and two-phase mixture in the tail pipe provides the driving force for the NC flow. Each loop has 113 parallel channels (only three are shown in Fig. 5.23B) connected between the ring header and the steam drum and four downcomers connecting the steam drum to the inlet header. Each loop serves one quarter of the core. Several peripheral channels are orificed and the orifice loss coefficients used in the calculation are shown in Table 5.2. As can be seen, most of the central channels are not orificed (loss coefficient equals 0). The channel power distribution at FP is shown in Table 5.3. The table also shows how the channels are named. A transient simulation was carried out considering all 113 channels in a loop. For this simulation, a cold start-up at 2% FP and 70 bar was simulated from stagnant initial conditions assuming the power distribution to be the same as at full power. Hence to obtain the channel powers at 2% FP, the value given for each channel in Table 5.3 is divided by 50. The predicted transients in all the channels are similar. Hence, typical predicted transients in the 10 channels of the K-row (the channels inside the dotted line in Table 5.2) are shown in Fig. 5.24A and B. It may be noted that

228

Single-phase, Two-phase and Supercritical Natural Circulation Systems

TABLE 5.2 Map for the inlet orifice loss coefficient.

13

12 14

11 15

10 16

9 17

8 18

Z

A

300

300

300

300

300

300

Y

B

SOR

300

300

300

300

300

7 19

6 20

5 21

4 22

3 23

2 24

1 25

300

X

C

0.0

0.0

0.0

SOR

0.0

0.0

0.0

0.0

W

D

AR

0.0

0.0

0.0

0.0

0.0

SOR

0.0

0.0

V

E

0.0

0.0

0.0

0.0

RR

0.0

0.0

0.0

0.0

0.0

U

F

0.0

0.0

SOR

0.0

0.0

0.0

0.0

SOR

0.0

0.0

0.0

T

G

SR

0.0

0.0

0.0

0.0

0.0

0.0

0.0

0.0

SOR

0.0

300

S

H

0.0

0.0

0.0

0.0

0.0

AR

0.0

0.0

0.0

0.0

0.0

300

R

J

0.0

0.0

0.0

0.0

0.0

0.0

0.0

0.0

RR

0.0

0.0

300

300

Q

K

SOR

0.0

0.0

SR

0.0

0.0

0.0

0.0

0.0

0.0

SOR

300

300

P O

300

0.0

0.0

0.0

0.0

0.0

0.0

0.0

SOR

0.0

0.0

0.0

300

300

M

0.0

0.0

0.0

0.0

0.0

0.0

0.0

0.0

0.0

0.0

0.0

300

300

N

SOR

0.0

0.0

SOR

0.0

0.0

SR

0.0

0.0

AR

0.0

SOR

300

L

Note: The symmetry lines go through the middle of channels in the N-th row and 13th column. RR, regulating rod; SR, shim rod; SOR, shut-off rod; AR, adjuster rod.

in all channels the NC flow initiates in single-phase condition. However, the flow rates are not equal as the power in each channel is different. The flow initiation in single-phase condition is accompanied by oscillations. The oscillations persist for a longer time in low-power channels. Subsequently, a stable marginally increasing single-phase flow is observed during the heat-up phase. Once the saturation values are achieved, a sharp jump in flow rate with a few oscillations is observed. The transition from single-phase to twophase NC is shown in detail in Fig. 5.24B. Subsequently, a stable steady state with flow rate more than thrice that of the initial single-phase NC is attained with two-phase NC. Fig. 5.25A shows the predicted downcomer flow rate variation during cold start-up at 2% FP. The flow rate variations of all four downcomers are almost identical. Although the downcomer flow rate remains entirely in single phase, a jump in flow rate is observed as the boiling occurs in the core outlet. Fig. 5.25B shows the sharp rise in flow during the phase change in the channels during the transient.

Steady-state and transient analysis Chapter | 5

229

TABLE 5.3 Typical channel power distribution in the AHWR core.

13

Z

A

Y

B

X

C

W

D

V

F

T

G

R

10

9

8

7

6

5

4

3

2

1

15

16

17

18

19

20

21

22

23

24

25

1.842

1.707

1.649

1.628

1.635

1.655

SOR

1.842

1.75

1.908

1.732

1.658

1.859

1.84

1.782

2.047

SOR

2.07

1.897

2.038

2.094

AR

1.779

1.977

2.163

2.045

2.163

SOR

2.147

2.01

1.855

1.979

2.224

2.103

RR

2.086

2.308

2.301

2.01

2.008

2.187

2.249

SOR

2.259

2.015

1.966

2.262

SOR

2.298

2.142

2.086

SR

2.333

2.254

2.109

1.914

1.845

1.984

2.259

2.303

SOR

2.028

1.846

2.338

2.152

2.064

2.076

1.895

AR

1.843

1.962

2.079

2.153

1.886

1.645

1.63

2.391

2.195

2.156

2.256

2.051

1.893

1.911

2.009

RR

2.034

2.055

1.716

1.606

H J

Q

K

P

L

O

11

14

E

U

S

12

SOR

2.422

2.383

SR

2.254

2.073

2.104

2.251

2.093

2.149

SOR

1.888

1.598

2.461

2.267

2.24

2.382

2.153

2.06

2.247

SOR

2.213

1.963

2.028

1.731

1.618

2.487

2.289

2.266

2.42

2.192

2.147

2.326

2.24

1.959

1.766

1.766

1.819

1.673

SOR

2.487

2.461

SOR

2.391

2.338

SR

2.187

1.855

AR

1.84

SOR

1.842

M N

Note: The symmetry lines go through the middle of channels in the N-th row and 13th column. Channel powers in MWth. RR, regulating rod; SR, shim rod; SOR, shut-off rod; AR, adjuster rod.

5.8 Closure A methodology for the steady-state analysis has been presented in this chapter along with derivation of explicit flow rate equations for both closed- and open-loop NCSs. A dimensionless equation for flow rate has also been derived. However, an explicit dimensionless flow rate equation can be obtained only for closed loops. For open loops, the dimensionless flow equation is a polynomial. If the acceleration term can be neglected then the explicit dimensionless equations for the closed and open loops become the same. Both the explicit dimensional and dimensionless flow rate equations are based on the homogeneous equilibrium model. The steady-state performance of the homogeneous model when compared with the two-fluid model showed only marginal deviation, which is attributed to the near thermal equilibrium that exists in low-velocity NC flows. In addition, the difference in vapor and liquid velocities is not large enough to cause significant deviation in the

230 (A)

Single-phase, Two-phase and Supercritical Natural Circulation Systems

6 Mass flow rate (kg/s)

FIGURE 5.24 Cold start-up at 2% power for a pressure tube-type natural circulation-based BWR. (A) Start-up transient for various channels in K-row; (B) transition from SPNC to TPNC.

7

K1 K2 K4 K5 K6 K7 K8 K9 K11 K12

5 4 3 2 1 0 0

5000

10,000

15,000

20,000

25,000

30,000

Time (s)

(B)

7 K1 K2 K4 K5 K6 K7 K8 K9 K11 K12

Mass flow rate (kg/s)

6 5 4 3

TPNC

2 SPNC 1 0 29,500

30,000

30,500

31,000

Time (s)

predictions for the loops considered for analysis. However, this may not remain valid for transients involving severe nonequilibrium effects like LOCA with cold water injection. The adequacy of the dimensionless equation derived has been tested with data from both uniform and nonuniform diameter loops which showed deviations within 6 40%. Using the explicit equations derived for the steady-state behavior, a parametric analysis of both open- and closed-loop NCSs has been carried out. The analysis covered the parametric influences of the operating parameters such as power, pressure, feed water enthalpy, as well as the geometric parameters like loop diameter. Unlike single-phase systems, the acceleration term has significant influence on the steady-state performance of two-phase NC systems. The influence of local loss coefficient is found to be location-specific, unlike single-phase NCSs. The local loss coefficient in the two-phase region has a stronger influence on the steady-state flow rate than the local loss coefficient in the single-phase region. An interesting characteristic of two-phase NCLs is the existence of buoyancy-dominant and

Steady-state and transient analysis Chapter | 5

160

70

140

Mass flow rate (kg/s)

60

Downcomer 1 Downcomer 2 Downcomer 3 Downcomer 4

120

50

100 40 80 30 60

Power (% FP)

(A)

231

FIGURE 5.25 Variation of the flow rate in the downcomer during cold start-up. (A) Downcomer flow variation; (B) sharp rise in flow due to the phase change.

20

40 Flow

10

20 Power 0

0

5000

10,000

15,000

20,000

25,000

30,000

0

Time (s) 5 Downcomer 1 Downcomer 2 Downcomer 3 Downcomer 4

Mass flow rate (kg/s)

80

4

60

3

Power

40

2

Power (% FP)

(B)

100

Flow

20

1

0 29,000

29,500

30,000

30,500

0 31,000

Time (s)

friction-dominant regimes which are not normally observed in single-phase closed loops except at very high elevations. Certain parameters like pressure have contradicting influences on the flow rate depending on the NC flow regime. A formulation for the transient analysis of two-phase NCLs has been presented based on the homogeneous equilibrium model. Using the formulation, transient analysis of simple open loops and a pressure tube-type boiling water reactor have been performed. Influences of some important parameters on the transient performance have also been explored. It may be noted that the transient formulation presented for two-phase NCSs is also applicable to single-phase NCSs.

Nomenclature A a

flow area, m2 dimensionless flow area (A/Ar)

232

Single-phase, Two-phase and Supercritical Natural Circulation Systems

b Cp D d f g Grm H h I II i K k L lk l N NG P p Q q R Ress s S Δs t Δt T v v V w x z Z

exponent in the friction factor equation specific heat, J/(kg K) diameter, m dimensionless hydraulic diameter (D/Dr) friction factor acceleration due to gravity, m/s2 modified Grashof number (D3r ρ2f β h gQh H=Ar μ3r Þ height, m heat transfer coefficient,W/m2K  dimensionless enthalpy I 5 ði 2ir Þ=Δi Hh  dimensionless enthalpy integral II 5 IdZ enthalpy, J/kg loss coefficient thermal conductivity, W/(mK) length, m length of kth control volume, m dimensionless length, (L/Lt) number of pipe segments geometric number defined by Eq. (5.33) pressure, Pa constant in friction factor equation power, W heat flux, W/m2 24 hydraulic resistance,  m  Reynolds number, Dr wss =Ar μr coordinate along the loop, m dimensionless coordinate around the loop (s/H) space step, m time, s time step, s temperature, K specific volume, m3/kg dimensionless specific volume, (v/vr) volume, m3 mass flow rate, kg/s quality elevation, m dimensionless elevation (z/H)

Subscripts 0 a B c cl DC e

reference ambient boiling cooler cold leg downcomer exit

Steady-state and transient analysis Chapter | 5 eff eq F f fg h hl i i, in j k o r s SD sp ss t tp tp,c w

effective equivalent feed saturated liquid gas minus liquid heater hot leg inside, ith segment, ith node inlet junction kth control volume outside reference secondary, steam, single-phase steam drum single-phase steady state total two-phase two-phase cooler wall

Superscripts 2 n

mean value nth time step

Greek symbols α βh βh γ Δ φ2LO θ μ ξ ρ

thermal diffusivity, m2/s thermal expansion coefficient based on enthalpy, (J/kg)21 thermal expansion coefficient based on temperature, K21 H parameter defined as Li =Ai , 1/m refers to difference two-phase friction multiplier angle between the flow direction and the horizontal dynamic viscosity, kg/(ms) perimeter, m density, kg/m3

List of acronyms used AHWR AR BWR CANDU DFM ESBWR EVET

advanced heavy water reactor adjustor rod boiling water reactor Canadian deuterium uranium reactor drift flux model economic simplified boiling water reactor equal velocity equal temperature

233

234

Single-phase, Two-phase and Supercritical Natural Circulation Systems

FDR FP GDR ID LHS LOCA MHTS NC NCL NCS PWR RHS RR SBLOCA SD SOR SPNC SR TDMA TFM TPNC UDL UVUT

friction-dominant regime full power gravity-dominant regime inside diameter left-hand side loss of coolant accident main heat transport system natural circulation natural circulation loop natural circulation system pressurized water reactor right-hand side regulating rod small break LOCA steam drum shut-off rod single-phase natural circulation shim rod tridiagonal matrix algorithm two-fluid model two-phase natural circulation uniform diameter loop unequal velocity unequal temperature

References Baroczy, C.J., 1966. A systematic correlation for two-phase pressure drop. Chem. Eng. Progr. Symp. Ser. 62, 232
249. Dimmick, G.R., Chatoorgoon, V., Khartabil, H.F., Duffey, R.B., 2002. Natural convection studies for advanced CANDU reactor concepts. Nucl. Eng. Design 215, 27
38. Dubey, P., Rao, G.S.S.P., Pilkhwal, D.S., Vijayan, P.K., Saha, D., 2004. Analysis of Experimental Data on Two-Phase Natural Circulation from the Flow Pattern Transition Instability Studies Facility at Apsara Reactor. Report No. BARC/2004/E/031, Bhabha Atomic Research Centre, Trombay, Mumbai. Duffey, R.B., 2000. Natural convection and natural circulation flow and limits in advanced reactor concepts, Natural circulation data and methods for advanced water cooled nuclear power plant designs. In: Proceedings of Technical Committee meeting, IAEA-TECDOC-1281, Vienna, pp. 49
65. Fukuda, K., Kobori, T., 1979. Classification of two-phase flow stability by density-wave oscillation model. J. Nucl. Sci. Technol. 16, 95
108. Furutera, M., 1986. Validity of homogeneous flow model for instability analysis. Nucl. Eng. Design 95, 65
74. Gartia, M.R., Vijayan, P.K., Pilkhwal, D.S., 2006. A generalized flow correlation for two-phase natural circulation loops. Nucl. Eng. Design 236, 1800
1809. Ishii, M., Kataoka, I., 1983. Scaling laws for thermal-hydraulic system under single phase and two-phase natural circulation. Nucl. Eng. Design 81, 411
425. Lee, S.Y., Lee, D.W., 1991. Linear analysis of flow instabilities in an open two-phase natural circulation loop. Nucl. Eng. Design 128, 317
330.

Steady-state and transient analysis Chapter | 5

235

Liles, D.R., Reed, Wm. H., 1978. A semi implicit method for two-phase fluid dynamics. J. Comput. Phys. 26, 390
407. Linzer, W., Walter, H., 2003. Flow reversal in natural circulation systems. Appl. Therm. Eng. 23, 2363
2372. Lockhart, R.W., Martinelli, R.C., 1949. Proposed correlation of data for isothermal two-phase, two-component flow in pipes. Chem. Eng. Prog. 45, 39
48. Martinelli, R.C., Nelson, D.B., 1948. Prediction of pressure drop during forced circulation boiling of water. Trans. ASME 70, 695
702. McAdams, W.H., Woods, W.K., Heroman Jr., R.H., 1942. Vaporisation inside horizontal tubes. II. Benzene-oil mixtures. Trans. ASME 64, 193. Mendler, O.J., Rathbw, A.S., Van Huff, N.E., Weiss, A., 1961. Natural circulation tests with water at 800 to 2000 Psia under non-boiling, local boiling and bulk boiling condition. J. Heat Transfer 261
273. Naveen, K., Rajalakshmi, R., Kulkarni, R.D., Sagar, T.V., Vijayan, P.K., Saha, D., 2000. Experimental Investigation in High Pressure Natural Circulation Loop. BARC/2000/E/002, Bhabha Atomic Research Centre, Trombay, Mumbai. Naveen, K., Doshi, J.B., Vijayan, P.K., Saha, D., 2010. Investigations on low power low pressure instability. In: 20th National and 9th international ISHMT-ASME Heat and Mass Transfer Conference, January 4
6, Mumbai. Naveen, K., Iyer, K.N., Doshi, J.B., Vijayan, P.K., 2014. Investigations on single-phase natural circulation loop dynamics Part 1: model for simulating start-up from rest. Prog. Nucl. Energy 76, 148
159. Naveen, K., Nayak, A.K., Rao, A.R., Vijayan, P.K., 2017a. NUTAN-Th(V1): A Multi-Channel Multi Physics Thermal-Hydraulic Analysis Code for Natural Circulation Systems. BARC Internal Report No. BARC/2017/I/006, Bhabha Atomic Research Centre, Mumbai. Naveen, K., Nayak, A.K., Rao, A.R., Vijayan, P.K., 2017b. AHWR Cold Start-Up Simulation Using Thermal Hydraulics Analysis Code. NUTAN-Th, BARC News Letter, ISSN 09762108, Jan-Feb 2017. Nayak, A.K., Vijayan, P.K., Saha, D., Venkat Raj, V., Aritomi, M., 1998. Linear analysis of thermo-hydraulic instabilities of the Advanced Heavy Water Reactor (AHWR). J. Nucl. Sci. Technol. 35, 768
778. Rao, G.S.S.P., Vijayan, P.K., Jain, V., Borgohain, A., Sharma, M., Nayak, A.K., Belokar, D.G., Pal, A.K., Saha, D., Sinha, R.K., 2002. AHWR integral test loop
Scaling Philosophy and System Description, BARC/2002/E/017, RP311983, Reactor Engineering Division, Bhabha Atomic Research Centre, Trombay, Mumbai. Sekoguchi, K., Saito, Y., Honda, T., 1970. JSME preprint No. 700-7, 83. The RELAP5 Development Team, 1995. RELAP5/MOD3 Code Manual. NUREG/CR-5535, INEL95/0174, Idaho National Engineering Laboratory, P.O. Box, 1625, Idaho Falls, ID 83415. Todreas, N.E., Kazimi, M.S., 1990. Nuclear Systems II, Elements of Thermal Hydraulic Design. Hemisphere Publications, New York. Van Bragt, D.D.B., Van Der Hagen, T.H.J.J., 1998. Stability of natural circulation Boiling Water Reactors: Part I
Description stability model and theoretical analysis in terms of dimensionless groups. Nucl. Technol. 121, 40
51. Wang, Q., Chen, X.J., Kakac, S., Ding, Y., 1994a. An experimental investigation of density-wave type oscillations in a convective boiling upflow system. Int. J. Multiphase Flow 15, 241
246. Wang, F.S., Hu, L.W., Pan, C., 1994b. Thermal and stability analysis of a two-phase natural circulation loop. Nucl. Sci. Eng. 117, 33
46.