Research Article  Open Access
Soroush Heidari Sangestani, Mohammad Rahgoshay, Naser Vosoughi, Mitra Athari Allaf, "Study of Fast Transient Pressure Drop in VVER1000 Nuclear Reactor Using Acoustic Phenomenon", Science and Technology of Nuclear Installations, vol. 2018, Article ID 7862847, 11 pages, 2018. https://doi.org/10.1155/2018/7862847
Study of Fast Transient Pressure Drop in VVER1000 Nuclear Reactor Using Acoustic Phenomenon
Abstract
This article aims to simulate the sudden and fast pressure drop of VVER1000 reactor core coolant, regarding acoustic phenomenon. It is used to acquire a more accurate method in order to simulate the various accidents of reactor core. Neutronic equations should be solved concurrently by means of DRAGON 4 and DONJON 4 coupling codes. The results of the developed package are compared with WIMS/CITATION and final safety analysis report of Bushehr VVER1000 reactor (FSAR). Afterwards, time dependent thermalhydraulic equations are answered by employing Single Heated Channel by Sectionalized Compressible Fluid method. Then, the obtained results were validated by the same transient simulation in a pressurized water reactor core. Then, thermalhydraulic and neutronic modules are coupled concurrently by use of producing group constants regarding the thermal feedback effect. Results were compared to the mentioned transient simulation in RELAP5 computer code, which show that mass flux drop is sensed at the end of channel in several milliseconds which causes heat flux drop too. The thermal feedback resulted in production of some perturbations in the changes of these parameters. The achieved results for this very fast pressure drop represent accurate calculations of thermoneutronic parameters fast changes.
1. Introduction
One of the most important aspects in design of different safety systems with sufficient preparation is simulation and analysis of transient states in reactor core. For these kinds of analysis basic equations in neutronic and thermalhydraulic modules have to be coupled. Coupling of neutronic and thermalhydraulic modules is different from each other, considering numerical solution methods and time and body meshing size. Thus, written codes for different transient states are mostly used for study of fuel and coolant temperature changes, power peak level, coolant pressure, stability time, and other parameters. Computer coding submits models with different degree of accuracy and validity. Most codes are not able to analyze coupled conditions of very fast transient (FT) states in very short time. This deficiency is associated with neutronic and thermalhydraulic calculations or both. Therefore designing a code which is capable of analyzing FT conditions is highly needed.
One of the efficient ways in analyzing FT is using waveform method. Chan (1991) studied asymptotic waveform evaluation in analysis of timedependent calculations [1]. Ooi et al. (2003) studied the finite element method (FEM) in thermal analysis that usually produces a formulation in the space/time domain. However, the sizes of the equations in FEM usually are large, and thus the conventional algorithms involve considerable computational time. The conventional methods have to take a very small time step size to avoid undesirable numerically induced oscillations or numerical instabilities. Thus, a new solution algorithm, named the asymptotic waveform evaluation scheme, was introduced by them to solve transient problems [2]. Liu et al. (2006) studied fast thermal simulation when power density increases by fast spectrum in frequency domain for computing steady state response. It indicates that the resulting thermal analysis algorithms lead to 10x–100x speedup over the traditional integrationbased transient analysis with small accuracy loss. The studied parameters minimum time order is 15 ms [3]. Ham and Bathe (2012) used FEM to solve timedependent twodimensional wave propagation problems [4]. Ranaa et al. (2014) presented the well condition asymptotic waveform evaluation to solve heat conduction problem in the frequency domain. The method is presented for timedependent problems [5]. In addition, Ishii et al. (2009) investigated the effect of acoustic phenomenon in steam dryer (in a Boiling water reactor) which indicates the process of pressure pulsations caused by hydroacoustic resonance propagation along the steam dryer [6]. Proskuryakov (2017) studied the effect of acoustic vibrations in the nuclear power plant (NPP) coolant such as VVER1000 reactor and created new scientific direction “diagnosis, prognosis and prevention of vibration  acoustic resonances in the NPP equipment” shows that the developed methods can be used to predict and prevent the occurrence of vibrationacoustic resonances in the NPP equipment [7].
Various kinds of thermalhydraulic and neutronic calculating models are put to use in transient calculating code development studies. Reducing costs and runtime and achieving required accuracy are three main purposes of them. Leung et al. (1981) studied acoustic impact techniques for increasing the accuracy of FT states modeling in CODA code [8]. However, regarding its very small meshing, acoustic methods are timeconsuming and are not used any more.
In order to decrease complicated solving parameters in using compressible fluid method, geometry of every fuel assembly could be turned into a single heated channel (SHC). Four different methods are used to solve SHC transient equations. They are channel integral (CI), single mass velocity (SV), momentum integral (MI), and sectionalized compressible fluid (SC). SC model considers both sound impact and thermal expansion, while the other three models only consider thermal expansion [9]. The SHC method is used in costanza code (1967) in order to analyze the dynamics of liquid cooled nuclear reactor [10]. It is also used in DynCo code (2011) which is intended for complex neutronphysic and thermalhydraulic dynamic calculation of the reactor core in 3D hexZ geometry. DynCo code is designed to simulate the dynamic behavior of the Russian 3000MWt pressurized water reactor (PWR) [11]. The coupling of CI and POWEXK code (2011) simulates power excursion in reactor of Budapest University of Technology and Economics Reactor [12]. Hosseini et al. (2015) developed a coupled 3D neutronic with 1D thermalhydraulic model in order to find reasonable power distribution. Neutronic module includes threedimensional diffusion equation in two energy groups which was solved using analytic nodal method. Thermalhydraulic module contains the conservation equations solved for 1D axial homogeneous downward flow through channel using SHC model (MI method). In conclusion, Xenon saturation transient analysis of a research reactor core was carried out [13].
The SHC compressible fluid method is one of the ways that make it possible to use wave propagation method and acoustic phenomenon. Shapiro (1953) considered the derivation and properties of the dynamics and thermodynamics of compressible fluid flow which is along with acoustic theory [14]. Meyer (1961) investigated that several approximations can be used to decouple the momentum and energy equations to facilitate solution of the transient problem. Additionally, the numerical solution of a transient problem would be particularly simplified if the compressibility of the fluid could be ignored [15]. Todreas and Kazimi (2001) investigated a rapid increase in the heat flux without change in the applied pressure drop in a PWR. The SC approach for a 10% heat flux step increase in the PWR channel shows that because the pressure begins to rise at internal channel points, a reduction in inlet flow rate and an accelerated exit flow rate occur [9].
The mechanism of compressible fluid method was published by BarMeir in 2007 [16]. Khola and Pandey (2013) studied the numerical simulation of transients in twophase flow which is crucial to simulate accidentlike conditions of nuclear reactors for safety analysis. In their work, a code for numerical computation of unsteady onedimensional twophase flow has been developed. The governing equations were obtained by the homogenous equilibrium mixture model and were decoupled and approximated using the SC and MI model [17].
Numerical considerations (i.e., the stability and/or accuracy) of the difference solution require that the time step of integration be less than the time interval for sonic wave propagation across the spatial grid points. Compared with the transport velocity, the fluid sonic velocity is large. It causes limitation of the time step in most numerical schemes to very small values. Acoustic phenomenon is used in accommodation of very short time and very small body meshing. This accommodation is determined by Courant’s criterion [8, 20]:
is time meshing period, is body meshing distance, is sound velocity in coolant, and is the mean transport velocity. The fluid sonic velocity is large comparing with the transport velocity. Therefore, limiting the time step in most numerical schemes to very small values leads to a computationally expensive solution of this problem. Regarding meshing size criterion, pressure wave propagation method mostly requires longterm computations for numeral stability. Therefore acoustic phenomenon features (despite the high accuracy) are not usually considered by researchers. However these features are very useful during core parameter vital changes.
FT pressure drop is a type of loss of coolant accident (LOCA). That is well known as the double ended guillotine break. When double ended guillotine break occurs (main coolant pipeline cold leg breaks at the reactor inlet), suddenly reactor coolant pressure decreases with leak coolant flow rate of 45000 Kg/s [18]. A pressure wave is also produced, while the large break LOCA occurs, which propagates across the channel. Therefore International Atomic Energy Agency (IAEA, 2003) studied the beforebreak vital moment (leak before break) [21].
The coolant fast pressure drop accident can lead the fluid to become twophased and thermoneutronic parameters to change. GonzalezSantalo and Lahey Jr. (1972) investigated this matter by study of pressure drop transient in twophase condition [22].
Calculation program of VVER1000 reactor core FT pressure drop by means of SC method and acoustic phenomenon was developed in this investigation. In order to use the mentioned method every one of the 28 fuel assemblies should be considered as one SHC. Fuel assembly conversion into SHC and meshing method toward axis are both shown in Figure 1.
2. Materials and Methods
2.1. An Overview of Neutronic and ThermalHydraulic Model
A computational program for simulation of VVER1000 reactor core FT state of sudden pressure drop has been developed in this study. It includes the two thermalhydraulic and neutronic basic models. Group constants of the fuel assemblies and reflectors are produced by DRAGON 4 code [23]. DRAGON 4 is a lattice code developed to solve Boltzmann transport equation in two and three dimensions, to apply selfshielding effects and to compute fewgroup macroscopic cross sections and diffusion coefficients. The mentioned program includes different models for simulation of fuel assembly behavior, such as interpolation of microscopic cross sections, resonance selfshielding calculation, different solvers for the Boltzmann transport equation with ability to take into account leakage effects, and calculation of condensed and homogenized parameters.
In this study, the integral transport equation is solved by the SYBILT module, selfshielding calculations are performed by the SHI module by means of generalized Stamm’ler method, and the CPO module is utilized for production of equivalent fuel assembly parameters in consistent format that can be used in forgoing calculation. After that, timedependent multigroup neutron diffusion equations are solved by DONJON 4 models [24]. DONJON 4 is a multigroup diffusion solver. A threedimensional simulation is performed by Thomas–Raviart–Schneider method using The TRIVAT module, the group constants of fuel assembly calculated by DRAGON 4 are recovered using the CRE module, and the FLUD module is to compute multiplication factor. The FEM is utilized for discretization of equation. Finally momentary amount of power distribution is achieved in 1/6 of reactor core, across every fuel assembly.
In thermalhydraulic model, conservation of mass and momentum and energy dependent equations are solved by applying considered transients and use of channel axial division in MATLAB software. The channel is regarded as a typical coolant subchannel inside an assembly that only receives coolant through its bottom inlet. The fuel and clad heat transport equations are solved separately from coolant equations. Onedimensional transient transport equations of the coolant with radial heat input from the clad surfaces are resolved. The flow area is assumed axially uniform but pressure loss due to local area changes is still regarded. Any axial position of flow area could be considered as control area in order to derive radially averaged coolant flow equations. Mentioned equation is envisioned for SHC state and SC method of every 28 fuel assemblies (1/6 core). From the four SHC different methods only SC model considers both sound impact and thermal expansion, while the other three models considering thermal expansion [9]. Therefore it is able to study very FT states. Some of the thermodynamic features of water are also regarded [25].
Several approximations can be used to decouple the momentum and energy equations to facilitate solution of this transient. Additionally, the numerical solution of the discussed transient problem in this study is particularly simplified following the work of Meyer [15] and specific applications from Todreas and Kazimi [9]. Therefore the lateral variations of properties in the flow channels can be neglected. For this condition, (2), (3), and (4) are directly applicable:where is enthalpy, is mass flux, is density, is fluid level passing, is of fuel assembly height, is time, is friction coefficient, is heat flux, is thermalhydraulic diameter, is speed sound, is pressure, and is heated surface. Equations (2), (3), and (4) and appropriate constitutive relations provide the solutions for , , and for the initial and boundary conditions. The heat flux , which is a DONJON output, is specified as an input for the above equations. The boundary conditions for solving the momentum and energy equations are to be specified as and at the inlet and at the outlet. Furthermore, constitutive equations for and are required to complete definition of the problem. The equation of state for the density, assumed differentiable with respect to and , is specified asThe friction factor can be specified as
In SC model, numerical solution approach involves a set of difference equations representing the differential transport equations, arranged to consider , , and , and state variables at multiple points along the channel. Every one of the 28 fuel assemblies is considered as a single heated channel. The term sectionalized reflects the need to divide a channel into 360 segments to execute the numerical solution, regarding the height of every fuel assembly (3.6 m). Using (5) we getwhereFrom (2) and (7) we getEquations (4) and (9) may be combined into two equations by eliminating and :where we have defined as
Equations (3) and (10) are partial differential equations in , , and (i.e., the density does not appear as a differentiated variable). These equations were written in pointwise difference form (forward implicit method) and solved in , , and .
Sound velocity is different in various spaces; also the fluid is running in a specified direction with a specified speed. Fluid velocity is accordingly dependent on its other parameters. Regime pressure drop transient numerical considerations (i.e., the stability and/or accuracy) of the difference solution need that the time step of integration be less than the time interval for sonic wave propagation across the spatial grid points, that is, (1). The time discretization that is used in most codes varies from the semiimplicit scheme and multistep schemes or the fully implicit discretization.
Different codes have developed various methods for predicting theses states. CATHARE code [26] used critical flow rate correlation method; that way the geometry of the flow restriction is simplified and a coarse meshing is used. TRACE [27] and RELAP5 [28] codes use characteristic velocity; that way the geometry of the flow restriction is also simplified and a coarse meshing is used. A sound characteristic velocity is set to zero and simplified equations are used to predict flow evolution from upstream to the sonic section. The standard option in RELAP5 is a semiimplicit scheme. Only those governing parameters are evaluated at new time values which are needed to maintain stability and accuracy of the method or to avoid Courant time step size limitations based on pressure waves. CATHARE code [26] in another method assumes that the flow from upstream to sonic section is precisely calculated by 2fluid equations using a very fine meshing in the vicinity of the throat. Such a method is only possible with an implicit numerical scheme which allows high velocity flow in small meshes without material Courant limitation.
All these codes, in which Courant criterion is not used, consider different ways in order to decrease the time of calculations but still they do not have the accuracy of the codes using Courant criterion, in time of critical flow conditions and very short times. Therefore Courant criterion was taken into account in the simulating codes of this study. Time mesh amount is obtained about 10 μs considering for solution stability and sound velocity of 900 m/s in coolant. DRAGON and DONJON/SHC developed code is regarded 1 μs for more accuracy.
2.2. Very FT Utility Analysis
Very FT coolant pressure drop in VVER1000 reactor core is explained in
It is highly mentioned that axial fluid motion is considered steady in SC model which eliminates cross flow effect [9]. By considering core symmetry, one of the six parts of core is modeled, which includes 28 fuel assemblies’ channels. Therefore, despite the ones where core is assumed as one channel, radial accuracy of thermalhydraulic and neutronic parameters is increased. The impact of thermal feedback on cross sections is utilized for coupling the two mentioned modules. Eventually standard deviation of predicted results for noted transient is estimated through accuracy evaluation of thermalhydraulic and neutronic codes. Some advantages of SHC method by SC model are increasing the sufficient accuracy in sudden transient state simulations and more accurate analysis of different transient states. As mentioned before, acoustic phenomenon increases the number of time meshes, and therefore the time of calculation becomes longer. On that account, performing in FT states is the optimum condition for using noted code (considering user demands).
2.3. Coupling Procedure
Coupling algorithm of the two thermalhydraulic and neutronic modules is as follows. Primarily every fuel assembly is modeled by DRAGON 4 cell calculating code according to Table 1. This code is used (considering fresh fuel) in order to determine the macroscopic cross sections sequences and other neutronic parameters of all reactor fuel assemblies. Instead of modeling the whole core, 28 fuel assemblies are modeled according to core symmetry. The next step is obtaining axial and radial power distributions by means of available triangular geometry module of DONJON 4 and numeral solving by means of FEM, also applying VVER1000 reactor power. Considering nonpoisonous fresh fuel, reactor relative power radial distribution is measured in order to evaluate the calculation algorithm authenticity of the two codes’ correlation and validated by comparing to WIMS/CITATION [19] results and Bushehr VVER1000 reactor FSAR [18].

