Numerical and experimental analysis of cavitation characteristics in safety valves of the nuclear power second circuit using a modified cavitation model

Abstract

Cavitation frequently arises in the safety valve of nuclear power plants’ secondary circuits operating under high pressure conditions. This study integrates valve flow characteristics and velocity strain rate corrections into the Zwart-Gerber-Belamri model to accurately simulate cavitation inside the valve, reducing the impact of physical empirical coefficient variations on cavitation length prediction. Subsequently, a visualisation test rig is developed to validate the accuracy of the numerical model, and experimental cavitation results are obtained using the grayscale detection method. The evaporation/condensation coefficients are optimised using the AES-MSI model and GA based on the experimental results. The accuracy of the constructed model is validated by comparing it with experimental results obtained under various operating conditions. Finally, the high-fidelity numerical model is employed to investigate the effects of pressure drop and valve openings on cavitation, elucidating the underlying mechanisms governing cavitation variations resulting from pressure drops. Furthermore, a comprehensive equation is derived to determine the effective flow area, aiding in the identification of cavitation locations and offering insights into the relationship between cavitation behaviour and valve openings. The modified cavitation model proposed in this study can be readily extended to investigate cavitation prediction in other valves or throttle elements.

KEYWORDS:

• Safety valve
• cavitation model
• nuclear power plant
• CFD
• visualisation experiments

1. Introduction

Nuclear power is a clean, efficient, and reliable form of energy that harnesses nuclear reactions to generate electricity (Bayramov, 2023). Safety valves play a crucial role in ensuring safe operation of nuclear power plants by serving as the last line of defense against dangerous pressure buildup in the steam, cooling water, and compressed air systems (Bjerre et al., 2017). Among these valves, the safety valve used in the second circuit is specifically designed to prevent excessive pressure generated by the water pump supplying water to the steam generator, as illustrated in Figure 1 (Li et al., 2023b). Under high-pressure conditions, the valve opens and reduces the flow-through area, causing a decrease in pressure near the valve disc. When the local pressure falls below the saturated vapour pressure of the fluid medium, cavitation can occur, which is a phase change process that transforms liquid into gas (Nozawa, 2023). Cavitation bubbles then rapidly expand and release a huge pressure pulse when they collapse, causing material failure in the cavitation zone (Wang et al., 2021a). Cavitation damage inside the valve has been identified as a major cause of valve failure (Sreedhar et al., 2017). Therefore, studying cavitation is essential for safe engineering applications in flow systems, especially for the safe operation of nuclear power systems.

Figure 1. Safety valve in second loop of nuclear power plant. The red dotted line indicates the designated installation position of the safety valve within the secondary circuit of a nuclear power plant. In the event that the pump pressure exceeds the system's maximum permissible pressure, the safety valve promptly opens, releasing the medium to restore the system pressure to its normal level.

The study of cavitation can be traced back to 1753 (Sun et al., 2022). When Euler discovered that the pressure in a water pipe could decrease to a negative value, causing the water to separate from the pipe wall and create a cavity through vaporisation of the liquid. This vacuum should be avoided to prevent damage to the system. In the late nineteenth century, the shipbuilding industry encountered problems when propeller speed was increased continuously, leading to a decrease in a ship's speed driven by a steam turbine. This was due to the impact of cavitation, which led scholars to start recognising its effects, leading to increased attention towards the study of cavitation. In the available literature, the earliest investigation of valve cavitation dates back to the work of MacLellan et al. in 1960 (MacLellan et al., 1960). In their study, they detected the cavitation phenomenon during fluid flow through the spool valve orifice. Following this, numerous scholars focused on the study of valve cavitation. Advanced testing technology and improved computer performance have led to experimental testing and numerical simulation becoming the primary methods for studying the phenomenon in valves. Due to the intricate structure of the valve flow channel and the highly nonlinear nature of cavitation in fluid dynamics, it is impractical to rely solely on theoretical models to derive mathematical equations for the cavitation phenomenon in valves.

Studies on valve cavitation using experimental techniques are divided into two main categories. The first category involves the indirect study of cavitation through testing the mass flow rate of the valve, as cavitation significantly impacts valve throttling. For instance, Amirante et al. (2014) conducted an experimental analysis to investigate the effect of cavitation on the performance of a proportional direct-acting hydraulic reversing valve. They used a constructed experimental rig and a high-precision gear flow metre to measure the mass flow rate under cavitation conditions. Valdés et al. (2014) investigated the cavitation flow in ball check valves by measuring the mass flow to control various degrees of cavitation through pressure adjustments at the valve outlet. The characteristic flow coefficients of the valves were subsequently analyzed. However, while flow measurement provides some insight into cavitation, it does not constitute a comprehensive study of valve cavitation. The cavitation phenomenon should be more critically evaluated by considering the cavitation length, as the morphology of the cavitation bubble is directly related to the degree of corrosion experienced by the valve's key components. Consequently, visualisation experiments have become an essential approach to studying cavitation. Oshima and Ichikawa (1985) observed the cavitation phenomenon of converging and diverging flows in the poppet valve using the experimental half-cut test model. They investigated the effects of cavitation on the mass flow rate, noise, pressure distribution, and cavitation onset boundary of the water hydraulic poppet valve.

Kumagai et al. (2016) developed a visualisation system to analyze the mechanism of cavitation-induced spool vibration. They used a high-speed camera to capture cavitation inside a transparent cone valve (with the body made of PMMA). Li et al. (2013) developed a visualisation experimental system to observe the cavitation phenomenon between the flapper and nozzle in an electro-hydraulic servo valve. They determined the location of the cavitation source and explored the cavitation evolution with respect to Reynolds number. Wang et al. (2018) utilised a camera and a transparent throttle valve made of PMMA material to capture the inner bubble morphology of the valve. They employed convolutional neural networks for rapid detection of spool features and developed intelligent algorithms based on two-dimensional detection and three-dimensional reconstruction of the bubbles. Lu et al. (2022) employed 46# anti-wear hydraulic oil as the working medium, and the valve body of the spool valve was fabricated using highly transparent Plexiglas material. They captured cavitation images inside the valve using a high-speed camera and adjusted the valve opening by manipulating bolts to obtain cavitation images at different openings. Yi et al. (2015) observed the cavitation phenomenon inside the relief valve by using a camera and a transparent valve body. They experimentally investigated the interaction between the spool vibration characteristics and the cavitation characteristics.

Capturing cavitation cavities through high-speed cameras is a straightforward approach for conducting cavitation experiments. However, this method requires a transparent valve body, which limits the range of temperature and pressure of the medium during experimental tests. Moreover, the extensive time and economic costs associated with this method hinder the optimal valve design. As computer performance has improved, numerical simulation methods based on computational fluid dynamics (CFD) have become crucial tools for cavitation research. For example, Zong et al. (2022) conducted a numerical simulation of heat-fluid-solid coupling for explosion-proof valves with complex flow paths using the CFD method. They proposed a reduced-dimensional computational fluid dynamics modelling method to simulate the barrier effect of the flame retardant sheet. Steady and transient simulations were conducted, leading to the conclusion that the computational accuracy of the equivalent volumetric method proposed in their study was higher than that of the equivalent porous method. Within the computational framework based on the theory of CFD, two models are commonly adopted for gas–liquid two-phase flow: the two-fluid model and the homogeneous equilibrium flow model (HEM) (Tahmasebi et al., 2017). HEM has gained widespread use in practical engineering due to its superior convenience and relatively mature method in practical applications (Salvador et al., 2017; Webber, 2011). Under the assumption of HEM, two mainstream cavitation models, the barotropic equation model (BEM), and the transport equation model (TEM), have been formed based on the determination of different density fields (Chen et al., 2016). However, researchers found that the BEM has certain limitations in predicting the convection and transport of cavitation as it does not reflect the phase change process in the physical model. As the understanding of cavitation improved, researchers turned to the TEM to capture the phase change process in cavitation with greater accuracy. The most representative models based on the assumption of TEM are the Schnerr-Sauer cavitation model (Schnerr-Sauer model), Singhal fully cavitation model (Singhal model), and Zwart-Gerber-Belamri cavitation model (Z-G-B model) (Schnerr & Sauer, 2001; Singhal et al., 2002; Zwart et al., 2004). Liu et al. (2019) developed a CFD cavitation model for the regulating valve using the Schnerr-Sauer model and the standard k-ε turbulence model. They analyzed the changes in the location, shape, area, and density of cavitation in both the axial and radial directions at different inlet pressure conditions. Yaghoubi et al. (2018) performed a numerical simulation of a globe control valve to examine the effects of orifice diameter and various trim types on cavitation delay performance. They found that using three trims resulted in the best economic and technical performance in terms of delaying cavitation. In addition, some scholars have also adopted new cavitation models or methods for the simulation of cavitation flow. For example, Cheng et al. (2021) propose a novel Eulerian-Lagrangian cavitation model that considers the influence of non-condensable gases. The global flow field is solved using the Eulerian method, while the migration of non-condensable bubbles towards the vortex core is tracked employing the Lagrangian method. The mass source term in the original Schnerr-Sauer cavitation model is adjusted through modelling the relationship between the local gas concentration and its impact on cavitation. The model has been successfully applied to investigate the cavitation on NACA-16020 hydrofoils, and its effectiveness has been verified by comparing it with experimental results. Wang et al. (2021b) employed a multiscale Eulerian-Lagrangian method to simulate cavitation flow around a Clark-Y hydrofoil. The method combines large vortex simulation with the fluid volume method and tracks the microscale Lagrangian bubbles by solving the compressible Rayleigh-Plesset equation and the equations of motion for bubbles. A two-way coupled simulation of the interaction between the flow field and bubbles is performed using Gaussian kernel functions. The validity of the method is finally confirmed through experimental verification.

