# Theoretical Study on the Characteristics of Critical Heat Flux in Rectangular Channel of Natural Circulation under Motion Conditions

Authors: Mengmeng Xi; Rong Cai; Xiao Chu; Fan Yang; Xu Ran; Zhifang Qiu; Jian Deng; Wenxi Tian
Journal: Science and Technology of Nuclear Installations
 (2020)

---

## Abstract

With the wide application of sea-based reactors, the impact of ocean conditions on the safety performance of reactors has gradually attracted attention. In this paper, by establishing the thermal hydraulic transient analysis model and the critical heat flux (CHF) model of natural circulation system, the CHF characteristics in the rectangular channel of natural self-feedback conditions under ocean conditions are studied. The results show that the additional acceleration field generated by ocean conditions will affect the thermal hydraulic parameters of the natural circulation system, that is, the external macroscopic thermal hydraulic field. On the other hand, the boiling crisis mechanism will be affected, that is, the force on the bubble and the thickness of the liquid film. Within the parameters of the study, ocean conditions have a great impact on CHF of natural circulation, and the maximum degradation of CHF is about 45%. The obtained analysis results are significant to the improvement of design and safety operation of the reactor system.


---

## Body

## 1. Introduction

Critical heat flux (CHF) is one of the important safety guidelines for reactor thermal hydraulic design. The CHF is the heat flux when the heat transfer on the surface of the nuclear fuel element deteriorates. It is the result of the transition of the coolant flow boiling mechanism. According to the parameters such as the enthalpy and vapor content of the fluid, the DNB type and the Dryout type of boiling crisis phenomena may appear in the tube. If the boiling crisis occurs on the surface of the fuel element of the reactor, the surface temperature of the fuel element will be too high, causing the heated wall to burn out, and the radioactive material leaks, causing serious operational accidents, directly affecting the safety of the reactor. Compared with conventional energy-powered propulsion ships, the performance of nuclear-powered propulsion ships has been greatly improved, and the maneuverability, concealability, and endurance of nuclear-powered ships have been greatly improved. In order to meet maneuverability, ship reactors need to change power frequently, and if necessary, natural circulation is a reliable way to remove heat without any driving components such as pumps.With the application of sea-based reactors such as marine nuclear power platforms, the impact of ocean conditions on the operation of the reactor has gradually attracted attention. Isshiki [1] studied the effects of heaving and heeling on the CHF under normal pressure and low flow conditions. The results show that the CHF decreases linearly with the increase of the acceleration. The experimental observation of the ocean conditions causes the bubble distribution changes, and a semiempirical relationship is obtained by introducing the fluctuation influence factor and the heaving influence factor. OTSUJI et al. [2–4] used Freon 113 to study the effect of oscillating acceleration on the CHF of forced convection. It has been found that the acceleration of oscillating generally reduces the CHF, and the cavitation and flow fluctuation caused by oscillating acceleration will cause CHF to appear in advance, and OTSUJI believes that the mechanism of this CHF phenomenon is the evaporation of the liquid film under the steam block. There were only a few cases of CHF experiments under ocean conditions [5]. Due to the complexity of ocean conditions, researchers usually calculate CHF by modifying empirical relationships.There are many CHF prediction methods including empirical correlations, look-up tables, and various mechanistic models. Although CHF in boiling flow pipe was investigated by lots of researchers [6–10], they are all under forced circulation system and stationary condition. Umekawa et al. [11] experimentally and numerically studied Dryout CHF in a boiling channel under oscillatory flow condition. The numerical simulation represents well the transient behavior of the Dryout under the oscillatory flow condition. Geping et al. [12] experimentally studied CHF and Dryout point in narrow annuli pipe, and based on experiment data, a new correlation is developed. Celata et al. [13] presents an analytical model for the prediction of the CHF in water saturated flow boiling in round vertical and uniformly heated pipes. The model shows a quite good CHF predictive capability for selected experimental points. Du et al. [6] investigated the CHF in a vertical narrow rectangular channel. The influences of system pressure, inlet mass flux, and channel size on CHF in rectangular channels are analyzed. The CHF can be predicted well, compared with experimental data. Jayanti and Valette [14] presented a set of closure relations of a one-dimensional three-phase calculation methodology for the prediction of Dryout and post-Dryout heat transfer at high pressure conditions and developed new correlations for the interfacial friction factor between the liquid film and the gas phase, for the droplet diameter and for the transition to annular flow. Su et al. [10] developed a theoretical three-fluid model for predicting annular upward flow in vertical narrow annuli with bilateral heating. The present model can be used to calculate the CHF and critical quality in narrow annular gap.It can be seen from the above-mentioned research status at home and abroad that natural circulation, as an important form of circulation in nuclear power systems, is of great significance in the operation of reactors. The influence of ocean conditions on natural circulation is more significant than that on forced circulation. Under ocean conditions, the value of the critical heat flux has decreased significantly compared to the value under land-based conditions, and it will bring serious challenges to the reactor safety. The mechanism of CHF of natural circulation under ocean conditions is more complicated, but related studies are rarely published.This paper takes the natural circulation system under ocean conditions as the research object, establishes a reasonable mathematical physics model for the natural circulation system structure and dynamic self-feedback characteristics, and develops the analysis program using the Fortran 90 programming language. The thermal hydraulic characteristics of the natural circulation system and the CHF characteristics under natural circulation conditions are calculated and analyzed.


## 2. Mathematical Model

The schematic structure of the natural circulation loop is shown in Figure1. The natural circulation loop includes a heated section, a condenser, a preheater, a pressurizer, and pipes.Figure 1
Natural circulation schematic diagram.The basic field model for natural circulation is based on fundamental conservation principles: the mass, momentum, and energy conservation equations. With the assumption of one-dimensional flow, these equations, including single-phase and two-phase conservative equations, are introduced in detail as follows.


### 2.1. Thermal Hydraulic Model Equations

In the natural circulation loop, the calculation of transient thermal hydraulic parameters is solved by the basic conservation equations, including mass conservation equation, momentum conservation equation, and energy conservation equation.Mass conservation equation:(1)∂ρ∂t+1A∂W∂z=0.Momentum conservation equation:(2)∂W∂t+∂∂zW2ρA=−A∂p∂z−ρgA−fW22ρDeA.Energy conservation equation:(3)ρ∂H∂t+WA∂H∂z=qUhA.


### 2.2. Pressurizer Model

It is very important to accurately simulate the dynamic characteristics of the pressurizer to improve the simulation accuracy of the system. The natural circulation system adopts nitrogen pressurizer. According to the different phases and enthalpies of the fluid in the pressurizer, the pressurizer is divided into two regions, namely, the water region and the nitrogen region. And the following simplification assumptions are given:(1)
The two regions share the same pressure at the same time.(2)
The condensation and evaporation of the water at the interface of the gas-liquid interface are ignored.(3)
The surge water through the surge pipe is well proportioned with the water in the pressurizer.(4)
The same fluid has the same thermodynamic parameters at the same area and time.The detailed equations of the pressurizer can be found in the published original paper [15].


### 2.3. Pipe and Plenum Model

The adiabatic assumption was applied to the pipe and plenums in the natural circulation system. According to the thermal delay model, transient characteristic of coolant enthalpy in pipes and plenums was given as follows:(4)Mdhdt=Whin−h.


### 2.4. Ocean Motion Model

A sea-based reactor will be affected by ocean conditions (waves, wind, and tide) in working conditions. The ocean conditions mainly include heeling, heaving, rolling, and pitching motions (shown in Figure2). Coupled motions, which are coupling of more than one type of ocean condition, are also considered. In heeling motion, there is no additional acceleration or inertial. The forces’ effect on the flow is similar to that in steady state, except the gravity. In rolling and pitching motions, only the additional force in flowing direction should be considered. In this paper, it is assumed that the ocean motions are in sinusoidal order:(5)θ=θmsin2πtT,where θ is rolling or pithing angle, and θm and T are rolling or pitching amplitude and period, respectively. t is time. Taking derivatives, the angular velocity and angular acceleration for a fluid volume under rolling motion can be written as(6)ω=2πθmTcos2πtT,(7)β=−4π2θmT2sin2πtT.Figure 2
Schematic of ocean motions. (a) Heeling. (b) Heaving. (c) Rolling. (d) Pitching.
(a)(b)(c)(d)The heaving motion is similar to the law of motion.(8)d=dmsin2πtT,where d is heaving displacement, and dm and T are heaving amplitude and period, respectively. Taking derivatives, the velocity and acceleration for a fluid volume under heaving motion can be written as(9)v=2πdmTcos2πtT,a=−4π2dmT2sin2πtT.In the natural circulation system, pipes are three-dimensional. When the system is under rolling motion, the additional force on the fluid in the tube, which is perpendicular to the rolling plane, is omitted.In the progress of system calculation, the influence of ocean motion means adding additional terms in the momentum equation. According to [16], the momentum equation can be written as follows:(10)∂ρur∂t+∂ρurur∂z=−∂p∂z−fDeρur22+F+ρf⋅k,where F+ρf⋅k is the part that is influenced by ocean motions including the additional pressure and variant elevation pressure. The ocean motion additional pressure can be expressed as follows:(11)F=−ρa0+2ω×u+ω×ω×r+dωdt×r,where −2ω×ur is the Coriolis force, −ω×ω×r is the centrifugal force, −dω/dt×r is the tangential force, and a0 is additional acceleration.