In order to evaluate the authenticity of thermalhydraulic module, the mentioned transient is modeled in PWR core fuel assembly (conditions are noted in Table 2) and its obtained results are compared to the Todreas and Kazimi [9] study results. Core is composed of 28 fuel assemblies; but only Hot Fuel Assembly (HFA) coupling results are noted in this assay. Thermalhydraulic properties of coolant fast pressure drop accident in PWR core is simulated by SC method. The amounts of time periods which are caused by environmental sound velocity limitations are also considered for current simulation. After evaluating the authenticity of thermalhydraulic module, coupling results of the two neutronic and thermalhydraulic modules of coolant fast pressure drop accident in VVER1000 reactor core could be studied for validation which is compared to RELAP5 (advancedthermalhydraulic code, coupled with point kinetic model to simulate neutronic behavior) [28]. The thermalhydraulic and neutronic coupling module algorithm is as follows: the procedure starts with obtaining the axial power profile of neutronic module in steady state. The mentioned axial power profile and other core dependent parameters (sudden pressure drop transient, mass flux, pressure, temperature changes, etc.) are the entering data of thermalhydraulic module program. Axial thermal changes affect the nuclear cross sections and the new power profile in neutronic module is produced. Then the new axial power profile is applied on thermalhydraulic module. It is continued until the difference between axial thermal changes in every time step is less than a specified amount of convergence criteria. Neutronic and thermalhydraulic coupling calculations algorithm are shown in Figure 2.