In the study of cavitation length and mechanisms of cavitation evolution, Wang et al. (2023) proposed a novel multiscale modelling approach to investigate numerically the cavitation formation mechanism of hydrofoils. In their study, a mixture model and a Lagrangian model were used to capture the macroscopic vapour structure and the behaviour of microbubbles, respectively. The study provides a comprehensive explanation of the cloud cavitation and its evolution process, as well as the changing law of bubble dynamics. Gao et al. (2021) conducted numerical simulations of cage control valves using the Z-G-B cavitation model to compare the flow rate and cavitation strength of different valve cage structures. The comparison shows that the improved cage structure has the smallest cavitation strength, and the paper also deduces that the cavitation strength and flow coefficient gradually decrease with the increase in the number of cages and the height of grooves. Yuan et al. (2019) conducted a numerical simulation of cavitation caused by a jet inside the poppet valve using the Schnerr-Sauer model. They analyzed the inhibitory effect of cavitation growth on the vortex and explained the interaction mechanism between cavitation and the vortex. Jin et al. (2020) explored the effect of four typical valve core shapes on cavitation in valves using the Schnerr-Sauer model. They aimed to reduce cavitation in regulating valves and concluded that the cavitation strength was higher with ellipsoidal and cylindrical valve cores. Additionally, they studied the relationship between valve opening and vapour volume, finding that the vapour volume tends to increase and then decrease with the increase in valve core displacement. Liu et al. (2021) investigated the characteristics of transient cavitation and pressure fluctuation inside a regulating valve using the Schnerr-Sauer model. They analyzed the effect of variations in the length-to-radius ratio of the valve spool on cavitation characteristics. Furthermore, other scholars have conducted numerous studies on valve cavitation using numerical simulation methods (Habibnejad et al., 2022; He et al., 2023; Qian et al., 2023).

From the literature summary above, it is evident that there is a wealth of literature on valve cavitation from experimental and numerical perspectives. However, there are still several shortcomings.

1. Currently, research on valve cavitation primarily focuses on regulating/control valves, with less attention paid to safety valves. The few studies on safety valves have primarily focused on two-dimensional (2-D) or simplified three-dimensional (3-D) modelling of the flow core region only, instead of utilising full-size three-dimensional modelling of the valve. These modelling approaches may result in inaccurate capture of cavitation compared to full-size 3-D modelling, as geometric boundaries can significantly impact fluid flow characteristics, particularly pressure drop, which is relevant to simulating cavitation.

2. Most researchers use cavitation models developed for hydrofoils and inducers, which may not be appropriate for simulating safety valve cavitation due to specific values of empirical coefficients. Adjustments to these coefficients are often necessary to account for specific operating conditions. For instance, Zwart et al. (2004) recommended setting the evaporation and condensation coefficients to 50 and 0.01, respectively, while in another study he set the evaporation and condensation coefficients to 0.4 and 0.001, respectively. Therefore, it is necessary to propose a robust method capable of eliminating the effect of the variation of these empirical coefficients on the results or a method capable of efficiently determining these coefficients.

3. Moreover, although many researchers have provided general guidelines for the influence of certain parameters on valve cavitation, the definition of cavitation length on safety valve remains unclear, and the explanation for the cavitation evolution mechanism has not been widely accepted. Therefore, further research is required to address these issues.

Based on the analysis presented above, we have modified the mass transfer term in the Z-G-B model by incorporating the flow characteristics of the valve and the velocity strain rate. This modification effectively eliminates the effects of the bubble diameter and nucleation site volume fraction on the mass transfer during the phase change process. Subsequently, we combined an adaptive ensemble of surrogate models with the genetic algorithm to determine the optimal evaporation and condensation coefficients, based on the cavitation length observed in the visual experiment. Finally, the validated high-fidelity cavitation model is employed to explain the evolution mechanism of cavitation under different pressure drops and valve openings.

The paper is organised as follows: Section 2 presents the mathematical models used in the modified cavitation model and the adaptive ensemble of surrogate models. Section 3 provides a detailed account of the construction process of CFD cavitation models, including simplification of the fluid computational domain, verification of mesh independence, and setting boundary conditions. Additionally, we discuss the validation of the physical dimensionality reduction of the modified cavitation model and the optimisation of evaporation/condensation coefficients based on visualisation experiments. Section 4 explores the impact of essential parameters, such as pressure drop between the inlet and outlet and valve opening, on cavitation in safety valves and provides insightful explanations from a physical standpoint to elucidate this relationship. Finally, in Section 5, we summarise our findings and suggest future research considerations.

2. Mathematical models

In this section, a detailed description of the conservation equation for multiphase flow and the turbulence equation is provided to establish the theoretical foundation for simulating two-phase flow in the relief valve. Subsequently, the theory behind the Z-G-B model is presented, which is based on the Rayleigh-Plesset equation and describes mass transfer between the liquid and vapour phases. Furthermore, the modified cavitation model is explained in depth, incorporating the valve's internal flow characteristics and the velocity strain rate to account for the evaporation term. Finally, the construction concept and relevant theory of the adaptive ensemble of surrogate models are presented, which facilitates the establishment of the mapping relationship between evaporation/condensation coefficients and cavitation lengths.

2.1. Governing equations of multiphase flow

In multiphase flow, the conservation equations are fundamental in describing the behaviour of each phase. These equations encompass the continuity, momentum, and energy equations, and their interdependence results in coupling. Solving the conservation equations provides crucial information on the velocity field, phase distribution, and thermodynamic properties of each phase. Therefore, achieving accurate modelling and numerical simulation of multiphase flow is contingent on the proper formulation and thorough comprehension of the conservation equations. In this study, the hydraulic fluid constitutes the primary phase, while the secondary phase comprises the vapour formed through the vaporisation of the hydraulic fluid. As there is minimal temperature change during the phase transition process, the energy equation is omitted. The continuity equation and the momentum equation are presented below. (1) $\begin{array}{rl}& \frac{\mathrm{\partial }{\rho }_m}{\mathrm{\partial }t}+\mathrm{\nabla }\cdot \left({\rho }_m{\overrightarrow{v}}_m\right)=0\end{array}$(1) (2) $\begin{array}{rl}& \frac{\mathrm{\partial }\left({\rho }_m{\overrightarrow{v}}_m\right)}{\mathrm{\partial }t}+\mathrm{\nabla }\cdot \left({\rho }_m{\overrightarrow{v}}_m{\overrightarrow{v}}_m\right)\\ & =-\mathrm{\nabla }P+\mathrm{\nabla }\cdot \left[\right({\mu }_t+{\mu }_m\left)\right(\mathrm{\nabla }{\overrightarrow{v}}_m+\mathrm{\nabla }{{\overrightarrow{v}}_m}^T\left)\right]\end{array}$(2) Here t is time; ${\rho }_m$ is the mixture density; ${\overrightarrow{v}}_m$ is the mass-averaged velocity; P is the pressure; ${\mu }_m$ is the mixture viscosity; and ${\mu }_t$ represents the turbulence viscosity, modelled using the SST $k-\omega$ model in our study. The mixture density ${\rho }_m$, viscosity ${\mu }_m$, and ${\overrightarrow{v}}_m$ are calculated as follows: (3) $\begin{array}{rl}{\rho }_m& =\sum _{i=1}^N{\alpha }_i{\rho }_i\end{array}$(3) (4) $\begin{array}{rl}{\mu }_m& =\sum _{i=1}^N{\alpha }_i{\mu }_i\end{array}$(4) (5) $\begin{array}{rl}{\overrightarrow{v}}_m& =\frac{\sum _{i=1}^N{\alpha }_i{\rho }_i{\overrightarrow{v}}_i}{\sum _{i=1}^N{\alpha }_i{\rho }_i}\end{array}$(5)

Here ${\alpha }_i$, ${\rho }_i$, ${\mu }_i$ and ${\overrightarrow{v}}_i$ are the volume fraction, density, viscosity, and velocity of the $i^{\mathrm{t}\mathrm{h}}$ phase, respectively.

The turbulence viscosity ${\mu }_t$ and the turbulent kinetic energy k and specific dissipation ω are presented as follows (Mustafa et al., 2023; Viitanen et al., 2022). (6) $\begin{array}{rl}{\mu }_t& =\frac{{\rho }_m{a}_{1}k}{max\left(a_1\omega ,S{F}_2\right)}\end{array}$(6) (7) $\begin{array}{rl}\frac{\mathrm{\partial }\left({\rho }_{m}k\right)}{\mathrm{\partial }t}+\mathrm{\nabla }\cdot \left({\rho }_{m}kU\right)& ={\tilde{P}}_k-{\text{β}}^{\star }{\rho }_m\omega k\\ & +\mathrm{\nabla }\cdot \left[\left(\mu +{\sigma }_k{\mu }_t\right)\mathrm{\nabla }k\right]\end{array}$(7) (8) $\begin{array}{rl}\frac{\mathrm{\partial }\left({\rho }_m\omega \right)}{\mathrm{\partial }t}+\mathrm{\nabla }\cdot \left({\rho }_m\omega U\right)& =\gamma {\rho }_m{S}^2-\text{β}{\rho }_m{\omega }^2\\ & +\mathrm{\nabla }\cdot \left[\left(\mu +{\sigma }_{\omega }{\mu }_t\right)\mathrm{\nabla }\omega \right]\\ & +\left(1-F_1\right)C{D}_{k\omega }\end{array}$(8)