### 2.5. DNB Type CHF Model

For the prediction of DNB type CHF, the calculation model based on the microliquid layer evaporation mechanism is adopted. The model is improved based on the Lee and Mudawar model [17]. The model assumes that, in the near wall area, a small block of steam is formed due to the combination of small bubbles, and the small bubble polymerization increases only the length of the block in the axial direction (Figure 2). There is a very thin layer of liquid between the block and the wall, called the microliquid layer. Boiling crisis occurs when the microliquid layer below is completely depleted due to evaporation during the time the steam is swept:(12)qCHF=ρlδhfgUBLB.In the above formula,δ, UB, and LB represent the distance between the steam block and the heating wall, the moving speed of the steam block, and the length of the steam block, respectively. Solving these three parameters, we can find the CHF value qCHF under the corresponding working conditions.The computational model of this study still uses the basic idea of the Lee and Mudawar model.δ and UB are calculated by the radial and axial force balance of the bubble, while the LB is calculated by the majority of the microliquid layer model. That is, the Helmholtz wavelength is used. The calculation of δ and UB is given below.The velocityUB of the steam block in the direction of fluid flow is affected by buoyancy and drag forces in addition to the local liquid phase velocity.(13)FBa+FD=0,FBa=π4DB2LBρl−ρggcosθ+aat,FD=−12ρlCDUB−UBL2πDB24,where FBa is the buoyancy of the bubble, the unit N. FD is the drag force of the bubble, the unit N. ρl − ρg is the difference of vapor-liquid two-phase density, the unit kg/m3. UBL is the liquid velocity of the local radial position of the bubble, the unit m/s. g is the gravitational acceleration, the unit m/s2. θ is the angle of the fluid channel inclination, the unit °. aa is the axial component of the additional acceleration, the unit m/s2. CD is the drag coefficient. From equations (5)–(7), the steam block speed UB can be solved.The radial forces experienced by the steam block include the evaporation forceFE, the side lift force FL, the wall lubrication force FWL, and the radial buoyancy FBr. The evaporation force is formed by the vapor phase evaporating from the microliquid layer to the steam block. The lifting force is generated by the liquid boundary layer velocity gradient and the relative speed difference between the two phases. The wall lubrication force pushes the bubble near the wall toward the central region.The reason for radial buoyancyFBr is that the motion of the fluid channel generates an additional acceleration field, and the additional acceleration field has a component field in the radial direction (vertical to the flow direction). The two-phase fluid in the channel flows in this radial acceleration field. The density difference between the two phases causes the bubble to undergo radial buoyancy. The direction of the buoyancy depends on the direction of the additional radial acceleration field and is opposite to it.The effect of the bubble on the distance from the wall is very small. The force process is shown in Figure3. The expressions of the various forces are as follows. When the radial force is balanced, the distance δ between the block and the wall can be determined.(14)FE+FWL+FBr−FL=0,FE=DBLBq2ρghfg2,FL=CLρlUB−UBL∂UL∂yπ4DB2LB,FWL=LBπDB24CWLρlUB−UBL2DB,FBr=ρl,sat−ρgartLBπDB24.Figure 3
Conceptual sketch of liquid sublayer Dryout mechanism.In the above formulas,CL is the lifting force coefficient, (∂UL)/∂y represents the mainstream velocity gradient of the position of the steam block, which is solved by the Karman velocity equation, and CWL is the wall lubrication coefficient.For the DNB type CHF solution model, other constitutive relations are specifically described in the literature [18].


### 2.6. Dryout Type CHF Model

In an annular flow regime of narrow rectangular channels, the liquid flows upward adhered to the channel walls forming a continuous annular liquid film. The vapor with entrained droplets flows along the center of the channel and forms a continuous vapor core. Mass is continuously exchanged between the liquid film and droplets entrained in the vapor core in the flow direction: droplets entrained in the vapor core deposit to the liquid film, and the vapor core of high velocity flow cause droplets entrainment in the liquid film.


#### 2.6.1. Mass Conservation Equations

(1)
The mass conservation equation of the liquid filmIn the annular region, there exists mass transfer at the liquid film-vapor core interface by the complicated effects on film evaporation, droplets entrainment, and deposition. The mass conservation equation of the liquid film is shown as follows:(15)∂∂τρfαf+∂ρfαfuf∂z=−qPrqAhfg+PrwDep−PrwEnh+PrqEnqA,whereρf is the density of the liquid film, αf the void fraction of the liquid film, uf is the average flow velocity of the liquid film, q is the heat flux, hfg is the latent heat of vaporization, A is the cross-sectional area of rectangular channel, Prw is the channel periphery, Prq is the heating periphery, Dep is the liquid droplets deposition rate, Enh is the liquid droplets entrainment rate caused by breakup of disturbance waves on the liquid film-vapor core interface, and Enq is the liquid droplets entrainment rate induced by wall heat flux.(2)
The mass conservation equation of dropletsDroplets come from the liquid film evaporation, droplets entrainment, and deposition,(16)∂∂τρdαd+∂ρdαdud∂z=−PrwDep−PrwEnh+PrqEnqA,whereρd is the density of droplets, αd is the void fraction of droplets, and ud is the average flow velocity of droplets.(3)
The mass conservation equation of the vapor coreThe vapor core of the annular flow region mainly comes from the liquid film evaporation,(17)∂∂τρgαg+∂ρgαgug∂z=qPrqAhfg,where ρg is the density of the vapor core, αg is the void fraction of the vapor core, and ug is the average flow velocity of the vapor core.


#### 2.6.2. Momentum Conservation Equations

(1)
The momentum conservation equation of the liquid filmAccording to the law of momentum conservation on the liquid film control volume in which frictional force between the liquid film and droplets is ignored, the momentum conservation equation is deduced as follows:(18)∂ρfαfuf∂τ+∂ρfαfuf2∂z−DepPrwAud+qPrqAhfguf+EnhPrw+EnqPrqAuf=−αf∂p∂z+αfF+ρff⋅k−Mwf+Mfg,whereMwf is the wall friction, and Mfg is the vapor-liquid interface friction.(2)
The momentum conservation equation of dropletsThe forces of droplets include gravity, pressure, and frictional force between the vapor core and droplets. According to the forces of droplets and momentum exchange between droplets and the liquid film, the momentum conservation equation of droplets is deduced as follows:(19)∂ρdαdud∂τ+∂ρdαdud2∂z+DepPrwAud−EnhPrw+EnqPrqAuf=−αd∂p∂z+αdF+ρdf⋅k+Mgd,whereMgd is the friction between the vapor core and droplets.(3)
The momentum conservation equation of the vapor coreThe forces of the vapor core include gravity, pressure, frictional force between the vapor core and droplets, and frictional force between the vapor core and the liquid film. According to the force of the vapor core and momentum exchange between the vapor core and the liquid film, the momentum conservation equation of the vapor core is deduced as follows:(20)∂ρgαgug∂τ+∂ρgαgug2∂z−qPrqAhfguf=−αg∂p∂z+αgF+ρgf⋅k−Mfg−Mgd.


#### 2.6.3. The Energy Conservation Equation

To simplify the three-fluid model, it is assumed that the liquid film, droplets, and the vapor core of annular upward flow are in saturated condition.For the Dryout type CHF solution model, other constitutive relations are specifically described in the literature [19].


### 2.7. Auxiliary Model