As a result of time period selection limitations (short time periods), which leads to increase in the number of nodes, the rise in parameter numbers of defined study case and complicated equations or scales should not be considered. Therefore program solving time is lowered.
3. Results and Discussion
3.1. Results of Neutronic Calculations
Relative radial power distribution in 1/6 of VVER1000 reactor core is shown in Table 3. Power Peaking Factor is compared to WIMS/CITATION [19] and Bushehr VVER1000 reactor FSAR value [18], as observed. Obtained results show an unremarkable calculation error for the two neutronic coupled codes.

The maximum coupling calculation errors of DRAGON 4/DONJON 4 is 2.8% while it is 5.64% in WIMS/CITATION [19] codes. The obtained results of DRAGON 4 and DONJON 4 codes coupling shows higher accuracy than WIMS/CITATION codes as observed. Then relative axial power distribution in 1/6 of VVER1000 reactor core is analyzed, which is shown in Figure 3. As observed, the peak of relative axial power distribution (in beginning of cycle) occurs lower than the middle height of core. Relative axial power distribution is more coordinated to Bushehr VVER1000 reactor FSAR [18] when using DRAGON 4/DONJON 4 coupling rather than WIMS/CITATION coupling [19].
3.2. Results of ThermalHydraulic Calculations
After the defined model is validated, thermalhydraulic module preparation is started based on SHC by SC method. In order to evaluate the authenticity of newly made module, coolant fast pressure drop accident in PWR core fuel assembly is simulated by SC method (PWR core features are noted in Table 2). This simulation also considers different amounts of time periods which are caused by environmental limitations of sound velocity. Severe decreases in entrance pressure (also neutronic changes stability) cause mass flux reduction which is distributed in the channel during spent time [9, 29]. The outputs of thermalhydraulic program, which are normalized mass flux changes during axial side of the channel, are shown in Figure 4. Sound velocity is calculated 900 m/s, causing channel transient state to last for 4.1 ms. The mass flux fast reduction of ms would not reach the end of the channel; therefore upper sides will not be affected by pressure wave.
For ms, the reflected pressure wave affects the mass flux. It happens because of stable outer pressure theory, which means when opposite pole wave (with similar amplitude) moves in opposite direction, the diluted entrance wave is reflected like compressed wave in outer border. Pressure profile in ms and longer times shows the effect of reflection wave progress all the way from exit to entrance. In that short amount of time, average mass flux rate only decreases from 1 to 0.98. Figure 4 shows the results of program performance and has an acceptable accommodation with reference diagram [9], as observed.
The written program accuracy in neutronic states is compared to Bushehr VVER1000 reactor FSAR data [18]. Also it is compared to PWR data [9] in thermalhydraulic states. After the authenticity evaluation of neutronic and thermalhydraulic modules, correlated FT states and changes of various parameters during time are studied.
3.3. Results of Coupling Calculations
Changing diagrams of thermalhydraulic parameters by using sound effect in pressure propagation, for very FT coolant pressure drop in VVER1000 reactor HFA, are being studied. Pressure drop causes the mass flux reduction; therefore temperature increases and power drops in every node. The transfer of each node’s thermal data to DRAGON 4 module leads to the decrease of cross section moderation. Mentioned data and axial power profile are used as DONJON 4 module entrance. Finally lowered power distribution and pressure drop transient are applied as thermalhydraulic module and thermal distribution of every one of 360 nodes is obtained again. The process continues until the convergence criteria of each node thermal difference are available. Considering the central pressure of each node and reactor core pressure changes, which is 0.3 MPa, alterations of pressure per time for the first node is shown in Figure 5; that reveals pressure drop of 0.0012 MPa per 1 ms.
Figure 6 shows mass flux changes during fuel assembly in three different times. Increasing the time from 10 μs to 2 ms causes a mass flux drop up to more than half of the fuel assembly. Mass flux drop in 10 μs, which is 0.35 Kg/m^{2}s, occurs up to the beginning of fuel assembly. However, these changes are not sensed by the upper nodes. After 1 ms, mass flux drop arrives at the middle of fuel assembly, which increases from 0.35 Kg/m^{2}s to 12 Kg/m^{2}s for the first node and decreases to 25 Kg/m^{2}s after 2 ms. The channel transient time of sonic wave is about 4 ms. The rapid decrease of the mass flux has not yet reached the end of the channel in this period of time; that means the upstream region is not yet affected by the pressure wave.
Mass flux changes in 5 ms indicate the wave reflux from end of the channel. From now on, for a few milliseconds, the reversible wave causes more mass flux drop in the end of channel than the beginning of it. Mass flux changes in 5 and 9 ms, which are calculated by RELAP5 [28], are shown in the diagram. RELAP5 does not sense pressure drop transient in times before 4.2 ms, and therefore it has some delay comparing to DRAGON 4 and DONJON 4/SHC codes and lesser time periods are not possible to be compared.
The observed difference between DRAGON 4 and DONJON 4/SHC diagram (in 5 ms) and RELAP5 diagram (in 5 ms) is because of the difference between neutronic and thermalhydraulic modules accuracy. Also the effects of thermal feedback cause perturbation in mass flux alterations, in DRAGON 4 and DONJON 4/SHC code. RELAP5 [28] results only show the average mass flux drop caused by pressure drop. Average mass flux changes with steady slope have four reasons:(1)Driving force drop caused by inlet pressure drop that reduces the channel average mass flux.(2)The fluid thermal expansion due to heating which causes a small slope change; that is, becomes less than .(3)The fluids being considered uncompressible (i.e., ).(4)Weak neutronic model and not using sufficient time meshes sizes for the mentioned transient.
The obtained results from mass flux changes simulation in primary, middle, and ending points of relative length of fuel assembly, in 9 ms, in DRAGON 4 and DONJON 4/SHC and RELAP5 code are shown in Figure 7. Mass flux in the first node initially faces a severe drop which continues for 2 ms. When this drop becomes more visible in upper nodes, diagrams slope turns more gradual along with average drop of 0.068 (Kg/m^{2}s). Mass flux in central point of fuel assembly starts dropping after whereabouts 1.2 ms. However severe drop for ending points occurs after 4 ms. The profiles at ms and more show the effect of the reflected wave progressing toward the channel exit from the inlet. During this short period of time, the mass flux only decreases 60 Kg/m^{2}s. The net mass flux is the superposition of the forward wave and the reflected wave.
RELAP5 results show mass flux profile changes, with steady slope, in primary, middle, and ending points and does not present adequate description of the passing pressure wave.
HFA heat flux changing per time for beginning, middle, and ending nodes of fuel assembly, in “0–3 ms” period of time, is shown in Figure 8. The heat flux of first node drops whereabouts 40 KW/m^{2} in that amount of time. However, mass flux of fuel assembly middle node drops whereabouts 30 KW/m^{2} in 1.7 ms. It also occurs after 2.3 ms for ending node. Temperature changes per time for the first node are shown in Figure 9, which face drop of 0.01°C in 0.06 ms after the beginning of the accident. Mentioned drop is inevitable, based on constant axial power decrease.
4. Conclusions
Simulation of VVER1000 reactor core FT pressure drop, by using acoustic phenomenon, was studied in this article. Coupling of the two DRAGON 4 and DONJON 4 neutronic modules, along with thermalhydraulic module by means of SC method, was utilized in this simulation. It is the first time that neutronic and thermalhydraulic modules are externally coupled, assuming that the coolant fluid is compressible, in order to obtain the best simulation for this transient.
In this transient, at the channel inlet, the mass flux begins to decrease because of the pressure reduction; it was observed in calculations that a perturbation was produced and propagates into the channel as time elapses. After pressure waves arrive at the end of the channel, reflected pressure waves affect the mass flux profile for short time. The incoming rarefaction wave will be reflected as a compression wave, due to the assumption of constant outlet pressure, at the exit boundary. Therefore a wave travels in the opposite direction with the same amplitude but opposite sign. The reverse wave causes more mass flux drop in upper nodes and as time passes, it decreases the difference between upper nodes mass flux drop and the ones in beginning of the channel. Thus it is concluded that, in such transients, mass flux drop rate does not follow stable conditions.
This transient (double ended guillotine break) is very fast, and therefore DRAGON 4 and DONJON 4 codes are used in neutronic module (comparing to WIMS/CITATION) considering their high speed and accuracy in analyzing the problem. Also SC method is utilized in thermalhydraulic module for solutions of conservation equations by forward implicit method and Courant’s criterion in a SHC. Thus the very short time step (μs grades) was used in analyzing transient. The neutronic module was validated by WIMS/CITATION and FSAR. Thermalhydraulic module was also validated by simulation of a PWR (comparing to Todreas and Kazimi [9] simulated one). Finally the simulated pressure drop transient (by DRAGON and DONJON/SHC code) was also simulated in RELAP5 code and results were compared. RELAP5 code was not capable of sensing the mentioned transient in times before 4.2 ms and, unlike DRAGON and DONJON/SHC code, coolant fluid is assumed uncompressible in it. So changes of density per pressure is considered zero. Of course this assumption is highly reasonable for classes of core transients that do not have the problem of losing a great amount of coolant. Therefore mass flux changes slope is almost steady along the channel and mass flux changes are only dependent on pressure changes in channel inlet and outlet. Thus it became completely obvious that, despite the fact that RELAP5 code uses “nearly implicit” scheme based on a “fractional step” approach for providing a sufficient degree of implicitness to eliminate the material Couranttype stability limit, but it does not give an adequate description of the mentioned studied transient.
The obtained results from DRAGON and DONJON/SHC for the times that a great volume of coolant is lost in a short time, assuming coolant fluid is compressible, are highly acceptable and they show its ability to effectively and accurately calculate the fast transient. On account of that, the results are not desirable if the sonic effect in coolant is not regarded.
Abbreviations
FT:  Fast transient 
FEM:  Finite element method 
PWR:  Pressurized water reactor 
CI:  Channel integral 
SV:  Single mass velocity 
MI:  Momentum integral 
SC:  Sectionalized, compressible fluid 
NPP:  Nuclear power plant 
FSAR:  Final safety analysis report 
LOCA:  Loss of coolant accident 
SHC:  Single heated channel 
HFA:  Hot fuel assembly 
ms:  Millisecond 
μs:  Microsecond 
IAEA:  International Atomic Energy Agency. 
Conflicts of Interest
All authors have no conflicts of interest to declare.
Acknowledgments
This study was performed as part of longterm research into nuclear safety supported by the Nuclear Engineering Department of Science and Research Branch, Islamic Azad University. The valuable researches into various accidents of VVER1000 reactor led us into performing this study.
References
 P. K. Chan, “Comments on “Asymptotic Waveform Evaluation for Timing Analysis”^{1},” IEEE Transactions on ComputerAided Design of Integrated Circuits and Systems, vol. 10, no. 8, pp. 10781079, 1991. View at: Publisher Site  Google Scholar
 C. Ooi, K. Seetharamu, Z. Alauddin, G. Quadir, K. Sim, and T. Goh, “Fast transient solutions for heat transfer [FEM],” in Proceedings of the IEEE TENCON 2003. Conference on Convergent Technologies for the AsiaPacific Region, pp. 469–473, Bangalore, India, 2003. View at: Publisher Site  Google Scholar
 P. Liu, H. Li, L. Jin, W. Wu, S. X.D. Tan, and J. Yang, “Fast thermal simulation for runtime temperature tracking and management,” IEEE Transactions on ComputerAided Design of Integrated Circuits and Systems, vol. 25, no. 12, pp. 2882–2892, 2006. View at: Publisher Site  Google Scholar
 S. Ham and K. Bathe, “A finite element method enriched for wave propagation problems,” Computers & Structures, vol. 9495, pp. 1–12, 2012. View at: Publisher Site  Google Scholar
 M. S. Ranaa, K. Jeevanb, R. Harikrishnanc, and A. W. Rezad, “A wellcondition asymptotic waveform evaluation method for heat conduction problems,” Advanced Materials Research, vol. 845, pp. 209–215, 2014. View at: Publisher Site  Google Scholar
 K. Ishii, T. Nagayoshi, S. Takahashi, Y. Wada, and N. Tanikawa, “Advanced simulation technologies for nuclear power plants,” Hitachi Review, vol. 58, no. 2, pp. 94–101, 2009. View at: Google Scholar
 K. N. Proskuryakov, “Scientific Basis for Modeling and Calculation of Acoustic Vibrations in the Nuclear Power Plant Coolant,” Recent Advances in Petrochemical Science, vol. 2, no. 1, 2017. View at: Google Scholar
 J. C. M. Leung, K. A. Gallivan, and R. E. Henry, Critical Heat Flux Predictions During Blow down Transient, Argonne National Laboratory, Argonne, IL60439, U.S.A, 1981.
 N. E. Todreas and M. S. Kazimi, “Thermal Hydraulic Fundamentals,” Massachusetts Institute of Technology, 1990. View at: Google Scholar
 G. Forti and E. Vincenti, “The codes costanza for the dynamics of liquid cooled nuclear reactor, joint nuclear research center Ispra stablishment  Italy,” Reactor Physics Department Reactor Teory and Analysis, 1967. View at: Google Scholar
 A. V. Levchenko, I. N. Leonov, and S. V. Kovalchuk, “Verification of the three dimensional dynamic DYNCO code on the basis of international test cases for pressurized water reactors,” Simulation Systems Ltd, 2011. View at: Google Scholar
 S. M. Khaleda and F. Al Mutairic, “A Numerical Technique For Solving Coupled Thermal hydraulic and MultiEnergy Group Neutron Diffusion Equations,” A Department of Studies and Basic Sciences, Community College, University of Tabuk, Saudi Arabia, 2011. View at: Google Scholar
 M. Hosseini, H. Khalafi, and S. Khakshournia, “Two group, three dimensional diffusion theory coupled with single heated channel model: A study based on xenon transient simulation of Tehran research reactor,” Progress in Nuclear Energy, vol. 85, pp. 108–120, 2015. View at: Publisher Site  Google Scholar
 A. H. Shapiro, The Dynamics and Thermodynamics of Compressible Fluid Flow, vol. 1, Ronald Press Co., New York, NY, USA, 1953.
 J. E. Meyer, “Hydrodynamic models for the treatment of reactor thermal transients,” Nuclear Science and Engineering, vol. 10, pp. 269–277, 1961. View at: Google Scholar
 G. BarMeir, “Fundamentals of Compressible Fluid Mechanics, 2007,” Version 0.4.4.2 aka 0.4.4.1j, May 21. Minneapolis, MN 554142411. View at: Google Scholar
 N. Khola and M. Pandey, “Numerical computation of onedimensional unsteady twophase flow using hem model and IAPWS IF97 equations of state,” in Proceedings of the 2013 21st International Conference on Nuclear Engineering, ICONE 2013, China, August 2013. View at: Publisher Site  Google Scholar
 AEOI., Reactor Final Safety Analysis Report VVER1000 Bushehr, Chapter 4, Atomic Energy Organization of Iran, 2005.
 A. H. Fadaei and S. Setayeshi, “Control rod worth calculation for VVER1000 nuclear reactor using WIMS and CITATION codes,” Progress in Nuclear Energy, vol. 51, no. 1, pp. 184–191, 2009. View at: Publisher Site  Google Scholar
 J. Plesek, R. Kolman, and D. Gabriel, “Estimation of The Critical Time Step for Explicit Integration,” in Proceedings of the 18th International Conference, Engineering Mechanics, Svratka, Czech Republic, 2012. View at: Google Scholar
 IAEATECDOC1361, Assessment and management of ageing of major nuclear power plant components important to safety Primary piping in PWRs, Scientific and Technical Publications, Vienna, Austria, 2003.
 J. M. GonzalezSantalo and R. T. Lahey Jr., An Exact Solution of Flow Decay Transients in TwoPhase Systems by The Method of Characteristics, The American Society of Mechanical Engineers, New York, NY, USA, 1972.
 G. Marleau, A. Hébert, and R. A. Roy, “A user guide for DRAGON Version 4,” Institute of Genius Nuclear, Department of Genius Mechanical, School Polytechnic of Montreal, 2011. View at: Google Scholar
 D. Sekki, A. Hébert, and R. Chambon, “A User Guide for Donjon 4 Version 4,” Institute of Genius Nuclear, Department of Genius Mechanical, School Polytechnic of Montreal, 2011. View at: Google Scholar
 W. Wagner and H.J. Kretzschmar, International Steam Tables, Faculty of Mechanical Engineering Chair of Thermodynamics, 2nd edition, 2007.
 G. Lavialle, “The CATHARE 2 V2.5 code: Main features, CATHARENEPTUNE International Seminar,” Grenoble, May 2004. View at: Google Scholar
 TRACEV5.0, Theory Manual, 2007, TRACE V5.0 User Manual, 2007, TRACE V5.0 Assessment Manual, 2007.
 RELAP5/SCDAP//MOD3.2 Code Manuals, A Computer Code for BestEstimate Transient Simulation of Light Water Reactor Coolant Systems During Severe Accidents, Prepared for the U.S. Nuclear Regulatory Commission, Idaho National Engineering and Environmental Laboratory, NUREG/CR6150, 1997.
 M. Rahgoshay, Encyclopedia of Advanced Subjects on Nuclear Heat Transfer, vol. 1, Chapter 3, Islamic Azad University Central Tehran Branch, Iran, 2016, ISBN: 9789641043898.
Copyright
Copyright © 2018 Soroush Heidari Sangestani et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.