Here where γ, σ_k, σ_ω, β, β^⋆ are model constants. ${\tilde{P}}_k$ represents the productions of turbulence kinetic energy, which can be written as follows: (9) $\{\begin{array}{l}{\tilde{P}}_k=min\left(P_k,10{\text{β}}^{\star }{\rho }_m\omega k\right)\\ P_k={\mu }_t{S}^2\end{array}$(9) In Equation (6), the value of parameter a₁ is 0.31. S represents the strain rate, which can be write as $S=\sqrt{2S_{ij}{S}_{ij}}$, and the $S_{ij}$ can be wirte as $S_{ij}=\frac{1}{2}\left(\frac{\mathrm{\partial }{u}_i}{\mathrm{\partial }{x}_j}+\frac{\mathrm{\partial }{u}_j}{\mathrm{\partial }{x}_i}\right)$. The remaining critical parameters in Equations (6)–(8) can be expressed as follows: (10) $\begin{array}{rl}{F}_1& =\tanh (\{min[max\left(\frac{\sqrt{}k}{{\text{β}}^{\star }\omega \mathrm{d}},\frac{500\nu }{{\mathrm{d}}^2\omega }\right),\\ & {\frac{4{\rho }_{\mathrm{m}}{\sigma }_{\omega 2}\mathrm{k}}{C{D}_{kw}{\mathrm{d}}^2}]\}}^4)\end{array}$(10) (11) $\begin{array}{rl}{F}_2& =\tanh \left({\left\{max\left(2\frac{\sqrt{}k}{{\text{β}}^{\star }\omega \mathrm{d}},\frac{500\nu }{{\mathrm{d}}^2\omega }\right)\right\}}^2\right)\end{array}$(11) (12) $\begin{array}{rl}C{D}_{kw}& =max\left(2{\rho }_m{\sigma }_{\omega 2}\frac{1}{\omega }\mathrm{\nabla }k\mathrm{\nabla }\omega ,{10}^{-10}\right)\end{array}$(12)

2.2. Theory of the cavitation model based on flow scale and velocity strain rate correction

The Z-G-B model, which is a frequently utilised cavitation model in commercial software for cavitation prediction, was originally introduced in 2000 by Zwart, Gerber, and Belamri (Jena et al., 2023). The firstly purpose of this model is to accurately capture cavitation resulting from the operation of a two-dimensional hydrofoil and subsequently also has been verified in other fields. In this study, An enhancement Z-G-B model is proposed by incorporating the flow scale associated with valve opening and velocity strain rate. The objective of this enhancement is to mitigate the impact of physical empirical parameters such as bubble diameter and nucleation volume fraction on the cavitation prediction. The modification of the Z-G-B model provides a more accurate prediction of cavitation inception, particularly in systems with varying flow scales. The Z-G-B model has been developed based on the Rayleigh-Plesset equation. Equation (13) presents the expression of the Rayleigh-Plesset equation, which explains the expansion and contraction of a single bubble due to the difference between the internal and external pressure (Geng et al., 2020; Snyder et al., 2016). This equation serves as the foundation for the development of the Z-G-B model. (13) $R_B\frac{d^2{R}_B}{d{t}^2}+\frac{3}{2}{\left(\frac{d{R}_B}{dt}\right)}^2+\frac{2\sigma }{R_B}=\frac{P_v-P}{{\rho }_l}$(13) Here, $R_B$ is the bubble radius, $\sigma$ is the surface tension coefficient, and $P_v$ is the vapour pressure.

By neglecting the second-order terms and the effect of surface tension, A simplified expression that shows the change rate of the bubble radius over time is obtained, which can be seen in Equation (14). In order to calculate the mass change rate during the phase change, the nucleation site volume fraction, the evaporation coefficient, and the condensation coefficient are introduced, which can be expressed as Equation (15). (14) $\begin{array}{rl}\frac{d{R}_B}{dt}& =\sqrt{\frac{2}{3}\frac{P_v-P}{{\rho }_l}}\end{array}$(14) (15) $\begin{array}{rl}\dot{m}& =\{\begin{array}{ll}{F}_{vap}{\displaystyle \frac{3r_{nuc}(1-r_v){\rho }_v}{R_B}}\sqrt{{\displaystyle \frac{2}{3}}{\displaystyle \frac{P_v-P}{{\rho }_l}}}& \mathrm{i}\mathrm{f}\mathrm{P}<{\mathrm{P}}_{\mathrm{v}}\\ F_{cond}{\displaystyle \frac{3r_v{\rho }_v}{R_B}}\sqrt{{\displaystyle \frac{2}{3}}{\displaystyle \frac{|P_v-P|}{{\rho }_l}}}& \mathrm{i}\mathrm{f}\mathrm{P}>{\mathrm{P}}_{\mathrm{v}}\end{array}\end{array}$(15) The parameters recommended for the above equation are $R_B=5\mathrm{E}-7\mathrm{m}$; $r_{nuc}=5\mathrm{E}-4$; $F_{vap}=50$; $F_{cond}=1\mathrm{E}-2$. These coefficients are determined empirically and are applicable only in specific engineering cases. The evaporation ($F_{vap}$) and condensation coefficients ($F_{cond}$) signify the relaxation time for the mutual transition process between the vapour and liquid phases. The relaxation time is derived in the mass transfer equation between the phases, which is based on bubble dynamics and obtained from the derivative of the bubble radius with respect to the time scale. The time scale is expressed as the ratio of the characteristic flow length to the fluid velocity. Consequently, a strong correlation exists between the evaporation and condensation coefficients and the local flow characteristics. Pope (2000) highlighted that the symmetric part of the velocity strain rate correlates directly with the shear stress. By analogy, this symmetric part can substitute for the local flow characteristics, enabling the linkage of the velocity strain rate with the evaporation and condensation coefficients. Due to the typically elevated flow velocity and significant pressure drop gradient in the cavitation region, the velocity strain rate correction facilitates the rupture and fragmentation of bubbles in areas with higher velocity strain rates, thereby diminishing the bubble diameter. Moreover, regions with higher velocity strain rates experience relatively high shear stress, impeding bubble nucleation and growth, and reducing the nucleation site volume fraction. Furthermore, Asnaghi et al. (2017) employed the velocity strain rate correction method in the Sauer cavitation model to examine the surface pressure distribution of a steady-state two-dimensional hydrofoil under cavitation conditions, successfully diminishing the impact of the empirical parameter (number of bubbles per unit volume, n) on the calculation outcomes. Likewise, a prior study by the authors (Li et al., 2023a) effectively eliminated the influence of bubble diameter and nucleation site volume fraction on the mass flow rate of a cavitation jet through an orifice plate, using a velocity strain rate correction within the Z-G-B model.

Building on this modified cavitation model, a characteristic length equation that considers the flow characteristics (valve opening) of the valve is derived. Figure 2 depicts the principle of the calculation of the characteristic length of the relief valve investigated in our study. When the valve opening (x_v) occurs, fluid exits the nozzle and forms a curtain surface shown in blue in Figure 2. This surface corresponds to the ring area with diameter d_v, which represents the characteristic length corresponds to valve opening. (16) $\begin{array}{rl}{F}_{vap-mod}& =F_{vap}\left[1+\frac{d_v}{U_{\mathrm{\infty }}}\left|\frac{1}{2}(\frac{\mathrm{\partial }{U}_i}{\mathrm{\partial }{x}_j}+\frac{\mathrm{\partial }{U}_j}{\mathrm{\partial }{x}_i})\right|\right]\end{array}$(16) (17) $\begin{array}{rl}{d}_v& =\sqrt{4\cdot d_{avg}\cdot x_v}\end{array}$(17) Here $d_v$ is the characteristic length, $U_{\mathrm{\infty }}$ is the velocity of the fluid, d_avg is the average of the inner and outer diameters of the valve disc. and x_v is the valve opening.

Figure 2. The schematic of the calculation of the characteristic length. When the valve is open at position x_v, the space between the valve disc and the nozzle forms a ring-shaped region through which the fluid flows, represented by the transparent blue section. This annular region can be considered as a circular area with a characteristic diameter d_v, which serves as the characteristic length for further analysis.

2.3. Theory of AES-MSI

The use of surrogate models is prevalent in optimising designs as it replaces complicated and time-consuming numerical analyses (Elmisaoui et al., 2023). These models, also known as approximate models or response surface models, significantly enhance efficiency and reduce optimisation difficulty. Common surrogate models include polynomial response surface (PRS), Kriging model (KRG), radial basis function (RBF), and support vector regression (SVR). However, selecting the most suitable surrogate model for a specific problem can be challenging. A model's performance depends not only on its nature but also on the problem and experiment design. There is no single surrogate model that fits all engineering applications. Consequently, using an ensemble of surrogate models is vital to improve predictive capabilities with effective ensemble strategies. In our previous works, the authors proposed an adaptive surrogate model called the Adaptive Ensemble of Surrogate Models by Minimum Screening Index (AES-MSI) model, which demonstrated accuracy and robustness when tested with several test functions (Zhang et al., 2022). Therefore, in this study, AES-MSI is utilised to map the relationships between critical parameters in the Z-G-B cavitation model and cavitation lengths. AES-MSI consists of two parts. In part (I), individual surrogate models in the model library are updated based on the screening index to obtain a model subgroup with relatively high accuracy. In part (II), a weight calculation strategy is developed to determine the weight of each individual surrogate model, and AES-MSI is constructed using the obtained weights. The framework of AES-MSI is illustrated in Figure 3.

Figure 3. The framework of AES-MSI.

2.3.1. Part (I): update the model library

The surrogate models in the original model library mainly consist of four commonly used ones: PRS, RBF, SVR, and KRG. Since the adaptability of these models to different engineering problems varies based on dimensionality, degree of nonlinearity, etc., it is necessary to evaluate their applicability to the research object of this paper. To achieve this, a screening index based on the correlation coefficient and CV is proposed, as shown in Equation (18). (18) $r_{ij}=\frac{\sum _{k=1}^n({{\hat{y}}_i}^k-{\mu }_i)({{\hat{y}}_j}^k-{\mu }_j)}{\sqrt{\sum _{k=1}^n{({{\hat{y}}_i}^k-{\mu }_i)}^2{({{\hat{y}}_j}^k-{\mu }_j)}^2}}$(18) where r_ij is the correlation coefficient between i-th and j-th individual surrogates. ${{\hat{y}}_i}^k$ and ${{\hat{y}}_j}^k$ are the predictions of i-th and j-th individual surrogates, respectively. ${\mu }_i$ and ${\mu }_j$ are the predictions mean value of i-th and j-th individual surrogates, respectively.

Corresponding to the model library constructed in this paper, the correlation coefficient between i-th individual surrogate and the other four individual surrogates can be written as: (19) $r_i=\left[r_{i1}{r}_{i2}{r}_{i3}{r}_{i4}\right]$(19) Then, construct the screening index as follow: (20) $S{I}_i=\frac{{\sum _{j=1}^{n}C{V}_i\left(x_j\right)}_i}{n\cdot max\left\{r_i\right\}}$(20) where SI_i is the screening index of i-th individual surrogates. $max\left\{r_i\right\}$ represents max correlation coefficient in r_i.

A larger $max\left\{r_i\right\}$ indicates the i-th surrogate is highly correlated with an individual surrogate in the model library, which means that the model has relatively better prediction accuracy. Conversely, a smaller MCV_i indicates better performance of the individual surrogate. Therefore, a smaller SI_i value represents a better performance of the individual surrogate. To compare the performance and exclude the influence of large SI_i values of individual surrogates, we calculate the normalised SI_i using the following formula: (21) $NS{I}_i=\frac{S{I}_i-S{I}_{min}}{S{I}_{max}-S{I}_{min}}$(21) where SI_max and SI_min are the maximum and minimum value of SI_i in the model library, respectively.

2.3.2. Part (II): calculation of adaptive weights

This methodology integrates global trends and local errors to determine the weights of each surrogate model in adaptive ensemble. The baseline surrogate model, with the smallest NSI, represents the global trend. The local error is calculated as the difference between the other surrogate models and the baseline surrogate model, as shown in Equation (22). A smaller value of e_ij indicates a greater weight of the individual surrogate in the ensemble. (22) $e_{ij}=|{\hat{y}}_{ij}-{\hat{y}}_{basej}|$(22) where e_ij represents local error at j-th sample point of i-th individual surrogate. ${\hat{y}}_{ij}$ and ${\hat{y}}_{basej}$ are the prediction value of the i-th individual surrogate and the baseline surrogate at the j-th sample point, respectively.

The normalised error is then calculated by combining the local error of each surrogate, including the baseline surrogate. It is noted that e is a vector, which includes the local error of each predicted value. Thus, there is a weight for all predicted values of each surrogate. (23) $l_{ij}=\frac{e_{ij}}{\sum _{m=1}^{M-1}\left(e_{mj}\right)}$(23) The weight of each surrogate in the ensemble is composed of two parts: a global trend weight represented by the baseline surrogate, and a weight based on the local errors. To balance these two parts, we set the weights of the baseline surrogate. (24) $\begin{array}{rl}{\omega }_{ij}& =\frac{\frac{1}{l_{ij}}}{2\sum _{m=1}^{M-1}\left(\frac{1}{l_{mj}}\right)}\end{array}$(24) (25) $\begin{array}{rl}{\omega }_{bj}& =0.5\end{array}$(25) where ${\omega }_{bj}$ and ${\omega }_{ij}$ are the weights of baseline surrogate and the rest of the surrogates at j-th sample point.

3. CFD cavitation simulation of the safety relief valve

In this section, we present a comprehensive description of the numerical simulation of the cavitation phenomenon in relief valves. Firstly, we discuss the simplification and meshing process of the computational domain of the relief valve. Next, we define the cavitation length as a criterion to quantify the cavitation phenomenon, taking into account the specific characteristics of the relief valve under study. Additionally, we construct a mesh-independent validation by examining this parameter, demonstrating that our proposed numerical method is not sensitive to the distribution of the computational grid. We also provide detailed information about the boundary conditions, which are set based on actual working conditions. Furthermore, we verify the effectiveness of the physical dimensionality reduction achieved through the use of the modified evaporation coefficient, which is corrected based on the velocity strain rate. The bubble diameter and nucleation site density serve as indicators to validate this reduction. Lastly, we establish a relationship between the evaporation/condensation coefficients and cavitation lengths, employing the AES-MSI method. Grayscale detection techniques are applied to visualise experiments, and the resulting cavitation lengths are utilised as targets for parameter optimisation using a Genetic Algorithm. Through these steps, we develop a high-fidelity cavitation model specifically tailored for relief valves.

3.1. Computational domain simplification and mesh generation

In our previous research, the computational domain modelling and mesh generation techniques were successfully employed to analyze safety valves in the second circuit of nuclear power plants for predicting fluid forces acting on the valve disc (Li et al., 2023b). In this study, we extend our previous work by increasing the valve inlet pressure and directing our attention to cavitation inside the valve. Since the only difference between the two operating conditions is the pressure, we can utilise the same fluid domain partitioning and mesh generation techniques used in the previous work. Figure 4(a) presents a 3D model of a direct-acting relief valve comprising five essential components: preload screws, spring, valve seat, valve body, and valve disc. In Figure 4(b), a simplified fluid domain model based on the 3D model depicted in Figure 4(a) is shown. To generate a pure hexahedral mesh, the fluid computational domain is divided into four crucial sections, namely the inlet, turbulent core, valve chamber, and outlet domains. The interfaces facilitate data transfer between different sub-fluid domains. Figure 4(c) illustrates the mesh that has been generated for the fluid domain. To enhance computational efficiency, all meshes have been partitioned into hexahedral meshes. Moreover, the O-block method has been employed to generate the mesh while ensuring that it satisfies the quality criteria for meshes with circular boundaries.

Figure 4. The calculation model configuration of the safety relief valve. (a) presents a semi-profile 3D model of the relief valve, offering a clear depiction of its key components. (b) presents the fluid computational domain obtained by employing the split computational domain processing method. Various sub-domains are highlighted using distinct colours, and data exchange between these sub-domains is facilitated through Interfaces. (c) showcases the resulting mesh division, where hexahedral meshes are employed throughout using the ‘O-block’ approach. This strategy not only optimises computational resources but also enhances computational convergence.

3.2. Independence verification of mesh and turbulence models

Based on the computational domain and mesh method described in the previous section, this section presents details of the grid independence verification and turbulence models. Numerical simulations’ outcomes are heavily impacted by these two factors. Two critical factors in the grid independence verification are the first layer height near the wall and the mesh density. The first layer height is crucial for accurate prediction of the flow near the wall, and its reasonableness is often evaluated using y⁺. In our study, the SST k-ω model is employed, and it requires the y⁺ value of the valve disc in the cavitation region to be close to 1. Therefore, a first layer height of 1E-6 m is selected for this region, and the distribution of y⁺ of the valve disc is presented in Figure 5. Subsequently, we performed grid independence verification using the cavitation length as the target, selecting five different grid sizes (2, 1, 0.5, 0.2, and 0.1 mm) based on the first layer height.

Figure 5. The distribution of y⁺ of the valve disc.

The cavitation region of the relief valve under investigation forms an approximate annular cavitation zone with a half-elliptical cross-section around the axis of the valve (as shown in Figure 6). In our study, we define two cavitation lengths: the radial cavitation length (x_length) and the axial cavitation length (y_length). The x_length is calculated as the difference between the length at the start and end of cavitation, using Equation (26). The y_length is computed from the bubble's volume, which is the volume of the bubble swept by half of the elliptical section around the centre of the section, using Equation (27). The mesh independence verification was conducted under specific boundary conditions, including 1.5 MPa inlet pressure, 0 MPa outlet pressure, and 0.5 mm valve opening, while keeping the parameters of the Z-G-B cavitation model at their default values. Table 1 compares the cavitation length under five different mesh densities. It shows that with increasing mesh density, both x_length and y_length slightly increase. However, when the mesh density reaches 0.2 mm, the two cavitation length stabilise with a variation rate of less than 1%. Hence, 0.2 mm is considered the optimal mesh density and is adopted as the final mesh density. (26) $\begin{array}{rl}{x}_{length}& ={r_v}_{apor-end}-r_{vapor-start}\end{array}$(26) (27) $\begin{array}{rl}{y}_{length}& =\frac{V_{vapor-volume}}{{\pi }^2\cdot x_{length}\cdot (r_{vapor-end}+r_{vapor-start})}\end{array}$(27) Turbulence models are crucial for accurately predicting cavitation in engineering applications. Accurate prediction of cavitation requires accounting for turbulence effects, particularly in modelling complex flows with vortices or high shear flows. Therefore, selecting and validating an appropriate turbulence model for specific flow conditions is essential to ensure accurate predictions of cavitation phenomenon. In this study, we discuss five different turbulence models (Standard k-ε, RNG k-ε, Realisable k-ε, Standard k-ω, and SST k-ω) to determine the most suitable turbulence model. Figure 7 shows the cavitation lengths (x_length and y_length) predicted by various turbulence models for the analyzed relief valves. Based on Figure 7, the cavitation length predictions by the turbulence models can be grouped into two categories: the k-ε turbulence model and the k-ω turbulence model. Among the turbulence models developed from the k-ε model, the RNG k-ε model predicts the longest cavitation length, while the Standard k-ε model predicts the smallest. As for the turbulence model based on the k-ω model, both the x_length are 3.5 mm, while the y_length are larger than those of the Standard k-ω model. According to the comparison with experimental results, the cavitation lengths predicted based on the k-ω model are more accurate. This finding will be further emphasised in the subsequent experimental section. Therefore, we selected the SST k-ω model to predict the cavitation length due to its accurate prediction in our study.

Figure 6. Represents the radial dimension (x_length) and axial dimension (y_length) of the cavitation area in the relief valve. The blue region represents a portion of the fluid domain surrounding the valve disc and seat. As cavitation exclusively transpires within the gap formed between the valve disc and seat, this specific area is magnified for clarity. Meanwhile, the red area denotes the cavitation region, taking the form of a ring-shaped region obtained by sweeping an ellipse as a cross-section with the ellipse's centre serving as the radius.

Figure 7. The effect of different turbulence models on cavitation length. The figure depicts the cavitation dimensions with different turbulence models, with the red bars representing the radial dimension, x_length, and aligned with the left axis. Conversely, the blue bars represent the axial dimension, y_length, corresponding to the right axis.

Table 1. The comparison of cavitation length with different numbers of meshes.

The size of mesh[mm]	Nodes	Elements	xlength[mm]	ylength[mm]	Error [xlength]	Error [xlength]
2	1,046,832	986,971	3.225	0.173	7.80%	12.63%
1	1,095,655	1,031,464	3.331	0.188	4.77%	5.05%
0.5	1,366,246	1,289,520	3.430	0.195	1.94%	1.52%
0.2	4,336,399	4,199,519	3.500	0.198	′0.06%	0
0.1	20,885,032	20,588,675	3.498	0.198	–	–

3.3. Boundary conditions

Boundary conditions define the behaviour of fluids at the edges or boundaries of the analyzed domain, including walls, inlets, outlets, and interfaces. They dictate the flow parameters, such as velocity, pressure, and species concentration, at these boundaries. Properly selecting and specifying boundary conditions is crucial to ensure the accuracy and reliability of CFD simulations. Inaccurate or inappropriate boundary conditions can lead to erroneous outcomes or convergence issues, rendering the simulation invalid. Therefore, this section presents a comprehensive description of the boundary conditions based on real working conditions.

Half of the model is used to reduce computational resources due to its symmetry, and a symmetric boundary is employed for the fluid domain model. Pressure boundary conditions are set for both the inlet and outlet. The mixture model describes the phase transition process in the cavitation, where the primary phase is oil, and the secondary phase is oil-vapour. Table 2 presents the physical properties of the two media. User Defined Functions (UDFs) are employed to compile the two-phase flow mechanism based on the Z-G-B cavitation model's mass transfer equations, with the mass transfer term modified by the velocity strain rate. Equations (15)–(17) present these equations. The acceleration in the y direction is set to 9.81 m/s² to account for the fluid's gravity and the formation of bubbles. The coupled solver is used due to its robustness and accelerated convergence, especially for high-pressure cavitation jets. Second-order upwind schemes are used for discretization of momentum, turbulence kinetic energy, and turbulence dissipation rate. To monitor convergence, the continuity equation, x-directional velocity, y-directional velocity, z-directional velocity, turbulent kinetic energy, and turbulent dissipation rate are considered. Convergence is achieved when the residual values of these physical variables fall below 1E-4. The solution algorithm and spatial discretization method are listed in Table 3.

Table 2. Material properties of the flow medium.

Object	Viscosity[kg/m·s]	Density[kg/m^3]
Oil	0.042	872
Oil-vapour	6E-6	0.13

Table 3. The solution algorithm and the spatial discretization method.

Object	Method
Solution Scheme	Coupled
Gradient term	Least Squares Cell-Based
Pressure term	Second Order Upwind
Momentum term	Second Order Upwind
Turbulence kinetic energy	Second Order Upwind
Convergence criterion	1E-4

3.4. Validation of the implementation effect of physical dimensionality reduction

As discussed in Section 2.2, implementing the velocity strain rate correction to the mass transfer term in the Z-G-B cavitation model, in conjunction with the valve flow length, can significantly reduce physical dimensionality. This section aims to validate the effectiveness of this approach by examining the cavitation length, specifically the bubble diameter (D_b) and nucleation site volume fraction (α_nuc). Figure 8 compares the cavitation length predicted by the modified and original cavitation models for the bubble diameter (D_b). The radial and axial cavitation lengths are represented by Figure 8(a and b), respectively. The red and blue dotted lines in Figure 8(a) represent the variation curves of x_length calculated by the original and modified cavitation models, respectively. The red and blue bars show the change rate of the original and modified cavitation models compared to the recommended bubble diameter (D_b = 1E-6 m), respectively. The results indicate that x_length decreases as D_b increases. The gradient of change in the cavitation length is more drastic when D_b is larger than the recommended value, and more moderate when it is smaller. The modified cavitation model effectively eliminates the effect of the D_b on the x_length, reducing the maximum error from 57% to 11.2%. In Figure 6(b), the axial cavitation lengths (y_length) predicted by both models decrease when the D_b is less than the recommended value. However, when the D_b is larger than the recommended value, the two models predict opposite trends. The error of the modified cavitation models are smaller than those of the original cavitation model, reducing the maximum error from 12.5% to 8.7%. Overall, the modified cavitation model effectively mitigates the influence of the bubble diameter (D_b) on the prediction of the cavitation length.

Figure 8. The effect of the modified cavitation model on the cavitation length compared with the original model for the bubble diameter. (a) and (b) represent the comparison between the original cavitation model and the modified cavitation model in terms of x_length and y_length for different bubble diameters, respectively. The bar charts in red and blue represent the cavitation lengths of the original and modified models, respectively, and are aligned with the left-hand axis. The dotted line plots in red and blue represent the cavitation lengths of the original and modified models, respectively, and are aligned with the right-hand axis.

Figure 9 illustrates a comparison between the cavitation length predicted by the modified and original cavitation models for the α_nuc. The dotted line plots and bar graphs represent the same features as in Figure 8, with the exception of the independent variable (Figure 9 for α_nuc). From Figure 9(a), it is evident that the x_length predicted by the original cavitation model gradually increases with an increase in α_nuc, and the change is most significant when α_nuc is less than the recommended value (α_nuc = 5E-4). In contrast, the modified cavitation model shows only a slight decrease with an increase in α_nuc. A comparison of the errors between the predicted and recommended values of the two cavitation models under different values of α_nuc reveals that the modified cavitation model effectively suppresses the influence of α_nuc on the x_length prediction and reduces the maximum error from 63% to within 5%. As shown in Figure 9(b), the prediction of y_length for both cavitation models increases with increasing α_nuc, but the change in y_length prediction for the original model is more pronounced when α_nuc is less than the recommended value (α_nuc = 5E-4). A comparison of the bars reveals that the errors calculated by the original model relative to the recommended values at different α_nuc are greater than those predicted by the modified cavitation model, and the maximum error is reduced from 22% to within 10%. In summary, the modified cavitation model effectively suppresses the effects of D_b and α_nuc variations on the prediction results of the cavitation length, and achieves a physical dimensionality reduction of the cavitation model.

Figure 9. The effect of the modified cavitation model on the cavitation length compared with the original model for the nucleation site volume fraction. (a) and (b) represent the comparison between the original cavitation model and the modified cavitation model in terms of x_length and y_length for different nucleation site volume fraction, respectively. The bar charts in red and blue represent the cavitation lengths of the original and modified models, respectively, and are aligned with the left-hand axis. The dotted line plots in red and blue represent the cavitation lengths of the original and modified models, respectively, and are aligned with the right-hand axis.

3.5. The construction of AES-MSI based on the evaporation and condensation coefficients

Based on the analysis in Section 3.4, the impact of D_b and α_nuc on the cavitation length is effectively reduced by correcting the evaporation term using the velocity strain rate correction. This section explores the relationship between the evaporation/condensation coefficients and cavitation lengths using the AES-MSI method. The evaporation coefficient's default value is 50. To explore a range of values, we varied it up or down by an integer factor from the default value, defining the sample space as [5E-1, 5E1]. The condensation coefficient has a default value of 1E-2, and we established a sample space of [1E-2, 1E0] by varying the coefficient up or down by an integer factor from its default value. To establish the relationship between evaporation/condensation coefficients and cavitation lengths, it is essential to effectively sample the design variables within the designated sample space. Subsequently, the cavitation lengths corresponding to the sample points are computed using CFD numerical models. The present study utilises the optimal Latin hypercube sampling (LHS) method to sample the evaporation and condensation coefficients. Since the design variables in this study are multiples of 10 and the ranges have an extreme difference of a factor of 100, the design space undergoes logarithmic transformation to achieve the highest possible sampling uniformity.

Figure 10 displays the distribution of the training points and test points in the sample space using the logarithmic coordinate system. The study uses a dataset of 34 points, including 30 uniformly distributed sample points and 4 boundary points, represented by black and blue circles, respectively. The training points are chosen based on the principle criteria proposed by Forrester and Jones, which require a number of training points exceeding ten times the number of dimensions. Additionally, 10 uniformly distributed validation points, represented by red circles, are employed to evaluate the surrogate model's accuracy. The CFD calculations for both the training and validation datasets utilise specific boundary conditions, namely, a valve opening of 0.5 mm and an inlet pressure of 1.5 MPa. The remaining settings align with those described in Section 3.2–3.3.

Figure 10. The distribution of the training points and validation points in the sample space. The black points in the figure represent the training points utilised to construct the surrogate model. The red points indicate the validation points employed to assess the accuracy of the surrogate model. To avoid extrapolation, the authors introduced four additional boundary points, illustrated by the blue points.

The relationship between evaporation/condensation coefficients and cavitation length is established using the AES-MSI method, as described in Section 2.3, based on the CFD numerical simulation results at specific sample points. The evaporation and condensation coefficients are treated as input variables, while the cavitation lengths are considered as the output variable. To assess the performance of AES-MSI, four commonly used engineering surrogate models (PRS, RBF, KRG, and SVR) are compared. Since relying on a single indicator for evaluation can be limiting, this study employs three evaluation indicators to assess the surrogate model's accuracy. These indicators are the coefficient of determination (R²), root mean square error (RMSE), and mean absolute percentage error (MAPE). R² is a widely used metric to visually evaluate the model fit, ranging between 0 and 1, with values closer to 1 indicating a better fit to the data. Both RMSE and MAPE visually represent the magnitude of the model's prediction error, with smaller values indicating better prediction performance. The results of the prediction accuracy evaluation for the five surrogate models are presented in Tables 4 and 5. The surrogate model constructed by AES-MSI have the highest R² value (closest to 1), the smallest RMSE and MAPE values, demonstrating the accuracy of AES-MSI in constructing the mapping relationship between the evaporation/condensation coefficient and cavitation lengths.

Table 4. Three evaluation indexes of five surrogate models for the prediction of x_length.

Object	R^2	RMSE	MAPE
KRG	0.910	0.030	0.535
PRS	0.902	0.031	0.623
RBF	0.869	0.036	0.883
SVR	0.513	0.145	3.365
AES-MSI	0.933	0.026	0.390

Table 5. Three evaluation indexes of five surrogate models for the prediction of y_length.

Object	R^2	RMSE	MAPE
KRG	0.876	0.007	0.753
PRS	0.901	0.006	0.564
RBF	0.918	0.006	0.592
SVR	0.652	0.032	3.792
AES-MSI	0.932	0.005	0.497

The constructed surrogate model provides a rapid way to access the cavitation length at specific evaporation and condensation coefficients, which is crucial for determining the optimal coefficients for both processes. Figure 11 show the surface correlation among the cavitation lengths, evaporation coefficient, and condensation coefficient, generated using AES-MSI. Figure 11(a and b) depict x_length and y_length as functions of evaporation and condensation coefficients, respectively. From Figure 11(a), it can be observed that the relationship between evaporation coefficient, condensation coefficient, and x_length is bowl-shaped. When either the evaporation or condensation coefficient remains constant while the other parameter gradually increases, x_length increases in a quadratic trend with a decreasing gradient of change, indicating a positive correlation between x_length and evaporation/condensation coefficients. The maximum x_length value appears in the upper left corner of Figure 11(a) (looking at the x_length – Lg(F_cond) plane), corresponding to the position where the evaporation/condensation coefficient are at their maximum. Conversely, the minimum x_length value appears in the lower right corner of Figure 11(a) (looking at the x_length – Lg(F_cond) plane), corresponding to the position where the evaporation/condensation coefficient are at their minimum. The relationship between evaporation and condensation coefficients and y_length shown in Figure 11(b) differs slightly from Figure 11(a). The maximum value appears at the same location as in Figure 11(a), although the y_length in the lower right corner is smaller. However, the global minimum value appears in the lower left corner. y_length tends to increase as the condensation coefficient decreases when the evaporation coefficient is small, while y_length tends to decrease as the condensation coefficient decreases when the evaporation coefficient is larger.

Figure 11. The surface correlation among the cavitation lengths, evaporation coefficient, and condensation coefficient. (a) and (b) represent the evaporation/condensation coefficients plotted against x_length and y_length, respectively, using a logarithmic scale for improved visualisation. The colour gradient ranging from red to blue indicates the variation of the coefficients, with red indicating the maximum value and blue indicating the minimum value. The accompanying colour bar provides a reference for the magnitude of the coefficients.

3.6. Visualisation experiments and genetic algorithms for determining evaporative condensation coefficients

Once the relationship between the evaporation/condensation coefficients and cavitation lengths has been established, it becomes crucial to determine these coefficients based on the cavitation lengths obtained under actual operating conditions. This section provides a detailed description of the experimental procedure and the independently built test rig used to measure cavitation phenomena. The cavitation lengths measured during the experiment are used as the target to find the optimal values for the evaporation/condensation coefficients of the relief valve using genetic algorithms. In our previous study (Li et al., 2023b), we conducted experimental measurements on the fluid forces exerted on the valve disc of a relief valve using a hydraulic bench. These experiments yielded high test accuracy and demonstrated good experimental repeatability. Building upon this hydraulic test rig setup, we further investigated the cavitation phenomenon occurring within the valve by integrating visualisation equipment. The hydraulic test rig are depicted in Figure 12.

Figure 12. The test rig for observing cavitation. (a) depicts the comprehensive layout of the experimental setup, with the control cabinet, represented by the white cabinet, serving as the central unit for oil pressure regulation. The oil supply system comprises the motor and oil tank, while the monitor enables the visualisation of the cavitation area surrounding the valve disc and valve seat. (b) represents the visualisation system, encompassing the CCD camera, industrial microscope, and the valve system under examination. (c) showcases the data acquisition system, which consists of the signal amplifier and an NI acquisition card. (d) depicts the valve opening adjustment mechanism, primarily comprises a linear motion system and a laser displacement sensor. (e) provides a detailed view of the transparent test valve and the valve opening adjustment mechanism.

The experimental setup consists of various components, such as a tank and motor-driven oil supply system, an oil pressure control system based on PID regulation, a data acquisition system using NI data acquisition cards, and a visualisation system that includes a transparent test valve and a high-definition (HD) industrial microscope. This study specifically focuses on the visualisation experimental system, which includes transparent test valves, an electronic magnifier system, a valve opening adjustment mechanism, and upstream and downstream pressure monitoring systems.

The transparent test valves are made of PMMA material, providing excellent visualisation capabilities to observe the flow conditions inside the valves. The visualisation system consists of a high-definition industrial microscope, a high-performance CCD camera, and brackets. The system allows adjustable optical magnification ranging from 21x to 135x, adjustable lens magnification from 0.7x to 4.5x, and can capture images with up to 4800 w pixels. The valve opening adjustment mechanism primarily consists of a rack-and-pinion-based linear motion mechanism and a laser displacement sensor. The linear motion mechanism features a scale bar with an adjustment accuracy of 0.01 mm. By rotating the adjustment knob, the valve opening can be adjusted within a maximum range of 5 mm, and the opening can be fixed using the locking nut. The laser displacement sensor measures the valve disc opening by projecting a laser onto the iron piece mounted on the spool, allowing for precise verification of the valve opening set by the linear mechanism.

The experimental procedure is exemplified using a valve opening of 0.5 mm and an inlet pressure of 1.5 MPa. Firstly, the power should be turned on, and a small pressure should be entered to allow the equipment to run for a while. This will ensure that there are no fluid leaks in the pipeline and that the pressure values displayed on the sensor match those from the data acquisition system. In the second step, the initial position of the adjustment mechanism relative to the scale when the valve disc is at 0 mm valve opening should be confirmed. The adjusting screw should then be rotated to adjust the displacement of the valve spool to 0.5 mm, and the laser displacement sensor value should be compared to confirm the accuracy of the valve opening. In the third step, the pressure sensor in front of the valve should be stable, and the cavitation image detected from the display should not change. Once this condition is met, the cavitation image should be recorded, marking the end of the cavitation test process. Figure 13(a) presents the cavitation region of the relief valve, observed through an electron magnifying lens. In this image, the yellow area represents oil, the upper portion of the bright gray area corresponds to the valve disc, the lower portion represents the nozzle, and the black area below the valve disc indicates the cavitation zone. To determine the cavitation length, our study employs digital image processing techniques based on grayscale detection, commonly used for identifying two-phase flows (Liu et al., 2019; Soyama & Hoshino, 2016). To reduce the number of pixels requiring calculation, images are cropped to appropriate regions. The average gray value (G) is determined by calculating the average value of the gray values along the same horizontal coordinate, as shown in Equation (28). And a threshold for defining the cavitation initiation and collapse termination points is defined follows the method proposed by Liu et al. (2019), namely, λ_s = (λ_max + λ_min)/2. If the gray value (λ) is greater than or equal to λ_s, the distance from the collapse to the starting boundary represents the cavitation length. (28) $\begin{array}{rl}& \lambda =\frac{\mathrm{G}}{256}\end{array}$(28) (29) $\begin{array}{rl}& s.t.\{\begin{array}{l}min\left(|x_{length}-3.67|\right)\\ min\left(|y_{length}-0.22|\right)\\ -1\le \mathrm{L}\mathrm{g}({\mathrm{F}}_{\mathrm{v}\mathrm{a}\mathrm{p}}/5)\le 1\\ -2\le \mathrm{L}\mathrm{g}({\mathrm{F}}_{\mathrm{c}\mathrm{o}\mathrm{n}\mathrm{d}})\le 0\end{array}\end{array}$(29) Figure 13(b) depicts the grayscale imagine obtained using the aforementioned method, where the white region represents the cavitation area determined through grayscale detection. In Figure 13(c), the grayscale values are plotted against the x-axis coordinates to illustrate their variation. Subsequently, we present the dimensional measurements of the cavitation region. It exhibits a radial length of 3.67 mm (referred to as x_length) and an axial length of 0.22 mm (referred to as y_length). Next, the genetic algorithm is employed to optimise the evaporation and condensation coefficients based on the relationship constructed by AES-MSI. The optimisation objectives and constraints are presented in Equation (29). Specifically, the optimisation of evaporation and condensation coefficients involves two objectives and two constraints. The two objectives are to minimise the difference between the experimental values and the predicted values of x_length and y_length, which are surrogate model outputs. The two constraints are to limit the range of values of evaporation and condensation coefficients. In this context, x_length represents the radial cavitation dimension that influences the cavitation length on the valve disc in the engineering problem. Moreover, y_length represents the axial cavitation length that directly affects the mass flow rate of the valve by obstructing the passage of fluid. Given these objectives and constraints, both objectives are weighted equally at 0.5 to select the optimal solution from the pareto solution set. Accordingly, the final evaporation and condensation coefficients are determined to be 5E0.757 and 1E-1.830, respectively. Figure 14 presents the vapour volume contour (the cavitation region) of the relief valve under seven different calculation conditions, including the Standard k-ε, Realisable k-ε, RNG k-ε, SST k-ω, Standard k-ω, MSST k-ω (a modified cavitation model corrected by velocity strain rate), and OMSST k-ω (an optimised cavitation model combining the flow characteristics of the valve, velocity strain rate correction, AES-MSI surrogate model, and genetic algorithm). Since cavitation only occurs in the valve disc region, Figure 14 displays only the valve disc region, with the cavitation region highlighted in red. The highlighted area represents an iso-surface with a gas-phase volume fraction greater than 0.5.

Figure 13. The method of cavitation length calculation. (a) is the cavitation region of the relief valve, observed through an electron magnifying lens. In this image, the yellow area represents oil, the upper portion of the bright gray area corresponds to the valve disc, the lower portion represents the nozzle, and the black area below the valve disc indicates the cavitation zone. (b) represents the grayscale map obtained using the grayscale detection method, where the bright white region represents the cavitation area. (c) represents the variation of the grayscale values with the x-axis coordinates to illustrate the cavitation length.

Figure 14. Comparison of cavitation regions for different turbulence models. (a) Standard k-ε. (b) Realisable k-ε. (c) RNG k-ε. (d) SST k-ω. (e) Standard k-ω. (f) MSST k-ω. (g) OMSST k-ω.

Table 6 shows the cavitation lengths predicted by these seven calculation methods. The results indicate that k-ω-based models predict cavitation lengths that are in better agreement with experimental values than k-ε-based models. Among the four calculation methods based on the k-ω turbulence model, our proposed cavitation prediction method exhibits the highest prediction accuracy, with prediction errors of x_length and y_length reduced to 0.81% and 1.72%, respectively. Furthermore, to further validate the accuracy of the simulation, we compared the cavitation lengths obtained from experiments and simulations under two operating conditions: valve openings of 0.7 and 0.9 mm. These evaporation/condensation coefficients were chosen based on the optimised coefficients mentioned above. The comparison results are presented in Figure 15.

Figure 15. Comparison of experimental and simulated cavitation lengths at different valve openings. (a1)-(a3) represent the simulated cavitation length, the experimental cavitation length, and grayscale curve at 0.7 mm opening, respectively. (b1)-(b3) represent the the simulated cavitation length, the experimental cavitation length, and grayscale curve at 0.9 mm opening, respectively.

Table 6. Three evaluation indexes of five surrogate models for the prediction of y_length.

Object	xlength	ylength	Error (xlength)	Error(ylength)
Standard k-ε	0.75	0.091	79.62%	58.7%
Realisable k-ε	1.10	0.114	70.11%	48.4%
RNG k-ε	2.10	0.121	42.93%	45.2%
SST k-ω	3.50	0.198	4.89%	9.94%
Standard k-ω	3.50	0.212	4.89%	3.73%
MSST k-ω	3.59	0.224	2.56%	1.76%
OMSST k-ω	3.65	0.224	0.81%	1.72%

4. Analysis of cavitation evolution based on modified cavitation model

In this section, we investigate the impact of two crucial operating parameters, namely the pressure difference between inlet and outlet (hereafter referred to as $\mathrm{\Delta }\mathrm{P}$) and valve opening, on cavitation in relief valves. Using the high-fidelity CFD cavitation model, we provide insights into the influence of $\mathrm{\Delta }\mathrm{P}$ on valve cavitation by examining factors such as fluid flow obstruction caused by vortices and turbulence. This analysis is supported by velocity contours and streamlines. Additionally, employing the same high-fidelity CFD cavitation model, we simulate the cavitation phenomenon for various valve openings and analyze the corresponding changes in the cavitation area. To facilitate this analysis, we propose a flow area calculation formula.

4.1. Effect of $\mathrm{\Delta }\mathrm{P}$ on the cavitation evolution

$\mathrm{\Delta }\mathrm{P}$ is a critical factor that affects the cavitation phenomenon of safety valves. Variations in $\mathrm{\Delta }\mathrm{P}$ can result in different fluid flow velocities through valve openings, thereby creating diverse local low-pressure zones and inducing varying degrees of cavitation. By conducting a comprehensive investigation of the impact of $\mathrm{\Delta }\mathrm{P}$ on the cavitation phenomenon, it is possible to optimise the design parameters and operating conditions of safety valves to ensure the secure operation of equipment. The investigation for various $\mathrm{\Delta }\mathrm{P}$ is conducted using a valve opening of 0.5 mm, while keeping the other parameter settings for the numerical simulations consistent with those presented above. Figure 16 displays the cavitation contours for various $\mathrm{\Delta }\mathrm{P}$ ranging from 0.5 to 1.5 MPa. For space constraints, only four of them (0.6, 0.8, 1.0, and 1.2 MPa) are presented in the main text, while the remainder will be presented in Appendix A. As evident from Figure 16 and Appendix A, the area of cavitation (highlighted in red) occurs in the gap formed between the valve disc and seat. The reduction in fluid flow area at this point causes an increase in flow velocity, leading to a local pressure drop below the saturation vapour pressure, ultimately resulting in cavitation. Observing these cavitation clouds reveals that the cavitation region gradually expands with an increase in $\mathrm{\Delta }\mathrm{P}$. However, the cavitation area does not change beyond a certain critical pressure, and the degree of cavitation differs between the left and right sides of the valve disc due to the valve's asymmetric structure, having only one side outlet. To gain a deeper understanding of the correlation between the cavitation region and pressure, the cavitation length defined in Section 3.2 is employed to measure the extent of cavitation at varying pressure drop, as presented in Figure 17.

Figure 16. Cavitation contours for various $\mathrm{\Delta }\mathrm{P}$. (a) 0.6 MPa. (b) 0.8 MPa. (c) 1.0 MPa. (d) 1.2 MPa.

Figure 17. Cavitation length for various $\mathrm{\Delta }\mathrm{P}$. The figure categorises the cavitation length into two groups based on their proximity to the outlet on either side of the valve disc, as indicated by the red and black colour. By examining the characteristics of the curves, the cavitation phenomenon can be classified into three distinct stages: the initial cavitation stage, the development cavitation stage, and the stable cavitation stage.

In previous literature, the valves cavitation process with varying pressure has been categorised into three stages: the no-cavitation stage, the initial cavitation stage, and the stable cavitation stage. In this study, the three stages mentioned above have been further refined based on the varying degrees of cavitation around the valve disc. Figure 17 considers the rate of change of the radial cavitation length (x_length) and the difference in x_length between the two sides of the valve disc (i.e. close/far away from the valve outlet) to classify cavitation into three stages: initial cavitation, developing cavitation, and stable cavitation. When $\mathrm{\Delta }\mathrm{P}$ is less than 0.7, if the x_length is less than 0.1 mm, which considered that x_length remains unchanged with $\mathrm{\Delta }\mathrm{P}$. Additionally, x_length is consistent on both sides of the valve disc within this range of $\mathrm{\Delta }\mathrm{P}$, thus defining it as the initial stage of cavitation. When $\mathrm{\Delta }\mathrm{P}$ is greater than 0.7 and less than 1.3, the x_length exceeds 0.1 mm, which exhibits a significant increasing gradient with $\mathrm{\Delta }\mathrm{P}$. At this time, x_length on the left and right sides of the valve disc starts to differ, and this difference progressively amplifies. Consequently, this range of $\mathrm{\Delta }\mathrm{P}$ is designated as the development stage of cavitation. When $\mathrm{\Delta }\mathrm{P}$ exceeds 1.3, the rate of change of x_length is less than 2%, which remains constant with $\mathrm{\Delta }\mathrm{P}$. Additionally, the x_length on the left and right sides of the valve disc tends to equalise. Hence, this range of $\mathrm{\Delta }\mathrm{P}$ is categorised as the cavitation stabilisation stage.

To investigate the reasons why cavitation lengths on both sides of the valve disc differ in the developing and stable cavitation stages, yet ultimately show consistency, and why the scales gradually increase and then stabilise with rising pressure, we analyzed the flow characteristics of the valve at various pressures based on the high-fidelity CFD cavitation model. Figure 18 displays the velocity contours, streamlines, pressure contours and vapour contours of the valve at 0.7, 1.2, and 1.5 MPa. From the velocity contours, it can be observed that the flow velocity at the three different $\mathrm{\Delta }\mathrm{P}$ follow a similar trend. As the fluid enters the nozzle, the flow area decreases, leading to an increase in velocity. As the fluid exits the nozzle and enters the valve disc groove region, the through-flow area increases, leading to a decrease in flow velocity. When the fluid reaches the bottom of the valve disc groove, it is divided into two pathways based on the centre of the valve disc (represented by the black centre line). The fluid is accelerated as it passes through the bottom of the valve disc groove (black solid line) and flows into the gap formed by the valve disc and seat (jet flow region) through the valve disc wall (purple solid line). In this region, the flow velocity increases due to the reduction in flow area, leading to a relatively large pressure drop (as shown in Figure 18(a3–c3)) and the formation of a cavitation area (as shown in Figure 18(a4–c4)). The velocity contours (as shown in Figure 18(a1–c1)) at various pressures with the primary difference being the size of the high-speed region (as indicated by the black line) that increases with higher pressure. At pressures of 0.7 and 1.5 MPa, the size of the high-speed region on both sides of the valve disc is equal, whereas at 1.2 MPa, the size of the high-speed region on the left side is smaller than that on the right side. A larger high-speed region corresponds to a greater pressure gradient drop, resulting in a larger cavitation length. This is consistent with Figure 17, which shows that the size of the cavitation region on either side of the valve disc is different at 1.2 MPa but the same at 0.7 and 1.5 MPa.

Figure 18. The velocity contours, streamlines, pressure contours and vapour contours of the valve for three kinds of $\mathrm{\Delta }\mathrm{P}$. where (a), (b) and (c) represent the three $\mathrm{\Delta }\mathrm{P}$ (0.7, 1.2, and 1.5 MPa), respectively. Numbers 1, 2, 3 and 4 represent velocity contours, streamlines, pressure contours and vapour contours, respectively.

The analysis above indicates that differences in flow characteristics are the primary cause of varying cavitation lengths on either side of the valve disc. To investigate the causes of the different high-speed zones on either side of the valve disc, Figure 18(a2–c2) shows the streamlines at three typical $\mathrm{\Delta }\mathrm{P}$. At a $\mathrm{\Delta }\mathrm{P}$ of 0.7 MPa, the fluid flows through the gap formed between the valve disc and the valve seat on the left side and then moves upwards along the wall. Due to the shear force of the wall and the viscous force of the fluid, the fluid gradually decelerates and forms a vortex, which subsequently flows to the outlet through the upper part of hole 1. As the $\mathrm{\Delta }\mathrm{P}$ increases to 1.2 MPa, the flow velocity of the fluid on the left side also increases, resulting in the formation of more vortices, and hence, obstruction of the flow (Martin & Bhushan, 2014). This obstructive effect leads to a shorter high-speed area on the left side of the valve disc. Meanwhile, Most of the fluid on the right side of the valve disc flows directly from hole 1, and only a small portion of it flows through hole 2. As a result, the high-speed area on the right side of the valve disc is larger than that on the left side, leading to different cavitation lengths on the two sides of the valve disc. With further increase in pressure, the fluid on the left side of the valve disc no longer moves upward along the wall due to the higher speed. Instead, it directly impacts the wall and flows into the lower part of hole 1, as shown in Figure 18(c2), resulting in an increased flow capacity due to the lack of vortex obstruction. At this stage, the left and right sides of the valve disc present almost the same high-speed area size, leading to the formation of the same cavitation length.

By employing our developed high-fidelity CFD cavitation model, together with velocity contours and streamlines, we have discovered a strong correlation between the length of cavitation and the extent of the high-speed zone in the gap region (formed by the valve disc and seat). Moreover, it has been determined that high-speed flow plays a critical role in initiating turbulence, wherein turbulence intensity serves as a valuable parameter for assessing and controlling turbulence dynamics. Considering the aforementioned perspective, it becomes feasible to establish a correlation between turbulence intensity and cavitation. Turbulence intensity is defined by the ratio of the root-mean-square of velocity fluctuations to the mean flow velocity, as expressed in Equation (30). (30) $I\equiv \frac{U^{\mathrm{\prime }}}{U_{avg}}=0.16(R{e}_{D_H}{)}^{-1/8}$(30) where $I$ is the turbulence intensity, $U^{\mathrm{\prime }}$ is the root-mean-square of the velocity fluctuations, $U_{\mathrm{a}\mathrm{v}\mathrm{g}}$ is the mean flow velocity.

Figure 19(a) depicts 9 circular faces intercepted in the cavitation region along the flow direction, with each face spaced 0.5 mm apart. Figure 19(b) shows the distribution of turbulence intensity on the 9 surfaces at different pressures. It can be observed from Figure 19(b) that the overall trend of turbulence intensity increases with the rise of pressure. When the pressure is constant, the turbulence intensity tends to increase and then decrease from face 1 to face 9, reaching a maximum between 1 and 2 mm, and then gradually decreasing as the pressure rises, causing the degree of cavitation to decrease. Two critical points of cavitation length are 0.7 and 1.3 MPa in Figure 15. From Figure 19(b) it is evident that the flow can be classified into three categories based on turbulence intensity. When the pressure is below 0.7 MPa, the turbulence intensity is almost flat, corresponding to the initial cavitation stage in Figure 17. The turbulence intensity begins to show a clear trend of increasing and then decreasing when the pressure is between 0.7 and 1.3 MPa. The turbulence intensity at the tail of the flow (distance of 3.5–4 mm) shows a parallel trend at all pressures in this range, corresponding to the developing cavitation stage in Figure 17. When the pressure is greater than 1.3 MPa, the maximum turbulence intensity gradually shifts backward and appears in the region of 1.5–2 mm, while at the tail of the flow (distance of 3.5–4 mm), the turbulence intensity appears to rise abruptly (no longer parallel trend). The effect of turbulence inhibits the further development of cavitation because the flow becomes turbulent, which corresponds to the stable cavitation stage in Figure 17.

Figure 19. Comparison of the turbulence intensity of a fluid passing through the specified cross-sections at different pressures. (a) illustrates the position of the cross-section through which the fluid flows. (b) compares the turbulence intensity of the fluid at various pressures.

4.2. Effect of valve opening on the cavitation evolution

The valve opening is a critical parameter in regulating fluid velocity and pressure distribution in a pipeline. Decreased valve opening results in reduced fluid flow and an increase in fluid velocity, which in turn alters the pressure distribution. These conditions may create localised pressure reductions, potentially leading to liquid vaporisation. Conversely, fully opening the valve allows for maximum fluid circulation through the valve cross-section, resulting in more uniform fluid velocity and pressure distribution, reducing the risk of vaporisation due to a more even pressure and temperature distribution. However, the complex flow dynamics of the relief valve nozzle and valve disc area complicate the relationship between valve opening and cavitation. Therefore, a more in-depth study is necessary to fully comprehend the impact of valve opening on cavitation. The investigation for various valve openings is conducted using the $\mathrm{\Delta }\mathrm{P}$ of 1.5 MPa, while keeping the other parameter settings for the numerical simulations consistent with those presented above. Figure 20 displays the vapour phase contour distribution at valve openings between 0.5 and 1.5 mm. However, due to spatial constraints, only partial valve openings are presented, while the vapour phase contours for all valve openings can be found in Appendix B. This paper divides the cavitation region into two sections: region 1 and region 2, which are determined based on the location of cavitation occurrence. Region 1 (as shown in Figure 20(a1) marked by red) is formed between the valve disc and the valve seat, while region 2 (as shown in Figure 20(a3) marked by red) refers to the nozzle region. Analysis of Figure 20 and Appendix B indicates that as the valve opening increases, the cavitation length in region 1 gradually shortens and disappears entirely. Additionally, with valve opening increases, cavitation emerges in region 2. The cavitation length of region 2 increases gradually before stabilising.

Figure 20. The vapour phase contour distribution at different valve openings. (a1) 0.6 mm. (a2) 0.8 mm. (a3) 0.9 mm. (a4) 1.0 mm. (a5) 1.2 mm. (a6) 1.5 mm.

To further quantify the variation in cavitation length with valve opening and investigate the cause of the shift in cavitation location, Figure 21 plots the radial cavitation length of region 1 (as shown the black dimension line in Figure 20(a1)) and the axial cavitation size of region 2 (as shown the black dimension line in Figure 20(a3)) as a function of valve opening. As depicted in Figure 21, the cavitation length in region 1 gradually diminishes as the valve opening increases. Moreover, once the valve opening reaches 0.8 mm, cavitation becomes apparent in region 2. This phenomenon can be attributed to the fluid's through-flow area. With an increase in valve opening, the through-flow area expands, resulting in a decrease in flow velocity and a subsequent reduction in pressure drop. Consequently, the cavitation length diminishes. However, when the valve opening reaches 0.8 mm, the through-flow area in region 1 continues to enlarge, surpassing the through-flow area of the fluid passing through the nozzle. This leads to the emergence of cavitation in region 2. Considering the hindering effect of cavitation bubbles on the flow in the vicinity of region 1, when the valve opening is x_v, the actual effective circulation opening can be determined by subtracting the axial cavitation length y_length from x_v. Referring to the method of calculating the effective circulation area in Figure 2, the formula for calculating the effective circulation area in region 1 is obtained as shown in Equation (31). (31) $S_{region1}=\pi \cdot d_{avg}\cdot \left[x_v-y_{length}\right]$(31) Similarly, in the vicinity of region 2, the obstruction of flow by cavitation bubbles is represented by the radial cavitation length x_length. The actual effective circulation diameter can be determined by subtracting x_length from the nozzle diameter. Therefore, the formula for calculating the effective circulation area for region 2 is obtained as shown in Equation (32). (32) $S_{region2}=\frac{\pi }{4}\cdot (d_v-2\cdot x_{length}{)}^2$(32) Combining Equations (26) and (27), the formula for calculating the final effective area of region 1 and region 2 can be obtained as shown in Equation (33). (33) $\{\begin{array}{l}{S}_{region1}=\pi \cdot d_{avg}\\ \cdot \left[x_v-{\displaystyle \frac{V_{vapor-volume}}{{\pi }^2\cdot ({r_{vapor-end}}^2-{r_{vapor-start}}^2)}}\right]\\ S_{region2}={\displaystyle \frac{\pi }{4}}\cdot (d_v-2\cdot x_{length}{)}^2\end{array}$(33) Figure 22 demonstrates the relationship between the area of region 1 and region 2 and the valve openings calculated using Equation (33), which indicates that the cavitation length of region 1 decreases as valve opening increases, resulting in a quadratic trend for the through-flow area of region 1 at valve openings between 0.5 and 1.3 mm. Once valve opening exceeds 1.3 mm, cavitation disappears completely, causing the through-flow area to increase linearly with valve opening. Conversely, as valve opening increases, the cavitation size of region 2 increases, leading to a decrease in the through-flow area, which stabilises once the cavitation region reaches a fixed size. By comparing the through-flow areas of region 1 and region 2 (i.e. values of S1 and S2), the cavitation location pattern can be inferred. When S1 is greater than S2, cavitation is likely to occur in region 2 due to the relatively small through-flow area of this region.

Figure 21. The cavitation length as a function of valve opening in region 1 and region 2. The cavitation region of the nozzle and the gap (between the valve disc and the valve seat) are represented by region1 and region2, respectively. The variation curves of these two regions with respect to the valve disc opening are depicted in red and black, respectively.

Figure 22. The relationship between the area of region 1 and region 2 and the valve openings. The location of cavitation can be determined by comparing the through-flow areas. The through-flow areas of the gap region and nozzle region, represented by S1 and S2 respectively, are calculated using Equation 33. These areas are depicted by the black and red curves, respectively. Based on the characteristics of these curves, they are categorised into three stages, each corresponding to different positions of cavitation.

5. Conclusion

This study conducts a comprehensive examination of cavitation characteristics in a safety valve located in the secondary circuit of a nuclear power plant. To analyze the phenomenon, a high-fidelity cavitation model is proposed. It integrates the valve characteristic length, velocity strain rate correction mass transfer term, adaptive ensemble of surrogate models, and genetic algorithm. The model is then compared with experimental data obtained from a self-built visualisation test rig, and a remarkable agreement is observed. This research primarily focuses on two key control parameters of the relief valve: $\mathrm{\Delta }\mathrm{P}$ and valve opening While investigating the variations in cavitation, the influence of $\mathrm{\Delta }\mathrm{P}$ on cavitation is elucidated through an analysis of vortex dynamics and turbulence intensity. Additionally, a formula is introduced to calculate the effective flow area under cavitation conditions, elucidating the impact of valve opening on cavitation. The principal findings of this study are summarised as follows.

1. The integration of the valve flow characteristic length and the velocity strain rate-based correction method for the evaporation term in the Zwart-Gerber-Belamri (Z-G-B) model effectively mitigates the influence of variations in bubble diameter (D_b) and nucleation site volume fraction (α_nuc) on the cavitation length. This approach results in a significant reduction in the maximum error of D_b from 12.5% to 8.7% and limits the maximum error of α_nuc from 22% to within 10%. As a result, it achieves a substantial reduction in dimensional complexity and ensures physical accuracy.

2. The visualisation test bench built in this study allows clear observation of the cavitation phenomenon in the valve, and the cavitation length can be accurately measured using the gray scale detection method. Furthermore, the AES-MSI model accurately establishes the relationship between the evaporation/condensation coefficients and the cavitation length. Experimental verification confirms that the above coefficients, when optimised by genetic algorithm, enable accurate prediction of the cavitation length.

3. Pressure drop and valve opening are two critical parameters that influence the cavitation phenomenon. Concerning the impact of pressure drop on cavitation, this study concludes that the cavitation length gradually increases with pressure until it reaches a steady state. This pressure-dependent cavitation process can be divided into three stages: initial cavitation stage, developing cavitation stage, and stable cavitation stage. During the developing cavitation stage, a significant difference in cavitation length exists between the left and right sides of the valve disc due to the presence of vortices within the fluid. With a further increase in pressure, the transition to the stable cavitation stage occurs, leading to changes in the flow pattern. As a result, the flow resistance on the left side of the valve disc decreases, and the cavitation length becomes consistent. Additionally, turbulence intensity can serve as a reliable indicator to distinguish the three cavitation stages.

4. Concerning the impact of valve opening on cavitation, this study concludes that the cavitation region around the valve disc decreases with increasing opening, and simultaneously shifts towards the nozzle region. This change can be attributed to the relative alteration in the effective flow area of the aforementioned region. The potential mechanism behind this change in cavitation area can be elucidated by adopting the effective flow area calculation formula proposed in this study.

5.1. Limit of this study

Although this study constructs a cavitation model for a nuclear power second circuit safety valve and explores the cavitation evolution law, there is a limitation in the visualisation test bench. As a result, the pressure set in this study deviates from the actual operating conditions of the valve to some extent. In future studies, using visualisation materials with higher pressure-bearing capacity or adopting more advanced testing methods for cavitation phenomenon will be considered to overcome the current pressure limitation. Furthermore, the valve disc shape studied in this study belongs to the classical cylindrical type, featuring a relatively simple structure. In the subsequent research, a valve disc with a more complex shape will be investigated to further validate the applicability of the improved cavitation model proposed in this study.

Author statement

The research was conceptualised and guided by Xueguan Song and Qingye Li. Shuai Zhang was responsible for constructing the surrogate model, MuChen Wang conducted the visualisation experiments and handled the image processing. Chaoyong Zong and Xuyang Li undertook the CFD modelling, calculations, and data analysis.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (No. 52075068). We would also like to thank the editors and the reviewers for their constructive comments and helpful suggestions.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Funding

This work was supported by the National Natural Science Foundation of China (No. 52075068).

Appendices

Appendix A. Cavitation contours of relief valve at different $\mathrm{\Delta }\mathrm{P}$

Appendix B. Cavitation contours of relief valve at different valve openings