The correct heat transfer and flow friction factor correlations are key aspects to the code. According to the operation conditions of the natural circulation system, the collier correlation is chosen for laminar heat transfer (Re < 2000), the Sieder–Tate correlation is chosen for the turbulent flow heat transfer (Re > 3000), and the linear interpolation is applied between the laminar and turbulent correlations for the transition flow heat transfer (2000 ≤ Re ≤ 3000). The Jens and Lottes correlation is adopted to calculate the subcooled boiling heat transfer, and the Chen correlation is chosen for saturated boiling heat transfer.The laminar flow friction factor is calculated byf=C/Re (C = 64 for circular tube, C = 96 for rectangular channel), the turbulent flow friction factor is calculated by Blasius correlation, and the transition flow friction is calculated by linear interpolation. The two-phase frictional multiplier is calculated by Lockhart–Martinelli correlation. The main heat transfer and flow conditions that might occur in the transient process of natural circulation system and corresponding optional correlations are also listed in Tables 1 and 2.Table 1
Flow regimes and heat transfer correlations.
Flow regimeCorrelationSingle-phase water:Laminar flowCollier correlationTransition flowInterpolationTurbulent flowSieder–Tate correlationTwo-phase:Subcooled boilingJens and Lottes correlationSaturated boilingChen correlationTable 2
Flow regimes and flow friction correlations.
Flow regimeCorrelationSingle-phase:Laminar flowf=C/Re (C = 64 for circular tube, C = 96 for rectangular channel)Transition flowInterpolationTurbulent flowBlasius modelTwo-phase:Mixture model (two-phase multiplier)Lockhart–Martinelli modelQuality and void fraction at the onset of annular flow should be given to determine the flow pattern. The typical two-phase flow patterns in a vertical tube are shown in Figure4. For a uniformly heated channel, the flow patterns contain single-phase flow, subcooled boiling flow, bubbly flow, churn flow, and annular flow. In the present work, the boundary between the subcooled and saturated flow boiling regions is the location where the equilibrium quality equals zero. And the onset of annular flow is the location where the mass quality is evaluated by the widely used Wallis [20] correlation:(21)xann=0.6+0.4gDeρf−ρgρf/G0.6+ρf/ρg.Figure 4
Schematic of flow patterns in a heated vertical tube.


### 2.8. Thermophysical Property Model

The thermophysical properties of water and steam are calculated using the correlations from the international standard IAPWS-IF97.


## 2.1. Thermal Hydraulic Model Equations

In the natural circulation loop, the calculation of transient thermal hydraulic parameters is solved by the basic conservation equations, including mass conservation equation, momentum conservation equation, and energy conservation equation.Mass conservation equation:(1)∂ρ∂t+1A∂W∂z=0.Momentum conservation equation:(2)∂W∂t+∂∂zW2ρA=−A∂p∂z−ρgA−fW22ρDeA.Energy conservation equation:(3)ρ∂H∂t+WA∂H∂z=qUhA.


## 2.2. Pressurizer Model

It is very important to accurately simulate the dynamic characteristics of the pressurizer to improve the simulation accuracy of the system. The natural circulation system adopts nitrogen pressurizer. According to the different phases and enthalpies of the fluid in the pressurizer, the pressurizer is divided into two regions, namely, the water region and the nitrogen region. And the following simplification assumptions are given:(1)
The two regions share the same pressure at the same time.(2)
The condensation and evaporation of the water at the interface of the gas-liquid interface are ignored.(3)
The surge water through the surge pipe is well proportioned with the water in the pressurizer.(4)
The same fluid has the same thermodynamic parameters at the same area and time.The detailed equations of the pressurizer can be found in the published original paper [15].


## 2.3. Pipe and Plenum Model

The adiabatic assumption was applied to the pipe and plenums in the natural circulation system. According to the thermal delay model, transient characteristic of coolant enthalpy in pipes and plenums was given as follows:(4)Mdhdt=Whin−h.


## 2.4. Ocean Motion Model

A sea-based reactor will be affected by ocean conditions (waves, wind, and tide) in working conditions. The ocean conditions mainly include heeling, heaving, rolling, and pitching motions (shown in Figure2). Coupled motions, which are coupling of more than one type of ocean condition, are also considered. In heeling motion, there is no additional acceleration or inertial. The forces’ effect on the flow is similar to that in steady state, except the gravity. In rolling and pitching motions, only the additional force in flowing direction should be considered. In this paper, it is assumed that the ocean motions are in sinusoidal order:(5)θ=θmsin2πtT,where θ is rolling or pithing angle, and θm and T are rolling or pitching amplitude and period, respectively. t is time. Taking derivatives, the angular velocity and angular acceleration for a fluid volume under rolling motion can be written as(6)ω=2πθmTcos2πtT,(7)β=−4π2θmT2sin2πtT.Figure 2
Schematic of ocean motions. (a) Heeling. (b) Heaving. (c) Rolling. (d) Pitching.
(a)(b)(c)(d)The heaving motion is similar to the law of motion.(8)d=dmsin2πtT,where d is heaving displacement, and dm and T are heaving amplitude and period, respectively. Taking derivatives, the velocity and acceleration for a fluid volume under heaving motion can be written as(9)v=2πdmTcos2πtT,a=−4π2dmT2sin2πtT.In the natural circulation system, pipes are three-dimensional. When the system is under rolling motion, the additional force on the fluid in the tube, which is perpendicular to the rolling plane, is omitted.In the progress of system calculation, the influence of ocean motion means adding additional terms in the momentum equation. According to [16], the momentum equation can be written as follows:(10)∂ρur∂t+∂ρurur∂z=−∂p∂z−fDeρur22+F+ρf⋅k,where F+ρf⋅k is the part that is influenced by ocean motions including the additional pressure and variant elevation pressure. The ocean motion additional pressure can be expressed as follows:(11)F=−ρa0+2ω×u+ω×ω×r+dωdt×r,where −2ω×ur is the Coriolis force, −ω×ω×r is the centrifugal force, −dω/dt×r is the tangential force, and a0 is additional acceleration.


## 2.5. DNB Type CHF Model

For the prediction of DNB type CHF, the calculation model based on the microliquid layer evaporation mechanism is adopted. The model is improved based on the Lee and Mudawar model [17]. The model assumes that, in the near wall area, a small block of steam is formed due to the combination of small bubbles, and the small bubble polymerization increases only the length of the block in the axial direction (Figure 2). There is a very thin layer of liquid between the block and the wall, called the microliquid layer. Boiling crisis occurs when the microliquid layer below is completely depleted due to evaporation during the time the steam is swept:(12)qCHF=ρlδhfgUBLB.In the above formula,δ, UB, and LB represent the distance between the steam block and the heating wall, the moving speed of the steam block, and the length of the steam block, respectively. Solving these three parameters, we can find the CHF value qCHF under the corresponding working conditions.The computational model of this study still uses the basic idea of the Lee and Mudawar model.δ and UB are calculated by the radial and axial force balance of the bubble, while the LB is calculated by the majority of the microliquid layer model. That is, the Helmholtz wavelength is used. The calculation of δ and UB is given below.The velocityUB of the steam block in the direction of fluid flow is affected by buoyancy and drag forces in addition to the local liquid phase velocity.(13)FBa+FD=0,FBa=π4DB2LBρl−ρggcosθ+aat,FD=−12ρlCDUB−UBL2πDB24,where FBa is the buoyancy of the bubble, the unit N. FD is the drag force of the bubble, the unit N. ρl − ρg is the difference of vapor-liquid two-phase density, the unit kg/m3. UBL is the liquid velocity of the local radial position of the bubble, the unit m/s. g is the gravitational acceleration, the unit m/s2. θ is the angle of the fluid channel inclination, the unit °. aa is the axial component of the additional acceleration, the unit m/s2. CD is the drag coefficient. From equations (5)–(7), the steam block speed UB can be solved.The radial forces experienced by the steam block include the evaporation forceFE, the side lift force FL, the wall lubrication force FWL, and the radial buoyancy FBr. The evaporation force is formed by the vapor phase evaporating from the microliquid layer to the steam block. The lifting force is generated by the liquid boundary layer velocity gradient and the relative speed difference between the two phases. The wall lubrication force pushes the bubble near the wall toward the central region.The reason for radial buoyancyFBr is that the motion of the fluid channel generates an additional acceleration field, and the additional acceleration field has a component field in the radial direction (vertical to the flow direction). The two-phase fluid in the channel flows in this radial acceleration field. The density difference between the two phases causes the bubble to undergo radial buoyancy. The direction of the buoyancy depends on the direction of the additional radial acceleration field and is opposite to it.The effect of the bubble on the distance from the wall is very small. The force process is shown in Figure3. The expressions of the various forces are as follows. When the radial force is balanced, the distance δ between the block and the wall can be determined.(14)FE+FWL+FBr−FL=0,FE=DBLBq2ρghfg2,FL=CLρlUB−UBL∂UL∂yπ4DB2LB,FWL=LBπDB24CWLρlUB−UBL2DB,FBr=ρl,sat−ρgartLBπDB24.Figure 3
Conceptual sketch of liquid sublayer Dryout mechanism.In the above formulas,CL is the lifting force coefficient, (∂UL)/∂y represents the mainstream velocity gradient of the position of the steam block, which is solved by the Karman velocity equation, and CWL is the wall lubrication coefficient.For the DNB type CHF solution model, other constitutive relations are specifically described in the literature [18].


