- Email: [email protected]
Development of a thermal-hydraulic analysis software for the Chinese advanced pressurized water reactor
Nuclear Engineering and Design 240 (2010) 112–122
Contents lists available at ScienceDirect
Nuclear Engineering and Design journal homepage: www.elsevier.com/locate/nucengdes
Development of a thermal-hydraulic analysis software for the Chinese advanced pressurized water reactor Y.W. Wu, G.H. Su ∗ , S.Z. Qiu, C.J. Zhuang State Key Laboratory of Multiphase Flow in Power Engineering, Department of Nuclear Science and Technology, Xi’an Jiaotong University, Xi’an City 710049, China
a r t i c l e
i n f o
Article history: Received 16 March 2009 Received in revised form 5 October 2009 Accepted 13 October 2009
a b s t r a c t A point reactor neutron kinetics model, a drift-flow U-tube steam generator model, a non-equilibrium three-region pressurizer model and other models were established and a transient analysis code with Visual Fortran 6.5 has been developed to analyze the thermal-hydraulic characteristics of the Chinese advanced pressurized water reactor (AC-600). Visual input, real-time processing and dynamic visualization output were achieved with Microsoft Visual Studio.NET 2003, which greatly facilitate applications in the engineering. The software were applied to analyze the transient thermal-hydraulic characteristics of the loss of feed-water accident, the double loops loss-of-flow accident, the reactivity insertion accident, the sudden increase of feed-water temperature accident and the loss of offsite power accident for the Qinshan nuclear power plant in China. The obtained analysis results are significant to the improvement of design and safety operation of the plant. © 2009 Elsevier B.V. All rights reserved.
1. Introduction Along with the rapid development of computing and computer technologies, many safety analysis codes for nuclear power system accident or hypothetical accident have been developed. At present, the best large thermal-hydraulic analysis codes for the light water nuclear reactor systems are RELAP series, RETRAN series, TRAC series, etc. Among them many have become indispensable and important tools for the design, evaluation and safety analysis of nuclear power system. These large thermal-hydraulic analysis codes have considered many factors, as well as have complex and detailed mathematical models and comprehensive functions. Therefore these large codes require users to have high professional technology and program technique relatively. In addition, it is very difficult to make modifications and supplements to large codes because they have complex structure, high data storage and transmission techniques. The operation speed of large codes is relatively slow and cannot realize the real-time simulation. The above factors have gravely limited the applications of large codes. In this study, a thermal-hydraulic analysis software named TTHAsoft 1.0 for the Chinese advanced pressurized water reactor (AC-600) with two loops has been developed. Compared with RELAP5 and other commercial codes, the TTHAsoft 1.0 needs
∗ Corresponding author. Tel.: +86 29 82663401; fax: +86 29 82663401. E-mail address: [email protected] (G.H. Su). 0029-5493/$ – see front matter © 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.nucengdes.2009.10.020
less computational time and can realize the real-time calculation. Unlike other existing real-time simulation tools with dynamic visualization (e.g. plant simulators) which were generally developed on the basis of large amounts of experimental data or theoretical calculation data, the TTHAsoft 1.0 is based on fundamental conservation principles: the mass, momentum and energy conservation equations. In comparison with the complicated preparation work on the input card of RELAP5, the pre-process is much more convenient; and in comparison with the output results of RELAP5, the postprocess is much more advanced for its real-time processing and dynamic visualization output. At the mean time, the TTHAsoft 1.0 is highly modularized using Fortran 90 language; therefore it is much easier for the user to add new mathematical models and physical models for further development. For instance, more accurate heat transfer or flow friction correlations could be easily adopted by the TTHAsoft 1.0 by adding the corresponding subroutine program module. The safety analysis of the Qinshan nuclear power plant transient thermal-hydraulic behaviors have been conducted using this software. All analyses are aiming at verifying whether the plant is safety in operation. 2. Mathematical models The object of this study is AC-600 which is a single reactor with two loops nuclear power system. The basic field model of AC-600 is based on fundamental conservation principles: the mass, momentum and energy conservation equations. With the assumption of
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
Nomenclature A c Ci DNBR F H h k L M N N–S p q T T¯ t u
v V W x z
cross sectional area (m2 ) specific heat (kJ/(kg K)) neutron precursor concentration of group i departure from nucleate boiling ratio surface area (W/(m2 K)) heat transfer coefficient (W/(m2 K)) enthalpy (kJ/kg) feedback coefficient water level (m) mass (kg) kinetic power (MW) or number of the control volume in the primary loop Nassi and Shneiderman pressure (MPa) heat flux (W/m2 ) temperature (K) average temperature (K) time (s) velocity (m/s) specific volume (m3 /kg) volume (m3 ) mass flow rate (kg/s) quality special coordinate (m)
Greek symbols ˛ void fraction ˇ total delayed neutron fraction ˇi delayed neutron fraction of group i neutron generation time decay constant of precursor group i i reactivity or density Subscripts afw auxiliary feed water b1 ascending bubble c coolant ccd coolant and clad cd clad ce interface net vaporization cs spray condensate water dc downcomer dg drift-flow ex externally introduced f fuel, saturated liquid or frictional fcd fuel and clad fs fuel surface fw feed water g saturated vapor in inlet loc local m tube wall out outlet p primary pr pressurizer re release valve rir riser s secondary sa safety valve sd steam space sp spray water stm steam su surge water ws wall condensate water
113
one-dimensional flow, these equations, including single-phase and two-phase conservative equations, can be easily found in reference (Collier and Thome, 1994). The characteristics of the mathematical models are introduced in detail as follows. 2.1. Core model The core fission power is calculated through the solution of the point kinetic equations with six groups of delayed neutron. The power distribution in axial is assumed to be known.
(t) − ˇ dN(t) i Ci (t) = N(t) + dt 6
(1)
i=1
dCi (t) ˇ = i N(t) − i Ci (t), dt
i = 1, 2, . . . , 6
(2)
The reactivity feedbacks caused by the temperature variations of the moderator and the fuel are specially considered.
(t) = 0 + ex (t) + kc T¯ c (t) − T¯ c (0) + kf T¯ f (t) − T¯ f (0)
(3)
In Eq. (3), on the right hand side, the first term (0 ) is the initial reactivity; the second term (ex ) is the externally introduced reactivity; and the remaining terms are the various feedback contributions of the moderator temperature and the fuel temperature. The parameter T is the mass weighted average temperature of the fuel or the coolant; k is the feedback coefficient. The single channel model is chosen for the reactor core thermalhydraulic calculation. The fuel element is cylinder and the axial heat conduction is ignored. Thus, to the fuel element, the heat transfer equation can be defined as dTf (t) −Hfcd · Ffcd (Tfs − Tcd ) + qf = Mf cf dt
(4)
To the clad, the heat transfer equation can be defined as H · F (T − Tcd ) − Hccd · Fccd (Tcd − Tc ) dTcd (t) = fcd fcd fs Mcd ccd dt
(5)
To the coolant, the variation of the temperature can be defined as H · F (T − Tc ) − Wc (∂h/∂z) dTc (t) = ccd ccd cd Mc cc dt
(6)
2.2. Pressurizer model It is very important to accurately simulate the dynamic characteristics of the pressurizer in order to improve the simulation accuracy of the whole nuclear power system. Many scholars have studied on the mathematical simulation for dynamic characteristics of the pressurizer since the 1960s. In one early model, the steam and water in the pressurizer were thought to be at saturation state with same thermodynamic characteristics and the pressurizer was thermodynamic analyzed as a closed system. The model was very simple and convenient; however, the prediction results were not satisfactory when it was applied in the rapid variation processes. Later, Redfoeld (1968), Baron (1973) and Abdallah et al. (1982), etc. proposed a two-region unbalanced model, namely, the pressurizer was divided into the steam region and water region, and the thermodynamic state was unbalanced. The unbalanced model was more close to the actual processes of the pressurizer than the balance model, but it came to its shortages when the stratification phenomenon appeared after the super-cooled water flowing into the pressurizer while there was a positive system pressure surge. Baggovra and Martin (1983) proposed a three-region model, namely, two liquid phase regions and a vapor phase region which considered vaporization and condensation dynamics phenomena. However, the disadvantage of this model was that the
114
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
where Wsu is the surge coolant mass flow rate and is equal to the total of the surge coolant mass flow rate of each control volume in N the primary loop, Wsu = − 1 Vi (∂i /∂t). The calculation modes of Wb1 , Wws and Wce could get from previous research (Wilson, 1962; Kang, 1984). According to the energy conservation equation, d(M1 h1 ) dp = Wsu hx − Wb1 hg + qb1 + V1 dt dt
(10)
dp d(M2 h2 ) = Wsp hf + Wb1 hg + Wcs hf + Wws hf − Wce hg + qb2 + V2 dt dt (11)
d(M3 h3 ) dp = Wce hg − Wcs hg − Wws h3 − Wre h3 − Wsa h3 + V3 dt dt (12)
hsu positive fluctuation . h1 negative fluctuation The volume of the pressurizer is constant. Thus,
where hx =
d(M v ) dV i i =0 = dt dt 3
(13)
i=1
The pressure equation of the pressurizer can be obtained from Eq. (13):
Fig. 1. Schematic diagram of the pressurizer.
3
effect on the pressure variation, caused by the heat transfer occurring on the interface between the liquid phase and the vapor phase, was ignored. Compared with other heat transfer models, the heat transfer on the interface was dominant when the water level in the pressurizer was lower than the immersion depth of the positive surge super-cooled water. Therefore, Beak (1986) proposed a comprehensive three-region model, in which a lot of important thermal-hydraulic phenomena were considered, such as the bubble rising and disappearance, the wall condensation, and the heat and mass transfer on the interface. Here, considering the states of each region at different conditions and all important thermalhydraulic phenomena occurring in the pressurizer, a more precise pressurizer model was proposed based on the Beak’s three-region model. According to the different phases and enthalpies of the fluid in the pressurizer, the pressurizer is divided into three regions (Fig. 1), namely, the surge water region, the main water region and the steam region. And the following simplification assumptions are given: (1) The three regions share the same pressure at the same time; (2) The same fluid has the same thermodynamic parameters at the same area and time; (3) The heat dissipation to the atmosphere from the pressurizer is ignored. It is assumed that the surge water region, the main water region and the steam region are denoted by 1, 2 and 3 regions, respectively. According to the mass conservation equation, dM1 = Wsu − Wb1 dt
(7)
dM2 = Wb1 + Wsp + Wcs + Wws − Wce dt
(8)
dM3 = Wce − Wcs − Wws − Wre − Wsa dt
(9)
− dp = dt
i=1
(Mi (∂vi /∂hi )(dhi /dt) + vi (dMi /dt))
3
i=1
(14)
(Mi (∂vi /∂p))
The water level equation of the pressurizer can be obtained from the variations of the whole volumes of the main water region and the surge water region: dL = dt
2
i=1
[Mi ((∂vi /∂hi )(dhi /dt) + (∂vi /∂p)(dp/dt)) + vi (dMi /dt)] Apr (15)
2.3. Steam generator model The one-dimensional two-fluid model with drift-flow model was adopted for the vertical U-type tube natural circulation steam generator and the following simplification assumptions are given: (1) There is no heat conduction along the axial direction of the Utype tube; (2) There is no boiling in the downcomer, and no heat exchange between the riser and downcomer; (3) Heat capacity in the steam generator except the U-type heat transfer tube is ignored; (4) The water and steam in the steam separator and the steam space are saturated. The steam generator are divided into ten modules: (1) the primary fluid channel; (2) the feedwater plenum; (3) the downcomer; (4) the heat transfer tube wall; (5) the subcooled section of the secondary loop; (6) the boiling section of the secondary loop; (7) the riser and the steam separator; (8) the inlet plenum of the primary fluid; (9) the outlet plenum of the primary fluid; (10) the steam space. Thus, by applying the basic equations of mass, momentum and energy to the steam generator, the following equations are obtained:
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
115
(1) The pressure equation of the secondary loop [hg (Ws − Wstm − Wre ) − Vrir (g hg − f hf )(d˛out /dt)] dp = dt {(Vsd + ˛Vrir )(hg (∂g /∂p) + g (∂hg /∂p) − 1000) + [Vfw + Vrir (1 − ˛out )](hf (∂f /∂p) + f (∂hf /∂p))}
(16)
(2) The enthalpy equation of the primary fluid dhp 1 =− f Ap dt
∂hp qp − Wp ∂z
(17)
(3) The enthalpy equation of the secondary fluid (a) the enthalpy equation of the secondary fluid in the subcooled section dhs 1 =− f As dt
∂hs qs − Ws ∂z
+
1 dp f dt
(18)
(b) the enthalpy equation of the secondary fluid in the boiling section ∂ ∂ [˛g hg + (1 − ˛)f hf ] + [˛g hg ug + (1 − ˛)f hf uf ] ∂t ∂z =
dp qs + As dt
(19)
table lookup or linear interpolation. In this way, the four-quadrant curves can be easily applied to the code. 2.5. Auxiliary models The correct heat transfer and flow frictional factor correlations, critical heat flux model and coolant thermophysical property model are also key aspects of the quality of the results. According to the operation conditions of the core, the Collier (1981) correlation is chosen for the laminar heat transfer (Re < 2000), the Petukhov correlation (Thomas, 1992) is chosen for the turbulent flow heat transfer (Re > 3000), linear interpolation is applied between the laminar and turbulent correlations for the transient (2000 ≤ Re ≤ 3000) Nusselt number calculation. The natural con-
(4) The void fraction equation of the secondary fluid {(qp /A) + (∂p/∂t) + h(∂f /∂t) − (Ws /h)(∂h/∂t)(f hf ) − ˛[(∂/∂t)(g hg − f hf ) − h(∂gf /∂t)]} d˛ = dt (g hg − f hf − hgf )
(20)
where h = (1 − ˛)f uf hf + ˛g ug hg gf = f − g (5) The wall temperature equation of the U-type tube qp − qs dTm = Mm cm dt
(21)
vection heat transfer is considered generally after Gr/Re2 ≥ 0.1, and the McAdams (1954) correlation is adopted to take account of this. The Jens and Lottes (1951) correlation is adopted to calculate the subcooled boiling heat transfer, and the Chen (1966) correlation is chosen for the saturated boiling heat transfer. The laminar flow friction factor is calculated by f = 96/Re, the turbulent flow friction
(6) The enthalpy equation of the feedwater plenum dhfw {[(1 − xrir )Ws hf + Wfw hfw + Wafw hafw − Ws hdc ](Wfw + Wafw − xrir Ws )hdc } = (dc Adc Ldc ) dt (7) The water level equation (Wfw + Wafw − xrir Ws )hdc L dL = − dc (dc Adc ) dt
∂dc dhdc ∂dc dp + ∂p dt ∂hdc dt
(23) (8) The recirculation flow rate equation of the secondary loop dWs = dt
−
gdz − pf − ploc − pdg
dz/A
(22)
(24)
2.4. Reactor coolant pump model The characteristics of the reactor coolant pump have great influence on the safety operation of nuclear power plants, thus the transient analysis model for the safety analysis code is very important. The main parameters describing the reactor coolant pump are pump head, torque, volume flow rate and angular velocity. There are certain relations and internal laws among these parameters. Usually the curves describing the corresponding relationships among these main parameters are called four-quadrant curves of the pump, which identify uniquely the relationships among the pump head, torque, volume flow rate and angular velocity and could be obtained from experiments. However, it is very difficult to apply these curves into the code development, thus these curves are converted into some simple analogy curves which are the functions of the angular velocity ratio and the volume flow ratio (the ratio of real value to rated value). The analogy curves could conveniently be input in the form of tables. Therefore, as the functions of the independent variables, the dependent variables are gained by way of
factor is calculated by Colebrook (1939) correlation, the transient flow friction factor is calculated by interpolating the laminar and turbulent correlations. The two-phase frictional multiplier is calculated by Chisholm (1967) correlation. W-3 equation is adopted for the calculation of critical heat flux. The thermophysical property correlations of water and vapor used in this study are chosen from the international standard (Rohsenow and Hartnett, 1973). Comparison between the model prediction and look-up table is conducted to verify and extend the applicability of those models under low pressure. 3. Software programming 3.1. Nodalization In order to numerically simulate the thermal-hydraulic characteristics of the AC-600, the whole system is divided into many control volumes and junctions according to different geometrical and heat transfer conditions. Fig. 2 schematically shows the nodalization of the AC-600. The primary circuit models used in the thermal-hydraulic analysis consist of a core, a steam generator, a pump, a pressurizer, a downcomer and plenums. The models of the secondary are mainly composed of a steam generator, pipes, plenums and valves. In the core and the steam generator, fine control volume division is used to properly predict the heat transfer phenomena. The number of nodes of the core, the exchanger and pipes can be changed according to the calculation requirement.
116
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
Fig. 2. The systematic diagram of nodalization of AC-600.
3.2. Numerical method and code description Once the nodalization of a nuclear reactor system was obtained, the transient analysis can be reduced to the solve of a set of ordinary differential equations with the initial conditions as shown in Eq. (25), where y is a vector including mass flow, enthalpy and other parameters of each control volume. Because of the specialty and complexity of nuclear reactor system, the corresponding differential equations are usually stiff differential equations (also called ill-conditioned equations).
dy៝ = f៝ (t, y៝ , y៝ ) dt y៝ (0) = y៝ 0
(25)
Gear method (Gear, 1971; Hind marah, 1974) equipped with Adams predictor–corrector method was adopted for better solution of stiff differential equations. Gear method is highly efficient for the solution of ill-conditioned problems for its good stability, high precision, etc. It is the first order upwind difference method; in the mean time, the code of Gear method has the merit of self-stability and self-comparatively. Adams method required less computation time for its simpler iterative procedure, however, may fail for strong ill-conditioned equations. In this study, the two methods are alternatively used according to the stiffness of the set of equations in the whole solving process. For example, in the original stage after shutdown, the core power sharply decreased because of the large negative reactivity caused by shutdown rods, while the fuel temperature, the coolant enthalpy and others might gently vary because of the relative delay of heat transfer. Thus Gear
method is used in this case. During the residual heat removal process, all parameters hardly change and therefore Adams method is employed to save running time. It is found in this study that the alternatively adoption of Gear method and Adams methods greatly reduce the computation time which is important for conducting long-term simulation. Based on the theoretical model and solution method, a simulation code was developed. The code is in Fortran 90 format and can be maintained in the PC/Windows environment. For convenient maintenance and readability of the code, modular programming techniques were adopted. The main function modules are data input module, initialization module, transient module, derivative module, numerical method module, auxiliary module, thermophysical property module, output module, etc. All involved heat transfer and friction coefficient correlations are supplied in the auxiliary module. Derivative module is used to calculate the right hand of the differential equations. Modification of these modules or addition of new modules for different calculation requirements can be easily done. The calling relationship of these modules and the flowchart of the code are shown in Figs. 3 and 4, respectively. 3.3. Software design Microsoft Visual Studio.NET is chosen as the visual software development tool. The program method based on control components is adopted in the software which could be fast installed and operated in common PC computers. The software calls the simulation code and graphics controls by the visual interface to achieve
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
Fig. 3. Calling relationships of the modules.
117
the dynamic real-time simulation for the thermal-hydraulic characteristics of the reactor system. The software, named TTHAsoft 1.0, is mainly composed of three parts: the main window interface, the input parameters setting interface and the help file. As shown in Fig. 5, the main window interface is mainly applied to show the temperature variation of each control volume in the reactor system reflected by the color variation and the real-time dynamic curves for key parameters (DNBR, the largest fuel temperature Tmax ). In addition, the realtime dynamic curves for the key parameters of each system will be observed by selecting the buttons of label controls in the framework of the main window interface, namely core (as shown in Fig. 6), steam generator I, steam generator II, pressurizer and reactor coolant pump. The input parameters setting interface is designed for the convenience of users’ modifying the nuclear power system simulation parameters. As shown in Fig. 7, sufficient notes have been given for the input parameters avoiding that users to look up the manual for parameter meaning. The help file is a CHM file which contains detailed operating rules and demonstrations for the beginners. 3.4. Comparison of the software For nuclear reactor systems, it is impossible or unpractical to conduct full-scale experimental investigations. Therefore it is very difficult to compare the computational result with the real experimental result. However, the comparison between two different codes is feasible. In this paper, RELAP5/MOD3 was used to validate the TTHAsoft 1.0 which was developed by us. The RELAP5/MOD3 was developed on the basis of RELAP5/MOD2. The general physical models used for system thermal-hydraulic analysis are as follows:
Fig. 4. N–S flowchart of the code.
(1) Non-equilibrium, two fluid models for hydrodynamics including transportation of non-condensable gases. (2) 2D/3D capability provided through ‘cross-flow’ options. (3) Convective and radioactive heat transfer.
Fig. 5. Main window interface.
118
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
Fig. 6. Window interface for the core.
(4) 1D heat conduction in system structures. (5) Point reactor kinetics. (6) External 3D kinetics provided through link to user supplied reactor kinetics packages. (7) Control system, trip logic, and special system components such as valves and pumps.
The RELAP5/MOD3 has been widely used in accidental analyses of different types of reactors including the PWR, BWR, research reactors and even some thermo hydraulic experimental facilities. As for the convective heat transfer models of RELAP5/MOD3, different heat transfer coefficient correlations are applied under forced turbulent convection, laminar convection and free natural
Fig. 7. Setting interface for input parameters.
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
119
Table 1 Main parameters of AC-600. Parameter
Value
Core Thermal power (MW) Mass flow rate in the primary loop (kg/s) Core inlet temperature (K) Core outlet temperature (K) Active zone diameter (m) Active zone height (m) Total heat transfer area (m2 ) Total flow area of the active zone (m2 )
1035 3333.3 589.3 561.1 2.6 2.9 2248.8 2.5
Pressurizer Pressure (MPa) Water level (m) Total volume (m3 ) Diameter (m) Height (m)
15.5 5.4 34.9 2.1 10.8
Steam generator Mass flow rate in the secondary loop (kg/s) Circulation ratio Inlet pressure (MPa) Inlet temperature (K) Water level (m) Heat transfer area of the primary side (m2 ) Heat transfer area of the secondary side (m2 )
280.6 3.4 5.7 493.15 10.5 3090.1 3468.4
Main pump Rated head (m) Rated torque (N m) Rated volume flow rate (m3 /s) Rated angular velocity (rad/s)
75 20,310 4.4 156
convection condition. Using the maximum value ensures a smooth transition between correlations. The flow friction models used in both the RELAP5/MOD3 and the TTHAsoft 1.0 are the same. More details of physical models used in RELAP5/MOD3 can be found in the code manual reference book. The major structural and initial parameters from the Qinshan AC-600 in China used for the theoretical analysis are shown in Table 1. Fig. 8(A and B) shows the change curves of the reactor thermal power, the primary loop mass flow rate, the pressure in the pressurizer, and the water level in the steam generator with time respectively in the loss of feed-water accident. In the initial stage (t < 12 s), when the feed-water is lost, the secondary side of the steam generator will generate more steam and then the water level of the steam generator decreases. If the primary coolant is not cooled down sufficiently, the primary pressure will increase because of the expansion of the coolant. When the water level in the steam generator decreases to the low-low water level at about 12 s, the reactor is shut down. Since that the steam generator releases steam continuously through the release valve after the emergency shutdown of the reactor, the water level will continue to decrease until the auxiliary feed-water system start to work (at about 60 s). The results show that the cooling of the primary loop could return to normal state in feed-water loss accident as long as the auxiliary feed-water system operating normally. The calculated results obtained through this present software are in good agreements with that by RELAP5. In addition, some typical accidents, such as the double loops loss-of-flow accident, the reactivity insertion accident, the sudden increase of feed-water temperature accident, the loss of offsite power accident are also verified and also show to be in good agreements with results obtained by RELAP5. This proves the correctness of the mathematical models, the validity of the numerical method and the reliability of the software. Also, this software is expected to be validated against experimental data.
Fig. 8. Comparison between the results calculated by TTHAsoft 1.0 and RELAP5 for the loss of feed-water accident. (A) Thermal power, (B) pressure and water level.
4. Calculation results and discussions 4.1. Double loops loss-of-flow accident Double loops loss-of-flow accident occurs when the two reactor coolant pumps idle at the same time due to mechanical or electrical failures. The protective measures will shut down the reactor immediately when there is a low-flow signal in any loops, and the main feed-water supply stops and the turbine trips at the same time. In the accident analysis, the bypass is not considered after the steam turbine trips, and the emergency shutdown will be only triggered by the low-flow signal. Fig. 9 presents the characteristic change curves of the typical parameters with time in the accident. During the initial stage of the accident, the mass flow rate in the primary loop decreases rapidly for the two reactor coolant pumps idle at the same time. However, the reactor power changes slightly and thus the ratio of the power to the mass flow rate increases. The heat leg exit quality and maximum pellet temperature increase while the DNBR decreases (to the minimum value, about 1.78 in about 5 s). The low-flow signal of the emergency shutdown is triggered at about 3.5 s and the power decreases rapidly which restrains the decrease of the DNBR. Since the decreasing speed of the power is larger than that of the mass flow rate, the hot leg exit quality and the maximum pellet temperature begin to decrease after reaching to the peak, the DNBR starts to increase, and then the accident has been mitigated. Therefore, the minimum DNBR not exceeding the safety limitation is the main safety criterion in the initial stage of the accident while the residual
120
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
Fig. 9. Curves for the double loops loss-of-flow accident.
heat in the core must be removed correctly in order to ensure the integrity of the core in long term. 4.2. Reactivity insertion accident Reactivity insertion accident may occur during the refueling, the unexpected withdrawing of the control rod or the improper regulation of boron. These events are considered the most severe accidents that could lead to core damage, because during such events the core become supercritical and the core power rises to level beyond the heat removal system’s capability. In this paper, the variations of the transient parameters are analyzed while a +0.11$ reactivity is inserted into the system as a result of the improper regulation of control rod or boron. Here, it is assumed that the accident does not cause the core shutdown or other reactor protection systems actions and the power is completely controlled by negative feedback effects caused by the Doppler action and moderator-density. As shown in Fig. 10, the power increases rapidly because of the sudden increase of the reactivity, then decreases gradually as a result of the negative feedback effect, and goes back to the original power level finally. For the external load does not change, the pressure in the secondary side of the steam generator increases and
the water level decreases a little. In the initial stage of the accident, the increase of the power makes the DNBR decrease. Under the present calculation conditions, the minimum DNBR obtained form the W-3 equation as mentioned in Section 2.5 drops to 1.71, which is still within the safety limitation. Therefore, the nuclear power system could adjust itself to a stable and safety state while a +0.11$ reactivity is inserted into the system. 4.3. Sudden increase of feed-water temperature accident All the protection, control and regulation systems are not put into operation except for the reactivity feedbacks in the sudden increase of the feed-water temperature accident. Fig. 11 shows the change curves of the reactor thermal power, the reactivity and the pressure in the pressurizer with time assuming a sudden increase of the feed-water temperature by 30 K. It is shown that the increase of the inlet feed-water enthalpy (temperature) of the secondary side of the steam generator will make the secondary side quality increase which would cause the increases of the pressure and the water level in the secondary side of the steam generator. The heat removal capability of the steam generator decreases. Thus, the system pressure and temperatures in the primary loop will increase. The negative feedback effect of the reactivity caused
Fig. 10. Curves for the reactivity insertion accident.
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
121
Fig. 11. Curves for the sudden increase of feed-water temperature accident.
by the increases of the system pressure and temperatures could make the power decrease to a stable power level finally. Generally speaking, the parameters of the primary and secondary loops change slightly under the disturbance of the feed-water temperature. 4.4. Loss of offsite power accident The electricity supply of the reactor coolant pumps is interrupted in the loss of offsite power accident and then the coolant mass flow rate decreases rapidly because the driving force of the coolant is replaced by the inertial rotation of the pump and the inertial flow of the coolant. At the same time, a signal triggered by the power loss would make the reactor shut down, the main feedwater stop and the main steam line close. The steam bypass is not considered after the main steam valve has been shut down and it is assumed that there is without delay among the shutdowns of the reactor, the pump, the main feed-water, the main steam line and the loss of offsite power. If the delay time of the reactor trip protection is too long, the decrease velocity of the coolant mass flow rate is larger than that of the heat flux on the fuel element surface. And it may result in the decrease of the safety margin and the increase of the fuel element surface temperatures which may lead to the clad meltdown. The water level in the steam generator decreases as a result of the steam evacuation in the secondary loop. If the feed-water is not supplied in time and the water level reaches to the tube sheet height, the heat transfer efficiency would reduce further which is worse for the primary loop heat removal. In this paper, the emergency feed-water is assumed to be supplied after 60 s to relieve the decrease of the water level. The thermal power of the reactor and the mass flow rate decrease rapidly after the loss of offsite power as shown in Fig. 12(A). However, in the initial stage, the decrease velocity of the coolant mass flow rate is less than that of the thermal power. Thus, the minimum DNBR increases as shown in Fig. 12(B). The minimum DNBR will decrease slightly when the coolant mass flow rate and the thermal power tend to be stable. As shown in Fig. 12(B), the water level in the steam generator decreases as a result of the steam evacuation in the secondary loop. However, the water level will increase and tend to be stable after the emergency feed-water is supplied. The consequences of the loss of offsite power accident are less severe than the double loops loss-of-flow accident.
Fig. 12. Curves for the loss of offsite power accident. (A) Thermal power and mass flow rate, (B) MDNBR and water level.
5. Conclusion An analysis software named TTHAsoft 1.0 has been developed using the Fortran 90 and Visual C++ languages to evaluate the transient thermal-hydraulic behaviors of the AC-600. The simulated results by TTHAsoft 1.0 were compared with those of the
122
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
RELAP5/MOD3 code. Although there were differences between some local parameters, which were caused by the adoption of different models in the two codes, a good agreement between the computational results by the two codes was obtained which proved the applicability and accuracy of this software. Also, this software is expected to be validated against experimental data. Because of the adoption of modular programming techniques, this analysis software is expected to be applied to other reactors by easily modifying the corresponding function modules. The software has been applied to analyze the transient thermal-hydraulic characteristics of the Qinshan nuclear power plant in China. The obtained analysis results can be used as reference of the improvement design of the AC-600 and the safety operation of the Qinshan nuclear power plant. References Abdallah, M.A., Ahmed, H.M., Mohammad, A.R., Mohammad, E.S.N., 1982. Pressurizer transients dynamic model. Nucl. Eng. Des. 73, 447–453. Baggovra, B., Martin, W.R., 1983. Transient analysis of the three-Mile island unit 2 pressurizer system. Nucl. Technol. 62, 407–416.
Baron, B.C., 1973. Digital simulation of a nuclear pressurizer. Nucl. Sci. Eng. 52, 283–291. Beak, S.M., 1986. A Non-equilibrium three region model for transient analysis of pressurized water reactor pressurizer system. Nucl. Technol. 74, 213–221. Chen, J.C., 1966. Process design development. Ind. Eng. Chem. 5, 322. Chisholm, D., 1967. Naval Electronics Laboratory Report, No. 310. Colebrook, C.F., 1939. Turbulent flow in pipes with particular reference to the transition region between the smooth and rough pipe laws. Proc. Inst. Civ. Eng. 11, 133. Collier, J.G., 1981. Convective Boiling and Condensation, second ed. McGraw-Hill. Collier, J.G., Thome, J.R., 1994. Convective Boiling and Condensation, 3rd ed. Oxford University Press. Gear, C.W., 1971. Numerical Initial Value Problems in Ordinary Differential Equation. Prentice-Hall, Englewood Cliffs, NJ. Hind marah, A.C., 1974. Gear: ordinary differential equation system solver, Lawrence Livarmore Laboratory, Report UCID-30001, Revision 3, December. Jens, W.H., Lottes, P.A., 1951. Argonne National Laboratory, Report No. 4627. Kang, S.W., 1984. Pool heat transfer in a simulated PWR pressurizer. Trans. Am. Nucl. Soc. 46, 845–854. McAdams, W.H., 1954. Heat Transmission, third ed. McGraw-Hill, New York. Redfoeld, J.A., 1968. Pressurizer performance during loss of load tests at shippingport: analysis and test. Nucl. Appl. 4, 189–196. Rohsenow, W.M., Hartnett, J.R., 1973. Hand Book of Heat Transfer. McGraw-Hill, New York. Thomas, L.C., 1992. Heat Transfer. Prentice Hall, New York. Wilson, J.F., 1962. The velocity of rising steam in a bubble two phase mixture. Am. Nucl. Soc. 5, 151–163.
Contents lists available at ScienceDirect
Nuclear Engineering and Design journal homepage: www.elsevier.com/locate/nucengdes
Development of a thermal-hydraulic analysis software for the Chinese advanced pressurized water reactor Y.W. Wu, G.H. Su ∗ , S.Z. Qiu, C.J. Zhuang State Key Laboratory of Multiphase Flow in Power Engineering, Department of Nuclear Science and Technology, Xi’an Jiaotong University, Xi’an City 710049, China
a r t i c l e
i n f o
Article history: Received 16 March 2009 Received in revised form 5 October 2009 Accepted 13 October 2009
a b s t r a c t A point reactor neutron kinetics model, a drift-flow U-tube steam generator model, a non-equilibrium three-region pressurizer model and other models were established and a transient analysis code with Visual Fortran 6.5 has been developed to analyze the thermal-hydraulic characteristics of the Chinese advanced pressurized water reactor (AC-600). Visual input, real-time processing and dynamic visualization output were achieved with Microsoft Visual Studio.NET 2003, which greatly facilitate applications in the engineering. The software were applied to analyze the transient thermal-hydraulic characteristics of the loss of feed-water accident, the double loops loss-of-flow accident, the reactivity insertion accident, the sudden increase of feed-water temperature accident and the loss of offsite power accident for the Qinshan nuclear power plant in China. The obtained analysis results are significant to the improvement of design and safety operation of the plant. © 2009 Elsevier B.V. All rights reserved.
1. Introduction Along with the rapid development of computing and computer technologies, many safety analysis codes for nuclear power system accident or hypothetical accident have been developed. At present, the best large thermal-hydraulic analysis codes for the light water nuclear reactor systems are RELAP series, RETRAN series, TRAC series, etc. Among them many have become indispensable and important tools for the design, evaluation and safety analysis of nuclear power system. These large thermal-hydraulic analysis codes have considered many factors, as well as have complex and detailed mathematical models and comprehensive functions. Therefore these large codes require users to have high professional technology and program technique relatively. In addition, it is very difficult to make modifications and supplements to large codes because they have complex structure, high data storage and transmission techniques. The operation speed of large codes is relatively slow and cannot realize the real-time simulation. The above factors have gravely limited the applications of large codes. In this study, a thermal-hydraulic analysis software named TTHAsoft 1.0 for the Chinese advanced pressurized water reactor (AC-600) with two loops has been developed. Compared with RELAP5 and other commercial codes, the TTHAsoft 1.0 needs
∗ Corresponding author. Tel.: +86 29 82663401; fax: +86 29 82663401. E-mail address: [email protected] (G.H. Su). 0029-5493/$ – see front matter © 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.nucengdes.2009.10.020
less computational time and can realize the real-time calculation. Unlike other existing real-time simulation tools with dynamic visualization (e.g. plant simulators) which were generally developed on the basis of large amounts of experimental data or theoretical calculation data, the TTHAsoft 1.0 is based on fundamental conservation principles: the mass, momentum and energy conservation equations. In comparison with the complicated preparation work on the input card of RELAP5, the pre-process is much more convenient; and in comparison with the output results of RELAP5, the postprocess is much more advanced for its real-time processing and dynamic visualization output. At the mean time, the TTHAsoft 1.0 is highly modularized using Fortran 90 language; therefore it is much easier for the user to add new mathematical models and physical models for further development. For instance, more accurate heat transfer or flow friction correlations could be easily adopted by the TTHAsoft 1.0 by adding the corresponding subroutine program module. The safety analysis of the Qinshan nuclear power plant transient thermal-hydraulic behaviors have been conducted using this software. All analyses are aiming at verifying whether the plant is safety in operation. 2. Mathematical models The object of this study is AC-600 which is a single reactor with two loops nuclear power system. The basic field model of AC-600 is based on fundamental conservation principles: the mass, momentum and energy conservation equations. With the assumption of
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
Nomenclature A c Ci DNBR F H h k L M N N–S p q T T¯ t u
v V W x z
cross sectional area (m2 ) specific heat (kJ/(kg K)) neutron precursor concentration of group i departure from nucleate boiling ratio surface area (W/(m2 K)) heat transfer coefficient (W/(m2 K)) enthalpy (kJ/kg) feedback coefficient water level (m) mass (kg) kinetic power (MW) or number of the control volume in the primary loop Nassi and Shneiderman pressure (MPa) heat flux (W/m2 ) temperature (K) average temperature (K) time (s) velocity (m/s) specific volume (m3 /kg) volume (m3 ) mass flow rate (kg/s) quality special coordinate (m)
Greek symbols ˛ void fraction ˇ total delayed neutron fraction ˇi delayed neutron fraction of group i neutron generation time decay constant of precursor group i i reactivity or density Subscripts afw auxiliary feed water b1 ascending bubble c coolant ccd coolant and clad cd clad ce interface net vaporization cs spray condensate water dc downcomer dg drift-flow ex externally introduced f fuel, saturated liquid or frictional fcd fuel and clad fs fuel surface fw feed water g saturated vapor in inlet loc local m tube wall out outlet p primary pr pressurizer re release valve rir riser s secondary sa safety valve sd steam space sp spray water stm steam su surge water ws wall condensate water
113
one-dimensional flow, these equations, including single-phase and two-phase conservative equations, can be easily found in reference (Collier and Thome, 1994). The characteristics of the mathematical models are introduced in detail as follows. 2.1. Core model The core fission power is calculated through the solution of the point kinetic equations with six groups of delayed neutron. The power distribution in axial is assumed to be known.
(t) − ˇ dN(t) i Ci (t) = N(t) + dt 6
(1)
i=1
dCi (t) ˇ = i N(t) − i Ci (t), dt
i = 1, 2, . . . , 6
(2)
The reactivity feedbacks caused by the temperature variations of the moderator and the fuel are specially considered.
(t) = 0 + ex (t) + kc T¯ c (t) − T¯ c (0) + kf T¯ f (t) − T¯ f (0)
(3)
In Eq. (3), on the right hand side, the first term (0 ) is the initial reactivity; the second term (ex ) is the externally introduced reactivity; and the remaining terms are the various feedback contributions of the moderator temperature and the fuel temperature. The parameter T is the mass weighted average temperature of the fuel or the coolant; k is the feedback coefficient. The single channel model is chosen for the reactor core thermalhydraulic calculation. The fuel element is cylinder and the axial heat conduction is ignored. Thus, to the fuel element, the heat transfer equation can be defined as dTf (t) −Hfcd · Ffcd (Tfs − Tcd ) + qf = Mf cf dt
(4)
To the clad, the heat transfer equation can be defined as H · F (T − Tcd ) − Hccd · Fccd (Tcd − Tc ) dTcd (t) = fcd fcd fs Mcd ccd dt
(5)
To the coolant, the variation of the temperature can be defined as H · F (T − Tc ) − Wc (∂h/∂z) dTc (t) = ccd ccd cd Mc cc dt
(6)
2.2. Pressurizer model It is very important to accurately simulate the dynamic characteristics of the pressurizer in order to improve the simulation accuracy of the whole nuclear power system. Many scholars have studied on the mathematical simulation for dynamic characteristics of the pressurizer since the 1960s. In one early model, the steam and water in the pressurizer were thought to be at saturation state with same thermodynamic characteristics and the pressurizer was thermodynamic analyzed as a closed system. The model was very simple and convenient; however, the prediction results were not satisfactory when it was applied in the rapid variation processes. Later, Redfoeld (1968), Baron (1973) and Abdallah et al. (1982), etc. proposed a two-region unbalanced model, namely, the pressurizer was divided into the steam region and water region, and the thermodynamic state was unbalanced. The unbalanced model was more close to the actual processes of the pressurizer than the balance model, but it came to its shortages when the stratification phenomenon appeared after the super-cooled water flowing into the pressurizer while there was a positive system pressure surge. Baggovra and Martin (1983) proposed a three-region model, namely, two liquid phase regions and a vapor phase region which considered vaporization and condensation dynamics phenomena. However, the disadvantage of this model was that the
114
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
where Wsu is the surge coolant mass flow rate and is equal to the total of the surge coolant mass flow rate of each control volume in N the primary loop, Wsu = − 1 Vi (∂i /∂t). The calculation modes of Wb1 , Wws and Wce could get from previous research (Wilson, 1962; Kang, 1984). According to the energy conservation equation, d(M1 h1 ) dp = Wsu hx − Wb1 hg + qb1 + V1 dt dt
(10)
dp d(M2 h2 ) = Wsp hf + Wb1 hg + Wcs hf + Wws hf − Wce hg + qb2 + V2 dt dt (11)
d(M3 h3 ) dp = Wce hg − Wcs hg − Wws h3 − Wre h3 − Wsa h3 + V3 dt dt (12)
hsu positive fluctuation . h1 negative fluctuation The volume of the pressurizer is constant. Thus,
where hx =
d(M v ) dV i i =0 = dt dt 3
(13)
i=1
The pressure equation of the pressurizer can be obtained from Eq. (13):
Fig. 1. Schematic diagram of the pressurizer.
3
effect on the pressure variation, caused by the heat transfer occurring on the interface between the liquid phase and the vapor phase, was ignored. Compared with other heat transfer models, the heat transfer on the interface was dominant when the water level in the pressurizer was lower than the immersion depth of the positive surge super-cooled water. Therefore, Beak (1986) proposed a comprehensive three-region model, in which a lot of important thermal-hydraulic phenomena were considered, such as the bubble rising and disappearance, the wall condensation, and the heat and mass transfer on the interface. Here, considering the states of each region at different conditions and all important thermalhydraulic phenomena occurring in the pressurizer, a more precise pressurizer model was proposed based on the Beak’s three-region model. According to the different phases and enthalpies of the fluid in the pressurizer, the pressurizer is divided into three regions (Fig. 1), namely, the surge water region, the main water region and the steam region. And the following simplification assumptions are given: (1) The three regions share the same pressure at the same time; (2) The same fluid has the same thermodynamic parameters at the same area and time; (3) The heat dissipation to the atmosphere from the pressurizer is ignored. It is assumed that the surge water region, the main water region and the steam region are denoted by 1, 2 and 3 regions, respectively. According to the mass conservation equation, dM1 = Wsu − Wb1 dt
(7)
dM2 = Wb1 + Wsp + Wcs + Wws − Wce dt
(8)
dM3 = Wce − Wcs − Wws − Wre − Wsa dt
(9)
− dp = dt
i=1
(Mi (∂vi /∂hi )(dhi /dt) + vi (dMi /dt))
3
i=1
(14)
(Mi (∂vi /∂p))
The water level equation of the pressurizer can be obtained from the variations of the whole volumes of the main water region and the surge water region: dL = dt
2
i=1
[Mi ((∂vi /∂hi )(dhi /dt) + (∂vi /∂p)(dp/dt)) + vi (dMi /dt)] Apr (15)
2.3. Steam generator model The one-dimensional two-fluid model with drift-flow model was adopted for the vertical U-type tube natural circulation steam generator and the following simplification assumptions are given: (1) There is no heat conduction along the axial direction of the Utype tube; (2) There is no boiling in the downcomer, and no heat exchange between the riser and downcomer; (3) Heat capacity in the steam generator except the U-type heat transfer tube is ignored; (4) The water and steam in the steam separator and the steam space are saturated. The steam generator are divided into ten modules: (1) the primary fluid channel; (2) the feedwater plenum; (3) the downcomer; (4) the heat transfer tube wall; (5) the subcooled section of the secondary loop; (6) the boiling section of the secondary loop; (7) the riser and the steam separator; (8) the inlet plenum of the primary fluid; (9) the outlet plenum of the primary fluid; (10) the steam space. Thus, by applying the basic equations of mass, momentum and energy to the steam generator, the following equations are obtained:
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
115
(1) The pressure equation of the secondary loop [hg (Ws − Wstm − Wre ) − Vrir (g hg − f hf )(d˛out /dt)] dp = dt {(Vsd + ˛Vrir )(hg (∂g /∂p) + g (∂hg /∂p) − 1000) + [Vfw + Vrir (1 − ˛out )](hf (∂f /∂p) + f (∂hf /∂p))}
(16)
(2) The enthalpy equation of the primary fluid dhp 1 =− f Ap dt
∂hp qp − Wp ∂z
(17)
(3) The enthalpy equation of the secondary fluid (a) the enthalpy equation of the secondary fluid in the subcooled section dhs 1 =− f As dt
∂hs qs − Ws ∂z
+
1 dp f dt
(18)
(b) the enthalpy equation of the secondary fluid in the boiling section ∂ ∂ [˛g hg + (1 − ˛)f hf ] + [˛g hg ug + (1 − ˛)f hf uf ] ∂t ∂z =
dp qs + As dt
(19)
table lookup or linear interpolation. In this way, the four-quadrant curves can be easily applied to the code. 2.5. Auxiliary models The correct heat transfer and flow frictional factor correlations, critical heat flux model and coolant thermophysical property model are also key aspects of the quality of the results. According to the operation conditions of the core, the Collier (1981) correlation is chosen for the laminar heat transfer (Re < 2000), the Petukhov correlation (Thomas, 1992) is chosen for the turbulent flow heat transfer (Re > 3000), linear interpolation is applied between the laminar and turbulent correlations for the transient (2000 ≤ Re ≤ 3000) Nusselt number calculation. The natural con-
(4) The void fraction equation of the secondary fluid {(qp /A) + (∂p/∂t) + h(∂f /∂t) − (Ws /h)(∂h/∂t)(f hf ) − ˛[(∂/∂t)(g hg − f hf ) − h(∂gf /∂t)]} d˛ = dt (g hg − f hf − hgf )
(20)
where h = (1 − ˛)f uf hf + ˛g ug hg gf = f − g (5) The wall temperature equation of the U-type tube qp − qs dTm = Mm cm dt
(21)
vection heat transfer is considered generally after Gr/Re2 ≥ 0.1, and the McAdams (1954) correlation is adopted to take account of this. The Jens and Lottes (1951) correlation is adopted to calculate the subcooled boiling heat transfer, and the Chen (1966) correlation is chosen for the saturated boiling heat transfer. The laminar flow friction factor is calculated by f = 96/Re, the turbulent flow friction
(6) The enthalpy equation of the feedwater plenum dhfw {[(1 − xrir )Ws hf + Wfw hfw + Wafw hafw − Ws hdc ](Wfw + Wafw − xrir Ws )hdc } = (dc Adc Ldc ) dt (7) The water level equation (Wfw + Wafw − xrir Ws )hdc L dL = − dc (dc Adc ) dt
∂dc dhdc ∂dc dp + ∂p dt ∂hdc dt
(23) (8) The recirculation flow rate equation of the secondary loop dWs = dt
−
gdz − pf − ploc − pdg
dz/A
(22)
(24)
2.4. Reactor coolant pump model The characteristics of the reactor coolant pump have great influence on the safety operation of nuclear power plants, thus the transient analysis model for the safety analysis code is very important. The main parameters describing the reactor coolant pump are pump head, torque, volume flow rate and angular velocity. There are certain relations and internal laws among these parameters. Usually the curves describing the corresponding relationships among these main parameters are called four-quadrant curves of the pump, which identify uniquely the relationships among the pump head, torque, volume flow rate and angular velocity and could be obtained from experiments. However, it is very difficult to apply these curves into the code development, thus these curves are converted into some simple analogy curves which are the functions of the angular velocity ratio and the volume flow ratio (the ratio of real value to rated value). The analogy curves could conveniently be input in the form of tables. Therefore, as the functions of the independent variables, the dependent variables are gained by way of
factor is calculated by Colebrook (1939) correlation, the transient flow friction factor is calculated by interpolating the laminar and turbulent correlations. The two-phase frictional multiplier is calculated by Chisholm (1967) correlation. W-3 equation is adopted for the calculation of critical heat flux. The thermophysical property correlations of water and vapor used in this study are chosen from the international standard (Rohsenow and Hartnett, 1973). Comparison between the model prediction and look-up table is conducted to verify and extend the applicability of those models under low pressure. 3. Software programming 3.1. Nodalization In order to numerically simulate the thermal-hydraulic characteristics of the AC-600, the whole system is divided into many control volumes and junctions according to different geometrical and heat transfer conditions. Fig. 2 schematically shows the nodalization of the AC-600. The primary circuit models used in the thermal-hydraulic analysis consist of a core, a steam generator, a pump, a pressurizer, a downcomer and plenums. The models of the secondary are mainly composed of a steam generator, pipes, plenums and valves. In the core and the steam generator, fine control volume division is used to properly predict the heat transfer phenomena. The number of nodes of the core, the exchanger and pipes can be changed according to the calculation requirement.
116
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
Fig. 2. The systematic diagram of nodalization of AC-600.
3.2. Numerical method and code description Once the nodalization of a nuclear reactor system was obtained, the transient analysis can be reduced to the solve of a set of ordinary differential equations with the initial conditions as shown in Eq. (25), where y is a vector including mass flow, enthalpy and other parameters of each control volume. Because of the specialty and complexity of nuclear reactor system, the corresponding differential equations are usually stiff differential equations (also called ill-conditioned equations).
dy៝ = f៝ (t, y៝ , y៝ ) dt y៝ (0) = y៝ 0
(25)
Gear method (Gear, 1971; Hind marah, 1974) equipped with Adams predictor–corrector method was adopted for better solution of stiff differential equations. Gear method is highly efficient for the solution of ill-conditioned problems for its good stability, high precision, etc. It is the first order upwind difference method; in the mean time, the code of Gear method has the merit of self-stability and self-comparatively. Adams method required less computation time for its simpler iterative procedure, however, may fail for strong ill-conditioned equations. In this study, the two methods are alternatively used according to the stiffness of the set of equations in the whole solving process. For example, in the original stage after shutdown, the core power sharply decreased because of the large negative reactivity caused by shutdown rods, while the fuel temperature, the coolant enthalpy and others might gently vary because of the relative delay of heat transfer. Thus Gear
method is used in this case. During the residual heat removal process, all parameters hardly change and therefore Adams method is employed to save running time. It is found in this study that the alternatively adoption of Gear method and Adams methods greatly reduce the computation time which is important for conducting long-term simulation. Based on the theoretical model and solution method, a simulation code was developed. The code is in Fortran 90 format and can be maintained in the PC/Windows environment. For convenient maintenance and readability of the code, modular programming techniques were adopted. The main function modules are data input module, initialization module, transient module, derivative module, numerical method module, auxiliary module, thermophysical property module, output module, etc. All involved heat transfer and friction coefficient correlations are supplied in the auxiliary module. Derivative module is used to calculate the right hand of the differential equations. Modification of these modules or addition of new modules for different calculation requirements can be easily done. The calling relationship of these modules and the flowchart of the code are shown in Figs. 3 and 4, respectively. 3.3. Software design Microsoft Visual Studio.NET is chosen as the visual software development tool. The program method based on control components is adopted in the software which could be fast installed and operated in common PC computers. The software calls the simulation code and graphics controls by the visual interface to achieve
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
Fig. 3. Calling relationships of the modules.
117
the dynamic real-time simulation for the thermal-hydraulic characteristics of the reactor system. The software, named TTHAsoft 1.0, is mainly composed of three parts: the main window interface, the input parameters setting interface and the help file. As shown in Fig. 5, the main window interface is mainly applied to show the temperature variation of each control volume in the reactor system reflected by the color variation and the real-time dynamic curves for key parameters (DNBR, the largest fuel temperature Tmax ). In addition, the realtime dynamic curves for the key parameters of each system will be observed by selecting the buttons of label controls in the framework of the main window interface, namely core (as shown in Fig. 6), steam generator I, steam generator II, pressurizer and reactor coolant pump. The input parameters setting interface is designed for the convenience of users’ modifying the nuclear power system simulation parameters. As shown in Fig. 7, sufficient notes have been given for the input parameters avoiding that users to look up the manual for parameter meaning. The help file is a CHM file which contains detailed operating rules and demonstrations for the beginners. 3.4. Comparison of the software For nuclear reactor systems, it is impossible or unpractical to conduct full-scale experimental investigations. Therefore it is very difficult to compare the computational result with the real experimental result. However, the comparison between two different codes is feasible. In this paper, RELAP5/MOD3 was used to validate the TTHAsoft 1.0 which was developed by us. The RELAP5/MOD3 was developed on the basis of RELAP5/MOD2. The general physical models used for system thermal-hydraulic analysis are as follows:
Fig. 4. N–S flowchart of the code.
(1) Non-equilibrium, two fluid models for hydrodynamics including transportation of non-condensable gases. (2) 2D/3D capability provided through ‘cross-flow’ options. (3) Convective and radioactive heat transfer.
Fig. 5. Main window interface.
118
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
Fig. 6. Window interface for the core.
(4) 1D heat conduction in system structures. (5) Point reactor kinetics. (6) External 3D kinetics provided through link to user supplied reactor kinetics packages. (7) Control system, trip logic, and special system components such as valves and pumps.
The RELAP5/MOD3 has been widely used in accidental analyses of different types of reactors including the PWR, BWR, research reactors and even some thermo hydraulic experimental facilities. As for the convective heat transfer models of RELAP5/MOD3, different heat transfer coefficient correlations are applied under forced turbulent convection, laminar convection and free natural
Fig. 7. Setting interface for input parameters.
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
119
Table 1 Main parameters of AC-600. Parameter
Value
Core Thermal power (MW) Mass flow rate in the primary loop (kg/s) Core inlet temperature (K) Core outlet temperature (K) Active zone diameter (m) Active zone height (m) Total heat transfer area (m2 ) Total flow area of the active zone (m2 )
1035 3333.3 589.3 561.1 2.6 2.9 2248.8 2.5
Pressurizer Pressure (MPa) Water level (m) Total volume (m3 ) Diameter (m) Height (m)
15.5 5.4 34.9 2.1 10.8
Steam generator Mass flow rate in the secondary loop (kg/s) Circulation ratio Inlet pressure (MPa) Inlet temperature (K) Water level (m) Heat transfer area of the primary side (m2 ) Heat transfer area of the secondary side (m2 )
280.6 3.4 5.7 493.15 10.5 3090.1 3468.4
Main pump Rated head (m) Rated torque (N m) Rated volume flow rate (m3 /s) Rated angular velocity (rad/s)
75 20,310 4.4 156
convection condition. Using the maximum value ensures a smooth transition between correlations. The flow friction models used in both the RELAP5/MOD3 and the TTHAsoft 1.0 are the same. More details of physical models used in RELAP5/MOD3 can be found in the code manual reference book. The major structural and initial parameters from the Qinshan AC-600 in China used for the theoretical analysis are shown in Table 1. Fig. 8(A and B) shows the change curves of the reactor thermal power, the primary loop mass flow rate, the pressure in the pressurizer, and the water level in the steam generator with time respectively in the loss of feed-water accident. In the initial stage (t < 12 s), when the feed-water is lost, the secondary side of the steam generator will generate more steam and then the water level of the steam generator decreases. If the primary coolant is not cooled down sufficiently, the primary pressure will increase because of the expansion of the coolant. When the water level in the steam generator decreases to the low-low water level at about 12 s, the reactor is shut down. Since that the steam generator releases steam continuously through the release valve after the emergency shutdown of the reactor, the water level will continue to decrease until the auxiliary feed-water system start to work (at about 60 s). The results show that the cooling of the primary loop could return to normal state in feed-water loss accident as long as the auxiliary feed-water system operating normally. The calculated results obtained through this present software are in good agreements with that by RELAP5. In addition, some typical accidents, such as the double loops loss-of-flow accident, the reactivity insertion accident, the sudden increase of feed-water temperature accident, the loss of offsite power accident are also verified and also show to be in good agreements with results obtained by RELAP5. This proves the correctness of the mathematical models, the validity of the numerical method and the reliability of the software. Also, this software is expected to be validated against experimental data.
Fig. 8. Comparison between the results calculated by TTHAsoft 1.0 and RELAP5 for the loss of feed-water accident. (A) Thermal power, (B) pressure and water level.
4. Calculation results and discussions 4.1. Double loops loss-of-flow accident Double loops loss-of-flow accident occurs when the two reactor coolant pumps idle at the same time due to mechanical or electrical failures. The protective measures will shut down the reactor immediately when there is a low-flow signal in any loops, and the main feed-water supply stops and the turbine trips at the same time. In the accident analysis, the bypass is not considered after the steam turbine trips, and the emergency shutdown will be only triggered by the low-flow signal. Fig. 9 presents the characteristic change curves of the typical parameters with time in the accident. During the initial stage of the accident, the mass flow rate in the primary loop decreases rapidly for the two reactor coolant pumps idle at the same time. However, the reactor power changes slightly and thus the ratio of the power to the mass flow rate increases. The heat leg exit quality and maximum pellet temperature increase while the DNBR decreases (to the minimum value, about 1.78 in about 5 s). The low-flow signal of the emergency shutdown is triggered at about 3.5 s and the power decreases rapidly which restrains the decrease of the DNBR. Since the decreasing speed of the power is larger than that of the mass flow rate, the hot leg exit quality and the maximum pellet temperature begin to decrease after reaching to the peak, the DNBR starts to increase, and then the accident has been mitigated. Therefore, the minimum DNBR not exceeding the safety limitation is the main safety criterion in the initial stage of the accident while the residual
120
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
Fig. 9. Curves for the double loops loss-of-flow accident.
heat in the core must be removed correctly in order to ensure the integrity of the core in long term. 4.2. Reactivity insertion accident Reactivity insertion accident may occur during the refueling, the unexpected withdrawing of the control rod or the improper regulation of boron. These events are considered the most severe accidents that could lead to core damage, because during such events the core become supercritical and the core power rises to level beyond the heat removal system’s capability. In this paper, the variations of the transient parameters are analyzed while a +0.11$ reactivity is inserted into the system as a result of the improper regulation of control rod or boron. Here, it is assumed that the accident does not cause the core shutdown or other reactor protection systems actions and the power is completely controlled by negative feedback effects caused by the Doppler action and moderator-density. As shown in Fig. 10, the power increases rapidly because of the sudden increase of the reactivity, then decreases gradually as a result of the negative feedback effect, and goes back to the original power level finally. For the external load does not change, the pressure in the secondary side of the steam generator increases and
the water level decreases a little. In the initial stage of the accident, the increase of the power makes the DNBR decrease. Under the present calculation conditions, the minimum DNBR obtained form the W-3 equation as mentioned in Section 2.5 drops to 1.71, which is still within the safety limitation. Therefore, the nuclear power system could adjust itself to a stable and safety state while a +0.11$ reactivity is inserted into the system. 4.3. Sudden increase of feed-water temperature accident All the protection, control and regulation systems are not put into operation except for the reactivity feedbacks in the sudden increase of the feed-water temperature accident. Fig. 11 shows the change curves of the reactor thermal power, the reactivity and the pressure in the pressurizer with time assuming a sudden increase of the feed-water temperature by 30 K. It is shown that the increase of the inlet feed-water enthalpy (temperature) of the secondary side of the steam generator will make the secondary side quality increase which would cause the increases of the pressure and the water level in the secondary side of the steam generator. The heat removal capability of the steam generator decreases. Thus, the system pressure and temperatures in the primary loop will increase. The negative feedback effect of the reactivity caused
Fig. 10. Curves for the reactivity insertion accident.
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
121
Fig. 11. Curves for the sudden increase of feed-water temperature accident.
by the increases of the system pressure and temperatures could make the power decrease to a stable power level finally. Generally speaking, the parameters of the primary and secondary loops change slightly under the disturbance of the feed-water temperature. 4.4. Loss of offsite power accident The electricity supply of the reactor coolant pumps is interrupted in the loss of offsite power accident and then the coolant mass flow rate decreases rapidly because the driving force of the coolant is replaced by the inertial rotation of the pump and the inertial flow of the coolant. At the same time, a signal triggered by the power loss would make the reactor shut down, the main feedwater stop and the main steam line close. The steam bypass is not considered after the main steam valve has been shut down and it is assumed that there is without delay among the shutdowns of the reactor, the pump, the main feed-water, the main steam line and the loss of offsite power. If the delay time of the reactor trip protection is too long, the decrease velocity of the coolant mass flow rate is larger than that of the heat flux on the fuel element surface. And it may result in the decrease of the safety margin and the increase of the fuel element surface temperatures which may lead to the clad meltdown. The water level in the steam generator decreases as a result of the steam evacuation in the secondary loop. If the feed-water is not supplied in time and the water level reaches to the tube sheet height, the heat transfer efficiency would reduce further which is worse for the primary loop heat removal. In this paper, the emergency feed-water is assumed to be supplied after 60 s to relieve the decrease of the water level. The thermal power of the reactor and the mass flow rate decrease rapidly after the loss of offsite power as shown in Fig. 12(A). However, in the initial stage, the decrease velocity of the coolant mass flow rate is less than that of the thermal power. Thus, the minimum DNBR increases as shown in Fig. 12(B). The minimum DNBR will decrease slightly when the coolant mass flow rate and the thermal power tend to be stable. As shown in Fig. 12(B), the water level in the steam generator decreases as a result of the steam evacuation in the secondary loop. However, the water level will increase and tend to be stable after the emergency feed-water is supplied. The consequences of the loss of offsite power accident are less severe than the double loops loss-of-flow accident.
Fig. 12. Curves for the loss of offsite power accident. (A) Thermal power and mass flow rate, (B) MDNBR and water level.
5. Conclusion An analysis software named TTHAsoft 1.0 has been developed using the Fortran 90 and Visual C++ languages to evaluate the transient thermal-hydraulic behaviors of the AC-600. The simulated results by TTHAsoft 1.0 were compared with those of the
122
Y.W. Wu et al. / Nuclear Engineering and Design 240 (2010) 112–122
RELAP5/MOD3 code. Although there were differences between some local parameters, which were caused by the adoption of different models in the two codes, a good agreement between the computational results by the two codes was obtained which proved the applicability and accuracy of this software. Also, this software is expected to be validated against experimental data. Because of the adoption of modular programming techniques, this analysis software is expected to be applied to other reactors by easily modifying the corresponding function modules. The software has been applied to analyze the transient thermal-hydraulic characteristics of the Qinshan nuclear power plant in China. The obtained analysis results can be used as reference of the improvement design of the AC-600 and the safety operation of the Qinshan nuclear power plant. References Abdallah, M.A., Ahmed, H.M., Mohammad, A.R., Mohammad, E.S.N., 1982. Pressurizer transients dynamic model. Nucl. Eng. Des. 73, 447–453. Baggovra, B., Martin, W.R., 1983. Transient analysis of the three-Mile island unit 2 pressurizer system. Nucl. Technol. 62, 407–416.
Baron, B.C., 1973. Digital simulation of a nuclear pressurizer. Nucl. Sci. Eng. 52, 283–291. Beak, S.M., 1986. A Non-equilibrium three region model for transient analysis of pressurized water reactor pressurizer system. Nucl. Technol. 74, 213–221. Chen, J.C., 1966. Process design development. Ind. Eng. Chem. 5, 322. Chisholm, D., 1967. Naval Electronics Laboratory Report, No. 310. Colebrook, C.F., 1939. Turbulent flow in pipes with particular reference to the transition region between the smooth and rough pipe laws. Proc. Inst. Civ. Eng. 11, 133. Collier, J.G., 1981. Convective Boiling and Condensation, second ed. McGraw-Hill. Collier, J.G., Thome, J.R., 1994. Convective Boiling and Condensation, 3rd ed. Oxford University Press. Gear, C.W., 1971. Numerical Initial Value Problems in Ordinary Differential Equation. Prentice-Hall, Englewood Cliffs, NJ. Hind marah, A.C., 1974. Gear: ordinary differential equation system solver, Lawrence Livarmore Laboratory, Report UCID-30001, Revision 3, December. Jens, W.H., Lottes, P.A., 1951. Argonne National Laboratory, Report No. 4627. Kang, S.W., 1984. Pool heat transfer in a simulated PWR pressurizer. Trans. Am. Nucl. Soc. 46, 845–854. McAdams, W.H., 1954. Heat Transmission, third ed. McGraw-Hill, New York. Redfoeld, J.A., 1968. Pressurizer performance during loss of load tests at shippingport: analysis and test. Nucl. Appl. 4, 189–196. Rohsenow, W.M., Hartnett, J.R., 1973. Hand Book of Heat Transfer. McGraw-Hill, New York. Thomas, L.C., 1992. Heat Transfer. Prentice Hall, New York. Wilson, J.F., 1962. The velocity of rising steam in a bubble two phase mixture. Am. Nucl. Soc. 5, 151–163.