## 2.6. Dryout Type CHF Model

In an annular flow regime of narrow rectangular channels, the liquid flows upward adhered to the channel walls forming a continuous annular liquid film. The vapor with entrained droplets flows along the center of the channel and forms a continuous vapor core. Mass is continuously exchanged between the liquid film and droplets entrained in the vapor core in the flow direction: droplets entrained in the vapor core deposit to the liquid film, and the vapor core of high velocity flow cause droplets entrainment in the liquid film.


### 2.6.1. Mass Conservation Equations

(1)
The mass conservation equation of the liquid filmIn the annular region, there exists mass transfer at the liquid film-vapor core interface by the complicated effects on film evaporation, droplets entrainment, and deposition. The mass conservation equation of the liquid film is shown as follows:(15)∂∂τρfαf+∂ρfαfuf∂z=−qPrqAhfg+PrwDep−PrwEnh+PrqEnqA,whereρf is the density of the liquid film, αf the void fraction of the liquid film, uf is the average flow velocity of the liquid film, q is the heat flux, hfg is the latent heat of vaporization, A is the cross-sectional area of rectangular channel, Prw is the channel periphery, Prq is the heating periphery, Dep is the liquid droplets deposition rate, Enh is the liquid droplets entrainment rate caused by breakup of disturbance waves on the liquid film-vapor core interface, and Enq is the liquid droplets entrainment rate induced by wall heat flux.(2)
The mass conservation equation of dropletsDroplets come from the liquid film evaporation, droplets entrainment, and deposition,(16)∂∂τρdαd+∂ρdαdud∂z=−PrwDep−PrwEnh+PrqEnqA,whereρd is the density of droplets, αd is the void fraction of droplets, and ud is the average flow velocity of droplets.(3)
The mass conservation equation of the vapor coreThe vapor core of the annular flow region mainly comes from the liquid film evaporation,(17)∂∂τρgαg+∂ρgαgug∂z=qPrqAhfg,where ρg is the density of the vapor core, αg is the void fraction of the vapor core, and ug is the average flow velocity of the vapor core.


### 2.6.2. Momentum Conservation Equations

(1)
The momentum conservation equation of the liquid filmAccording to the law of momentum conservation on the liquid film control volume in which frictional force between the liquid film and droplets is ignored, the momentum conservation equation is deduced as follows:(18)∂ρfαfuf∂τ+∂ρfαfuf2∂z−DepPrwAud+qPrqAhfguf+EnhPrw+EnqPrqAuf=−αf∂p∂z+αfF+ρff⋅k−Mwf+Mfg,whereMwf is the wall friction, and Mfg is the vapor-liquid interface friction.(2)
The momentum conservation equation of dropletsThe forces of droplets include gravity, pressure, and frictional force between the vapor core and droplets. According to the forces of droplets and momentum exchange between droplets and the liquid film, the momentum conservation equation of droplets is deduced as follows:(19)∂ρdαdud∂τ+∂ρdαdud2∂z+DepPrwAud−EnhPrw+EnqPrqAuf=−αd∂p∂z+αdF+ρdf⋅k+Mgd,whereMgd is the friction between the vapor core and droplets.(3)
The momentum conservation equation of the vapor coreThe forces of the vapor core include gravity, pressure, frictional force between the vapor core and droplets, and frictional force between the vapor core and the liquid film. According to the force of the vapor core and momentum exchange between the vapor core and the liquid film, the momentum conservation equation of the vapor core is deduced as follows:(20)∂ρgαgug∂τ+∂ρgαgug2∂z−qPrqAhfguf=−αg∂p∂z+αgF+ρgf⋅k−Mfg−Mgd.


### 2.6.3. The Energy Conservation Equation

To simplify the three-fluid model, it is assumed that the liquid film, droplets, and the vapor core of annular upward flow are in saturated condition.For the Dryout type CHF solution model, other constitutive relations are specifically described in the literature [19].


## 2.6.1. Mass Conservation Equations

(1)
The mass conservation equation of the liquid filmIn the annular region, there exists mass transfer at the liquid film-vapor core interface by the complicated effects on film evaporation, droplets entrainment, and deposition. The mass conservation equation of the liquid film is shown as follows:(15)∂∂τρfαf+∂ρfαfuf∂z=−qPrqAhfg+PrwDep−PrwEnh+PrqEnqA,whereρf is the density of the liquid film, αf the void fraction of the liquid film, uf is the average flow velocity of the liquid film, q is the heat flux, hfg is the latent heat of vaporization, A is the cross-sectional area of rectangular channel, Prw is the channel periphery, Prq is the heating periphery, Dep is the liquid droplets deposition rate, Enh is the liquid droplets entrainment rate caused by breakup of disturbance waves on the liquid film-vapor core interface, and Enq is the liquid droplets entrainment rate induced by wall heat flux.(2)
The mass conservation equation of dropletsDroplets come from the liquid film evaporation, droplets entrainment, and deposition,(16)∂∂τρdαd+∂ρdαdud∂z=−PrwDep−PrwEnh+PrqEnqA,whereρd is the density of droplets, αd is the void fraction of droplets, and ud is the average flow velocity of droplets.(3)
The mass conservation equation of the vapor coreThe vapor core of the annular flow region mainly comes from the liquid film evaporation,(17)∂∂τρgαg+∂ρgαgug∂z=qPrqAhfg,where ρg is the density of the vapor core, αg is the void fraction of the vapor core, and ug is the average flow velocity of the vapor core.


## 2.6.2. Momentum Conservation Equations

(1)
The momentum conservation equation of the liquid filmAccording to the law of momentum conservation on the liquid film control volume in which frictional force between the liquid film and droplets is ignored, the momentum conservation equation is deduced as follows:(18)∂ρfαfuf∂τ+∂ρfαfuf2∂z−DepPrwAud+qPrqAhfguf+EnhPrw+EnqPrqAuf=−αf∂p∂z+αfF+ρff⋅k−Mwf+Mfg,whereMwf is the wall friction, and Mfg is the vapor-liquid interface friction.(2)
The momentum conservation equation of dropletsThe forces of droplets include gravity, pressure, and frictional force between the vapor core and droplets. According to the forces of droplets and momentum exchange between droplets and the liquid film, the momentum conservation equation of droplets is deduced as follows:(19)∂ρdαdud∂τ+∂ρdαdud2∂z+DepPrwAud−EnhPrw+EnqPrqAuf=−αd∂p∂z+αdF+ρdf⋅k+Mgd,whereMgd is the friction between the vapor core and droplets.(3)
The momentum conservation equation of the vapor coreThe forces of the vapor core include gravity, pressure, frictional force between the vapor core and droplets, and frictional force between the vapor core and the liquid film. According to the force of the vapor core and momentum exchange between the vapor core and the liquid film, the momentum conservation equation of the vapor core is deduced as follows:(20)∂ρgαgug∂τ+∂ρgαgug2∂z−qPrqAhfguf=−αg∂p∂z+αgF+ρgf⋅k−Mfg−Mgd.


## 2.6.3. The Energy Conservation Equation

To simplify the three-fluid model, it is assumed that the liquid film, droplets, and the vapor core of annular upward flow are in saturated condition.For the Dryout type CHF solution model, other constitutive relations are specifically described in the literature [19].


## 2.7. Auxiliary Model

The correct heat transfer and flow friction factor correlations are key aspects to the code. According to the operation conditions of the natural circulation system, the collier correlation is chosen for laminar heat transfer (Re < 2000), the Sieder–Tate correlation is chosen for the turbulent flow heat transfer (Re > 3000), and the linear interpolation is applied between the laminar and turbulent correlations for the transition flow heat transfer (2000 ≤ Re ≤ 3000). The Jens and Lottes correlation is adopted to calculate the subcooled boiling heat transfer, and the Chen correlation is chosen for saturated boiling heat transfer.The laminar flow friction factor is calculated byf=C/Re (C = 64 for circular tube, C = 96 for rectangular channel), the turbulent flow friction factor is calculated by Blasius correlation, and the transition flow friction is calculated by linear interpolation. The two-phase frictional multiplier is calculated by Lockhart–Martinelli correlation. The main heat transfer and flow conditions that might occur in the transient process of natural circulation system and corresponding optional correlations are also listed in Tables 1 and 2.Table 1
Flow regimes and heat transfer correlations.
Flow regimeCorrelationSingle-phase water:Laminar flowCollier correlationTransition flowInterpolationTurbulent flowSieder–Tate correlationTwo-phase:Subcooled boilingJens and Lottes correlationSaturated boilingChen correlationTable 2
Flow regimes and flow friction correlations.
Flow regimeCorrelationSingle-phase:Laminar flowf=C/Re (C = 64 for circular tube, C = 96 for rectangular channel)Transition flowInterpolationTurbulent flowBlasius modelTwo-phase:Mixture model (two-phase multiplier)Lockhart–Martinelli modelQuality and void fraction at the onset of annular flow should be given to determine the flow pattern. The typical two-phase flow patterns in a vertical tube are shown in Figure4. For a uniformly heated channel, the flow patterns contain single-phase flow, subcooled boiling flow, bubbly flow, churn flow, and annular flow. In the present work, the boundary between the subcooled and saturated flow boiling regions is the location where the equilibrium quality equals zero. And the onset of annular flow is the location where the mass quality is evaluated by the widely used Wallis [20] correlation:(21)xann=0.6+0.4gDeρf−ρgρf/G0.6+ρf/ρg.Figure 4
Schematic of flow patterns in a heated vertical tube.


## 2.8. Thermophysical Property Model

The thermophysical properties of water and steam are calculated using the correlations from the international standard IAPWS-IF97.


## 3. Code Design and Calculation Procedure

### 3.1. Nodalization

To numerically simulate the thermal-hydraulic characteristics of natural circulation system according to different geometrical structures and physical boundary conditions, the heated section, preheater, and condenser were divided into different control volumes and junctions, as shown in Figure5. The node number of the system can be changed according to the calculation requirement.Figure 5
Schematic diagram of the nodalization of the system.


### 3.2. Numerical Method

After the nodalization of the natural circulation system, the transient analysis can be obtained to solve a set of ordinary differential equations with initial conditions, which can be summarized as follows:(22)dy⟶dt=ft,y→,y→',y⟶t0=y⟶0,where y⟶ is a vector including mass flow rate, enthalpy, and other parameters of each control volumes.Gear method equipped with Adams predictor-corrector method was adopted for better solution [21].


### 3.3. Code Design

Based on the natural circulation system model and CHF model and numerical methods mentioned above, the code has been developed using FORTRAN 90 language to analyze the CHF characteristics of natural circulation system under ocean conditions. For the purpose of convenient modification and secondary development, modular programming techniques were adopted. The analysis code consists of several functional modules: data input module, initialization module, numerical method module, derivative module, CHF model, and data output module. All the above function modules are called by main program (Figure6). Furthermore, the initialization module calls the thermophysical property module and the auxiliary model module. The derivative module calculates the right hand of the differential equations and call the thermophysical property module and the auxiliary model module. The flowchart of the code is also given (Figure 7).Figure 6
Call relationships of modules.Figure 7
Program flowchart.


### 3.4. The Calculation Procedure

In this paper, according to the CHF under natural circulation conditions of ocean conditions, according to the different mechanisms of CHF of DNB and Dryout, the CHF mechanism model of DNB and Dryout is established under the corresponding ocean conditions. The DNB type and Dryout type CHF prediction programs are compiled based on the model. Finally, the DNB type and Dryout type CHF prediction models are coupled with the natural circulation system program by judging the boiling crisis region. In the coupled model, the natural circulation system model is used as the main program framework, and the DNB type and Dryout boiling crisis mechanism models are used as subroutines in the main program framework. The ocean condition module provides boundary conditions for the natural circulation system and boiling crisis calculation. Natural circulation system provides external thermal hydraulic environment boundary for boiling crisis calculation. The boiling crisis feeds back the calculation results to the natural circulation system. The relationship between the models is shown in Figure8.Figure 8
Model coupling diagram.The flow of the CHF prediction coupling procedure under natural conditions of ocean conditions is as follows. Figure7 shows the flowchart of the calculation procedure.(1)
First, the natural circulation system program reads geometric parameters, control parameters, and initialization parameters.(2)
Set an initial heated section heat flux densityqm, start steady state calculation, and also adjust the preheating section, condenser parameters, and the resistance coefficient of the valve in the loop to reach the targeted system flow rate, the inlet temperature of heated section.(3)
After the steady state calculation is completed, the steady state result is used as the initial condition to start the transient calculation. The transient calculation obtains the mass quality of the outlet of the heated section for each time step, and the starting point mass quality of the annular flow under the current working condition. By comparing the two mass qualities, we can know whether there is an annular flow in the heated section.(4)
If there is no annular flow in the heated section, the DNB type boiling crisis mechanism model is called for calculation. If annular flow occurs in the heated section, the Dryout type boiling crisis is called to calculate.


## 3.1. Nodalization

To numerically simulate the thermal-hydraulic characteristics of natural circulation system according to different geometrical structures and physical boundary conditions, the heated section, preheater, and condenser were divided into different control volumes and junctions, as shown in Figure5. The node number of the system can be changed according to the calculation requirement.Figure 5
Schematic diagram of the nodalization of the system.


## 3.2. Numerical Method

After the nodalization of the natural circulation system, the transient analysis can be obtained to solve a set of ordinary differential equations with initial conditions, which can be summarized as follows:(22)dy⟶dt=ft,y→,y→',y⟶t0=y⟶0,where y⟶ is a vector including mass flow rate, enthalpy, and other parameters of each control volumes.Gear method equipped with Adams predictor-corrector method was adopted for better solution [21].


## 3.3. Code Design

Based on the natural circulation system model and CHF model and numerical methods mentioned above, the code has been developed using FORTRAN 90 language to analyze the CHF characteristics of natural circulation system under ocean conditions. For the purpose of convenient modification and secondary development, modular programming techniques were adopted. The analysis code consists of several functional modules: data input module, initialization module, numerical method module, derivative module, CHF model, and data output module. All the above function modules are called by main program (Figure6). Furthermore, the initialization module calls the thermophysical property module and the auxiliary model module. The derivative module calculates the right hand of the differential equations and call the thermophysical property module and the auxiliary model module. The flowchart of the code is also given (Figure 7).Figure 6
Call relationships of modules.Figure 7
Program flowchart.


## 3.4. The Calculation Procedure

In this paper, according to the CHF under natural circulation conditions of ocean conditions, according to the different mechanisms of CHF of DNB and Dryout, the CHF mechanism model of DNB and Dryout is established under the corresponding ocean conditions. The DNB type and Dryout type CHF prediction programs are compiled based on the model. Finally, the DNB type and Dryout type CHF prediction models are coupled with the natural circulation system program by judging the boiling crisis region. In the coupled model, the natural circulation system model is used as the main program framework, and the DNB type and Dryout boiling crisis mechanism models are used as subroutines in the main program framework. The ocean condition module provides boundary conditions for the natural circulation system and boiling crisis calculation. Natural circulation system provides external thermal hydraulic environment boundary for boiling crisis calculation. The boiling crisis feeds back the calculation results to the natural circulation system. The relationship between the models is shown in Figure8.Figure 8
Model coupling diagram.The flow of the CHF prediction coupling procedure under natural conditions of ocean conditions is as follows. Figure7 shows the flowchart of the calculation procedure.(1)
First, the natural circulation system program reads geometric parameters, control parameters, and initialization parameters.(2)
Set an initial heated section heat flux densityqm, start steady state calculation, and also adjust the preheating section, condenser parameters, and the resistance coefficient of the valve in the loop to reach the targeted system flow rate, the inlet temperature of heated section.(3)
After the steady state calculation is completed, the steady state result is used as the initial condition to start the transient calculation. The transient calculation obtains the mass quality of the outlet of the heated section for each time step, and the starting point mass quality of the annular flow under the current working condition. By comparing the two mass qualities, we can know whether there is an annular flow in the heated section.(4)
If there is no annular flow in the heated section, the DNB type boiling crisis mechanism model is called for calculation. If annular flow occurs in the heated section, the Dryout type boiling crisis is called to calculate.


## 4. Model Verification

To verify the Dryout type CHF model proposed in this paper, the CHF prediction correlation in circular tube proposed by Katto and Ohno [22] is employed. The correlation proposed by Katto is described as follows:(23)qCHF=0.234Ghfgρgρf0.513σρfG2L0.433L/De0.271+0.0031L/De.The CHF values predicted by our model and the correlation are compared in different equivalent diameters, channel lengths, system pressures, and mass flux. The detailed comparison between predicted results and correlation data is shown in Figure9. It can be seen that the predicted values show good agreement with the Katto correlation.Figure 9
Comparisons between predicted CHF and Katto correlation. (a) With respect to equivalent diameter. (b) With respect to channel length. (c) With respect to system pressure. (d) With respect to mass flux.
(a)(b)(c)(d)The verification of the DNB type CHF model can be found in reference [18].


## 5. Results and Discussion

### 5.1. Influence of Heeling on CHF

Figure10 is a schematic diagram of the CHF reduction under heeling condition. The flow characteristic curve of the natural circulation is determined by the loop structure. After the loop is heeled, the height difference between the hot and cold cores is changed, so the flow characteristic curve of the loop also changes. For the natural circulation loop studied in this paper, since the condenser and the heated section are not at a vertical height, the cold and hot core differences of the loop firstly increase with the positive slope of the loop and begin to decrease after reaching the maximum value, so the system flow rate firstly increases to the maximum and then decreases. As shown in Figure 10, the system mass flow rate reaches maximum value when the inclination angle is +5°, and the loop flow characteristic curve is the line corresponding to loop heeling +5°. The loop flow characteristic curve moves below the steady state line when the loop inclines −5°. For the heated section, the annular liquid film in the rectangular channel in the static state is evenly distributed. With the inclination of the circuit, the flow of the upper and lower liquid films in the channel is redistributed due to the action of gravity, and the upper liquid film thickness will be thinner than the lower liquid film thickness, and the boiling crisis will first occur in the upper pipe wall, so the boiling crisis trend line will move to the left. The difference in the flow rate of the upper and lower liquid film is related to the inclination angle (i.e., the component of gravity along the vertical flow direction) and the fluid flow velocity (i.e., inertial force). In this paper, the distribution of liquid film flow in the upper and lower sides of the channel at the beginning of the annular flow is calculated by the following formula [9].(24)Gf,t=Gf21−162+Fr2,where Gf,t is the mass flow velocity of the upper side of the channel at the beginning of the annular flow, the unit kg·m−2·s−1. Fr is Froude number, calculated by(25)Fr=GgsinθDeρgρf−ρg0.5.Figure 10
Influence of heeling.


#### 5.1.1. Influence of Heeling Angle

Firstly, the variation law of CHF is analyzed when the whole natural circulation loop is inclined. At this time, the heeling motion not only affects the thermal hydraulic parameters of the natural circulation system, but also affects the liquid film distribution in the rectangular channel. Figure11 shows the influence of the inclination angle on CHF. It can be seen from the figure that as the system inclines, the system hot and cold section position difference changes, the loop flow rate changes, and the inlet temperature of the heated section also changes. The inclination causes the distribution of the liquid film uneven. Under the combined action of the two factors, the heat required for the evaporation of the liquid film at the end of the annular flow is reduced. The research object in this paper, shown in Figure 1, is not asymmetric about the axis, so the deviation of system mass flow rate and CHF around zero degree is asymmetric. CHF decreases as the inclination angle increases from 0° to −45° or 45° and the maximum decrease of the CHF under heeling condition is 28.89% when the inclination angle is −45°.Figure 11
Influence of heeling angle on CHF.All in all, the heeling not only changes the thermal hydraulic parameters of the system, but also makes the liquid film thickness distribution in the channel of the heated section uneven. Under the combined influence of the two factors, the CHF of the natural circulation system changes.


### 5.2. Influence of Heaving Motion on CHF

Figure12 shows the influence of the heaving condition on the loop flow characteristic curve and the boiling crisis trend line. It can be seen from the figure that the flow rate of the system fluctuates above and below the static value under the conditions of heaving, and the left shift of the critical trend line of the boiling is mainly caused by the change of the pressure in the flow channel caused by the additional force. The pressure change will affect liquid deposition, entrainment, and evaporation of liquid film and result in a decrease in CHF. The intersection of the new boiling crisis trend line and the new loop flow characteristic curve line is the CHF under heaving condition. When boiling crisis occurs, the system flow rate is between the minimum value and the steady state value. The dotted black line is the loop flow characteristic curve of the system mass flow rate when boiling crisis occurs.Figure 12
Influence of heaving.


#### 5.2.1. Influence of Heaving Period

Figure13 shows the variation of the CHF value of the natural circulation system when the fixed heaving acceleration amplitude is 0.1 g and 0.3 g. It can be seen that after the amplitude of the heaving acceleration is fixed, when the heaving period is very small (less than 5 s), the CHF value increases with the increase of the period. After the period is greater than 10 s, the CHF value tends to be constant with the increase of the heaving period. At this time, the heaving motion not only affects the macroscopic thermal hydraulic parameters of the natural circulation system, but also affects the additional force caused by the heaving brings about the pressure change in the heated section, while the pressure affects the deposition, entrainment, and evaporation of the liquid film, and the liquid film will also be axial turbulence, and the CHF will change under the combined effect of the two factors.Figure 13
Influence of heaving period on CHF. (a)A = 0.1 g. (b) A = 0.3 g.
(a)(b)When the amplitude of the heaving acceleration is 0.1 g, the maximum increase of the CHF is 4.06% when the heaving period is changed. When the heaving acceleration amplitude is 0.3 g, the maximum increase of the CHF is 24.71% when the heaving period is changed. It shows that the greater the amplitude of the heaving acceleration, the more significant the impact of the heaving period.


#### 5.2.2. Influence of Heaving Acceleration Amplitude

Figure14 shows the variation of the CHF of the natural circulation system when the fixed heaving period is 3 s and 10 s. At this time, the heaving motion has an effect on the thermal hydraulic parameters of the natural circulation system and the liquid film distribution in the heated section, that is, the external thermal hydraulic environment at the critical point of boiling changes. It can be seen from the figure that, after the fixed heaving period, as the amplitude of the acceleration increases, the amplitude of the flow fluctuation of the natural circulation system increases, and the system mass flow rate at the critical point of boiling decreases, and the inlet temperature of the heated section increases, so the system CHF value also gradually decreases.Figure 14
Influence of heaving acceleration amplitude on CHF. (a)T = 3 s. (b) T = 10 s.
(a)(b)When the heaving period is 3 s, the maximum CHF decrease is 35.36% with the increase of the acceleration amplitude. When the heaving period is 10 s, the maximum CHF decrease is 27.58%. It indicates that the smaller the heaving period, the more significant the effect of the acceleration amplitude of the heaving.


### 5.3. Influence of Rolling Motion on CHF

Figure15 shows the influence of the rolling condition on the flow characteristic curve of the natural circulation system and the boiling crisis trend line. Since the rolling motion causes the heeling of the loop, and the centripetal drive head and the tangential drive head are simultaneously generated, the average flow rate of the system is reduced. The system flow fluctuates up and down around the average flow rate. The additional pressure drop generated by the rolling motion causes change of the droplet entrainment, deposition rate, and liquid film evaporation rate, and the inclination of the loop causes the distribution of the upper and lower liquid films in the channel to be uneven. Under the joint action of the two factors, the boiling crisis trend line under rolling condition moves to the left of the steady state. When boiling crisis occurs, the system flow rate is between the minimum value and the average value. The dotted black line is the loop flow characteristic curve of the system mass flow rate when boiling crisis occurs.Figure 15
Influence of rolling.


#### 5.3.1. Influence of Rolling Period

As shown in Figure16, when the natural circulation loop is rolling, the CHF of the heated section and the inlet temperature of the heated section and system mass flow rate at the critical boiling point are analyzed. It can be seen from the figure that the fixed rolling angular acceleration amplitude is 0.1 rad·s−2 and 0.6 rad·s−2. As the rolling period increases, the rolling angle and angular velocity increase, and the additional pressure drop of the system increases. The fluctuation amplitude of the thermal hydraulic parameters of the system increases. The mass flow rate of the system decreases, and the inlet temperature of the heated section increases when the boiling threshold occurs. The corresponding CHF decreases. When the amplitude of the rolling angle acceleration is 0.1 rad·s−2, the maximum decrease of the CHF is 44.9% with the increase of the rolling period. When the amplitude of the rolling angle acceleration is 0.6 rad·s−2, the maximum reduction in CHF is 13.5% with the increase of the rolling period.Figure 16
Influence of rolling period on CHF. (a)A = 0.1 rad/s2. (b) A = 0.6 rad/s2.
(a)(b)


#### 5.3.2. Influence of Rolling Angle Acceleration Amplitude

Figure17 shows the influence of the rolling angle acceleration on the CHF and the thermal parameters of the boiling crisis point system with a fixed rolling period of 4.5 s. It can be seen from the figure that as the amplitude of the acceleration of the rolling angle increases, the additional pressure drop caused by the rolling affects the thermal hydraulic parameters of the system. From the previous analysis, the fluctuation of the system flow rate increases with the increase of the rolling angle acceleration. The mass flow rate at the critical point of boiling decreases. The inlet temperature of the heated section increases. On the other hand, with the increase of acceleration of the rolling angle, the influence on the liquid film at the exit of the heated section increases. Under the combined influence of the two factors, the CHF reduces. When the rolling period is 4.5 s, the maximum decrease in CHF is 8.5% with the increase of the amplitude of the rolling angle acceleration.Figure 17
Influence of rolling angle acceleration amplitude on CHF,T = 4.5 s.


#### 5.3.3. Influence of Rolling Angle Amplitude

Figure18 shows the influence of the rolling angle amplitude on the system parameters of boiling crisis point and CHF. It can be seen that as the amplitude of the rolling angle increases, the additional pressure drop of the system increases, the amplitude of the system fluctuations such as the flow rate increases, and the influence on the liquid film is enhanced. As the rolling angle increases, the mass flow rate at the boiling crisis point decreases, the inlet temperature increases, and the CHF decreases. When the rolling period is 4.5 s, the maximum CHF decrease is 12.95% as the amplitude of the rolling angle increases.Figure 18
Influence of rolling angle amplitude on CHF,T = 4.5 s.


## 5.1. Influence of Heeling on CHF

Figure10 is a schematic diagram of the CHF reduction under heeling condition. The flow characteristic curve of the natural circulation is determined by the loop structure. After the loop is heeled, the height difference between the hot and cold cores is changed, so the flow characteristic curve of the loop also changes. For the natural circulation loop studied in this paper, since the condenser and the heated section are not at a vertical height, the cold and hot core differences of the loop firstly increase with the positive slope of the loop and begin to decrease after reaching the maximum value, so the system flow rate firstly increases to the maximum and then decreases. As shown in Figure 10, the system mass flow rate reaches maximum value when the inclination angle is +5°, and the loop flow characteristic curve is the line corresponding to loop heeling +5°. The loop flow characteristic curve moves below the steady state line when the loop inclines −5°. For the heated section, the annular liquid film in the rectangular channel in the static state is evenly distributed. With the inclination of the circuit, the flow of the upper and lower liquid films in the channel is redistributed due to the action of gravity, and the upper liquid film thickness will be thinner than the lower liquid film thickness, and the boiling crisis will first occur in the upper pipe wall, so the boiling crisis trend line will move to the left. The difference in the flow rate of the upper and lower liquid film is related to the inclination angle (i.e., the component of gravity along the vertical flow direction) and the fluid flow velocity (i.e., inertial force). In this paper, the distribution of liquid film flow in the upper and lower sides of the channel at the beginning of the annular flow is calculated by the following formula [9].(24)Gf,t=Gf21−162+Fr2,where Gf,t is the mass flow velocity of the upper side of the channel at the beginning of the annular flow, the unit kg·m−2·s−1. Fr is Froude number, calculated by(25)Fr=GgsinθDeρgρf−ρg0.5.Figure 10
Influence of heeling.


### 5.1.1. Influence of Heeling Angle

Firstly, the variation law of CHF is analyzed when the whole natural circulation loop is inclined. At this time, the heeling motion not only affects the thermal hydraulic parameters of the natural circulation system, but also affects the liquid film distribution in the rectangular channel. Figure11 shows the influence of the inclination angle on CHF. It can be seen from the figure that as the system inclines, the system hot and cold section position difference changes, the loop flow rate changes, and the inlet temperature of the heated section also changes. The inclination causes the distribution of the liquid film uneven. Under the combined action of the two factors, the heat required for the evaporation of the liquid film at the end of the annular flow is reduced. The research object in this paper, shown in Figure 1, is not asymmetric about the axis, so the deviation of system mass flow rate and CHF around zero degree is asymmetric. CHF decreases as the inclination angle increases from 0° to −45° or 45° and the maximum decrease of the CHF under heeling condition is 28.89% when the inclination angle is −45°.Figure 11
Influence of heeling angle on CHF.All in all, the heeling not only changes the thermal hydraulic parameters of the system, but also makes the liquid film thickness distribution in the channel of the heated section uneven. Under the combined influence of the two factors, the CHF of the natural circulation system changes.


## 5.1.1. Influence of Heeling Angle

Firstly, the variation law of CHF is analyzed when the whole natural circulation loop is inclined. At this time, the heeling motion not only affects the thermal hydraulic parameters of the natural circulation system, but also affects the liquid film distribution in the rectangular channel. Figure11 shows the influence of the inclination angle on CHF. It can be seen from the figure that as the system inclines, the system hot and cold section position difference changes, the loop flow rate changes, and the inlet temperature of the heated section also changes. The inclination causes the distribution of the liquid film uneven. Under the combined action of the two factors, the heat required for the evaporation of the liquid film at the end of the annular flow is reduced. The research object in this paper, shown in Figure 1, is not asymmetric about the axis, so the deviation of system mass flow rate and CHF around zero degree is asymmetric. CHF decreases as the inclination angle increases from 0° to −45° or 45° and the maximum decrease of the CHF under heeling condition is 28.89% when the inclination angle is −45°.Figure 11
Influence of heeling angle on CHF.All in all, the heeling not only changes the thermal hydraulic parameters of the system, but also makes the liquid film thickness distribution in the channel of the heated section uneven. Under the combined influence of the two factors, the CHF of the natural circulation system changes.


## 5.2. Influence of Heaving Motion on CHF

Figure12 shows the influence of the heaving condition on the loop flow characteristic curve and the boiling crisis trend line. It can be seen from the figure that the flow rate of the system fluctuates above and below the static value under the conditions of heaving, and the left shift of the critical trend line of the boiling is mainly caused by the change of the pressure in the flow channel caused by the additional force. The pressure change will affect liquid deposition, entrainment, and evaporation of liquid film and result in a decrease in CHF. The intersection of the new boiling crisis trend line and the new loop flow characteristic curve line is the CHF under heaving condition. When boiling crisis occurs, the system flow rate is between the minimum value and the steady state value. The dotted black line is the loop flow characteristic curve of the system mass flow rate when boiling crisis occurs.Figure 12
Influence of heaving.


### 5.2.1. Influence of Heaving Period

Figure13 shows the variation of the CHF value of the natural circulation system when the fixed heaving acceleration amplitude is 0.1 g and 0.3 g. It can be seen that after the amplitude of the heaving acceleration is fixed, when the heaving period is very small (less than 5 s), the CHF value increases with the increase of the period. After the period is greater than 10 s, the CHF value tends to be constant with the increase of the heaving period. At this time, the heaving motion not only affects the macroscopic thermal hydraulic parameters of the natural circulation system, but also affects the additional force caused by the heaving brings about the pressure change in the heated section, while the pressure affects the deposition, entrainment, and evaporation of the liquid film, and the liquid film will also be axial turbulence, and the CHF will change under the combined effect of the two factors.Figure 13
Influence of heaving period on CHF. (a)A = 0.1 g. (b) A = 0.3 g.
(a)(b)When the amplitude of the heaving acceleration is 0.1 g, the maximum increase of the CHF is 4.06% when the heaving period is changed. When the heaving acceleration amplitude is 0.3 g, the maximum increase of the CHF is 24.71% when the heaving period is changed. It shows that the greater the amplitude of the heaving acceleration, the more significant the impact of the heaving period.


### 5.2.2. Influence of Heaving Acceleration Amplitude

Figure14 shows the variation of the CHF of the natural circulation system when the fixed heaving period is 3 s and 10 s. At this time, the heaving motion has an effect on the thermal hydraulic parameters of the natural circulation system and the liquid film distribution in the heated section, that is, the external thermal hydraulic environment at the critical point of boiling changes. It can be seen from the figure that, after the fixed heaving period, as the amplitude of the acceleration increases, the amplitude of the flow fluctuation of the natural circulation system increases, and the system mass flow rate at the critical point of boiling decreases, and the inlet temperature of the heated section increases, so the system CHF value also gradually decreases.Figure 14
Influence of heaving acceleration amplitude on CHF. (a)T = 3 s. (b) T = 10 s.
(a)(b)When the heaving period is 3 s, the maximum CHF decrease is 35.36% with the increase of the acceleration amplitude. When the heaving period is 10 s, the maximum CHF decrease is 27.58%. It indicates that the smaller the heaving period, the more significant the effect of the acceleration amplitude of the heaving.


## 5.2.1. Influence of Heaving Period

Figure13 shows the variation of the CHF value of the natural circulation system when the fixed heaving acceleration amplitude is 0.1 g and 0.3 g. It can be seen that after the amplitude of the heaving acceleration is fixed, when the heaving period is very small (less than 5 s), the CHF value increases with the increase of the period. After the period is greater than 10 s, the CHF value tends to be constant with the increase of the heaving period. At this time, the heaving motion not only affects the macroscopic thermal hydraulic parameters of the natural circulation system, but also affects the additional force caused by the heaving brings about the pressure change in the heated section, while the pressure affects the deposition, entrainment, and evaporation of the liquid film, and the liquid film will also be axial turbulence, and the CHF will change under the combined effect of the two factors.Figure 13
Influence of heaving period on CHF. (a)A = 0.1 g. (b) A = 0.3 g.
(a)(b)When the amplitude of the heaving acceleration is 0.1 g, the maximum increase of the CHF is 4.06% when the heaving period is changed. When the heaving acceleration amplitude is 0.3 g, the maximum increase of the CHF is 24.71% when the heaving period is changed. It shows that the greater the amplitude of the heaving acceleration, the more significant the impact of the heaving period.


## 5.2.2. Influence of Heaving Acceleration Amplitude

Figure14 shows the variation of the CHF of the natural circulation system when the fixed heaving period is 3 s and 10 s. At this time, the heaving motion has an effect on the thermal hydraulic parameters of the natural circulation system and the liquid film distribution in the heated section, that is, the external thermal hydraulic environment at the critical point of boiling changes. It can be seen from the figure that, after the fixed heaving period, as the amplitude of the acceleration increases, the amplitude of the flow fluctuation of the natural circulation system increases, and the system mass flow rate at the critical point of boiling decreases, and the inlet temperature of the heated section increases, so the system CHF value also gradually decreases.Figure 14
Influence of heaving acceleration amplitude on CHF. (a)T = 3 s. (b) T = 10 s.
(a)(b)When the heaving period is 3 s, the maximum CHF decrease is 35.36% with the increase of the acceleration amplitude. When the heaving period is 10 s, the maximum CHF decrease is 27.58%. It indicates that the smaller the heaving period, the more significant the effect of the acceleration amplitude of the heaving.


## 5.3. Influence of Rolling Motion on CHF

Figure15 shows the influence of the rolling condition on the flow characteristic curve of the natural circulation system and the boiling crisis trend line. Since the rolling motion causes the heeling of the loop, and the centripetal drive head and the tangential drive head are simultaneously generated, the average flow rate of the system is reduced. The system flow fluctuates up and down around the average flow rate. The additional pressure drop generated by the rolling motion causes change of the droplet entrainment, deposition rate, and liquid film evaporation rate, and the inclination of the loop causes the distribution of the upper and lower liquid films in the channel to be uneven. Under the joint action of the two factors, the boiling crisis trend line under rolling condition moves to the left of the steady state. When boiling crisis occurs, the system flow rate is between the minimum value and the average value. The dotted black line is the loop flow characteristic curve of the system mass flow rate when boiling crisis occurs.Figure 15
Influence of rolling.


### 5.3.1. Influence of Rolling Period

As shown in Figure16, when the natural circulation loop is rolling, the CHF of the heated section and the inlet temperature of the heated section and system mass flow rate at the critical boiling point are analyzed. It can be seen from the figure that the fixed rolling angular acceleration amplitude is 0.1 rad·s−2 and 0.6 rad·s−2. As the rolling period increases, the rolling angle and angular velocity increase, and the additional pressure drop of the system increases. The fluctuation amplitude of the thermal hydraulic parameters of the system increases. The mass flow rate of the system decreases, and the inlet temperature of the heated section increases when the boiling threshold occurs. The corresponding CHF decreases. When the amplitude of the rolling angle acceleration is 0.1 rad·s−2, the maximum decrease of the CHF is 44.9% with the increase of the rolling period. When the amplitude of the rolling angle acceleration is 0.6 rad·s−2, the maximum reduction in CHF is 13.5% with the increase of the rolling period.Figure 16
Influence of rolling period on CHF. (a)A = 0.1 rad/s2. (b) A = 0.6 rad/s2.
(a)(b)


### 5.3.2. Influence of Rolling Angle Acceleration Amplitude

Figure17 shows the influence of the rolling angle acceleration on the CHF and the thermal parameters of the boiling crisis point system with a fixed rolling period of 4.5 s. It can be seen from the figure that as the amplitude of the acceleration of the rolling angle increases, the additional pressure drop caused by the rolling affects the thermal hydraulic parameters of the system. From the previous analysis, the fluctuation of the system flow rate increases with the increase of the rolling angle acceleration. The mass flow rate at the critical point of boiling decreases. The inlet temperature of the heated section increases. On the other hand, with the increase of acceleration of the rolling angle, the influence on the liquid film at the exit of the heated section increases. Under the combined influence of the two factors, the CHF reduces. When the rolling period is 4.5 s, the maximum decrease in CHF is 8.5% with the increase of the amplitude of the rolling angle acceleration.Figure 17
Influence of rolling angle acceleration amplitude on CHF,T = 4.5 s.


### 5.3.3. Influence of Rolling Angle Amplitude

Figure18 shows the influence of the rolling angle amplitude on the system parameters of boiling crisis point and CHF. It can be seen that as the amplitude of the rolling angle increases, the additional pressure drop of the system increases, the amplitude of the system fluctuations such as the flow rate increases, and the influence on the liquid film is enhanced. As the rolling angle increases, the mass flow rate at the boiling crisis point decreases, the inlet temperature increases, and the CHF decreases. When the rolling period is 4.5 s, the maximum CHF decrease is 12.95% as the amplitude of the rolling angle increases.Figure 18
Influence of rolling angle amplitude on CHF,T = 4.5 s.


## 5.3.1. Influence of Rolling Period

As shown in Figure16, when the natural circulation loop is rolling, the CHF of the heated section and the inlet temperature of the heated section and system mass flow rate at the critical boiling point are analyzed. It can be seen from the figure that the fixed rolling angular acceleration amplitude is 0.1 rad·s−2 and 0.6 rad·s−2. As the rolling period increases, the rolling angle and angular velocity increase, and the additional pressure drop of the system increases. The fluctuation amplitude of the thermal hydraulic parameters of the system increases. The mass flow rate of the system decreases, and the inlet temperature of the heated section increases when the boiling threshold occurs. The corresponding CHF decreases. When the amplitude of the rolling angle acceleration is 0.1 rad·s−2, the maximum decrease of the CHF is 44.9% with the increase of the rolling period. When the amplitude of the rolling angle acceleration is 0.6 rad·s−2, the maximum reduction in CHF is 13.5% with the increase of the rolling period.Figure 16
Influence of rolling period on CHF. (a)A = 0.1 rad/s2. (b) A = 0.6 rad/s2.
(a)(b)


## 5.3.2. Influence of Rolling Angle Acceleration Amplitude

Figure17 shows the influence of the rolling angle acceleration on the CHF and the thermal parameters of the boiling crisis point system with a fixed rolling period of 4.5 s. It can be seen from the figure that as the amplitude of the acceleration of the rolling angle increases, the additional pressure drop caused by the rolling affects the thermal hydraulic parameters of the system. From the previous analysis, the fluctuation of the system flow rate increases with the increase of the rolling angle acceleration. The mass flow rate at the critical point of boiling decreases. The inlet temperature of the heated section increases. On the other hand, with the increase of acceleration of the rolling angle, the influence on the liquid film at the exit of the heated section increases. Under the combined influence of the two factors, the CHF reduces. When the rolling period is 4.5 s, the maximum decrease in CHF is 8.5% with the increase of the amplitude of the rolling angle acceleration.Figure 17
Influence of rolling angle acceleration amplitude on CHF,T = 4.5 s.


## 5.3.3. Influence of Rolling Angle Amplitude

Figure18 shows the influence of the rolling angle amplitude on the system parameters of boiling crisis point and CHF. It can be seen that as the amplitude of the rolling angle increases, the additional pressure drop of the system increases, the amplitude of the system fluctuations such as the flow rate increases, and the influence on the liquid film is enhanced. As the rolling angle increases, the mass flow rate at the boiling crisis point decreases, the inlet temperature increases, and the CHF decreases. When the rolling period is 4.5 s, the maximum CHF decrease is 12.95% as the amplitude of the rolling angle increases.Figure 18
Influence of rolling angle amplitude on CHF,T = 4.5 s.


## 6. Conclusion

In this paper, a series of theoretical models including CHF mechanistic models and natural circulation system model have been developed, and the CHF in the rectangular channel of natural circulation system under ocean conditions is studied. The prediction model of DNB and Dryout boiling crisis under ocean conditions, the thermal hydraulic characteristics analysis model of natural circulation system, and the theoretical model of additional force of ocean conditions are established. The corresponding calculation program was developed through multimodel and cross-dimensional coupling, and the CHF characteristics of natural circulation under ocean conditions were studied. Within the scope of the study parameters, the influence of ocean conditions on the CHF is great, and the maximum degradation of CHF is about 45%.



---

*Source: 9390645-2020-09-02.xml*