Numerical Study of Perforated Plate Balanced Flowmeter Performance for Liquid Hydrogen

1. Introduction

In the production, storage, transportation, and trade of cryogenic fluids such as liquid nitrogen (LN₂) and liquid hydrogen (LH2), flow rate measurement is essential, which requires accurate, wide range, and reliable detection and maintenance. The orifice flowmeter realizes flow rate measurement according to the linear relationship between the pressure difference before and after the throttling orifice and the volume flow rate square. This kind of flowmeter has no moving parts, is simple, reliable, and low cost, and has been widely used. Compared with the traditional orifice flowmeter, the balanced flowmeter uses a perforated plate instead of a single orifice plate, which not only inherits the advantages of the latter but also can stabilize the flow field, reduce pressure loss, and delay cavitation occurrence, consequently, it shows a wider measuring range ratio and accuracy [1]. Therefore, it gradually replaces the latter to obtain more and more extensive applications.

Although the balanced flowmeter has been extensively studied with water as the working fluid [2,3,4]. For the application of the flowmeter to cryogenic fluids, the physical properties of cryogenic fluids differ greatly from that of water, as shown in Table 1 [5], and the temperature in the storage and transportation process of LN2 and LH2 is often close to the saturation temperature, so it is more prone to trigger cavitation [6]. The thermal effect of the cryogenic cavitation process cannot be ignored, thus, the LN2 and LH2 balanced flowmeter often presents more complex flow characteristics than that of water, requiring special optimization and calibration.

The application of a balanced flowmeter to cryogenic fluids is relatively late compared to the orifice type one. As early as 2006, Kelley et al. [7] invented a balanced plate for liquid oxygen measurement. They reported that compared with the orifice type, the upstream straight pipe length of the balanced one is required to be less than 0.5d, the pressure recovery rate is increased by 100%, the accuracy is improved by 10 times, and the noise intensity is reduced by 15 times [8]. In 2016, Liu et al. [9] studied the influence of fluid types on flow coefficient and static pressure loss coefficient through simulation. The results show that the lower limit of the Reynolds number of cryogenic fluids in the self-similar region is relatively close to that of water, but the cryogenic fluids have a higher upper limit of the Reynolds number. The pressure drop depends on the shape of the perforated plate and the physical properties of the fluid. Jin et al. [10] also carried out numerical research on the effects of hole distribution on measuring accuracy. They designed a plate with a central distribution of circular holes and pointed out that the plate with a slightly larger diameter central hole is more suitable for LH2 measurement than the plate with an equal hole diameter. In 2017, Shaaban et al. [11] designed a new balanced flowmeter with improved LH2 distribution by carrying out multi-dimensional and multi-objective optimization through numerical simulation. The results show that the flowmeter has improved flow coefficient, static pressure loss coefficient, and cavitation characteristics compared with the flowmeter designed by Jin et al. [10]. In 2018, Wang Jie [12] used LN2 as a working fluid to compare the performance of a balanced flowmeter under different inlet temperatures and outlet pressures by numerical method. The results show that a higher upper limit of Re number can be obtained in the self-similar region by reducing the inlet temperature or increasing the outlet pressure. While, for the experimental study of the cryogenic balanced flowmeter, few research results have been reported until now, except with the work by Tian [13] in 2016, who carried out an experimental study using LN2 as a working fluid.

To explore the flow characteristics of LH2 through the perforated plate in a balanced flowmeter, this paper established a numerical calculation framework based on the “single-phase flow” Mixture model, Realizable turbulence closure, and Schnerr- Sauer cavitation model. The thermal effect with cryogenic cavitation was considered. To evaluate the model, two balanced flowmeters and an LN2 flowrate experimental setup were designed and built to measure the flow coefficient and pressure loss coefficient. Then, based on the validated model, the flow process of the perforated plate under non-cavitation and cavitation conditions was investigated with water, LN2, and LH2 as the working fluids. The cavitation cloud, turbulence distribution, and temperature contours are present with the three fluids, and their effects on the measuring range are analyzed. The results can help to understand the mechanisms and characteristics of the LN₂ and LH₂ balanced flowmeter.

2. Numerical Model

Conservation Equations

When the pressure drops to the saturation value corresponding to the local temperature, cavitation will be triggered and the generated bubble cloud will detach from the wall and flow downstream of the perforated plate with the mainstream fluid. Simultaneously, the temperature will drop inside the cavitation regime. Thus, it is necessary to establish a mathematical model considering cavitation with thermal effects in the simulation. The calculation for cavitating flow in a perforated plate adopts the homogeneous model. It simplifies the "two-phase fluid" to "single-phase one" and ignores the slip velocity between the two phases, which means that the two phases share the same velocity. The inlet Re numbers for the working conditions are in the range of 10⁴~10⁶, which is a compressible turbulent flow. Thus, the conservation equations of continuity, momentum, energy, and vapor mass transport of the mixture are as follows:

1) Mass conservation qutation

$${\displaystyle \frac{\partial }{\partial \mathrm{t}}}\left({\rho }_m\right)+{\displaystyle \frac{\partial }{\partial x_j}}\left({\rho }_m{u}_j\right)=0$$

(1)

Here, u_j is the velocity and ρ_m and is the vapor fraction weighted average density of the mixture, which is calculated as and the same below:

$${\varphi }_m={\alpha }_v{\varphi }_v+{\alpha }_l{\varphi }_l$$

(2)

α is the phase fraction, and subscripts v and l are vapor and liquid phases, respectively.

2) Momentum equation

$${\displaystyle \frac{\partial }{\partial \mathrm{t}}}\left({\rho }_m{u}_i\right)+{\displaystyle \frac{\partial }{\partial x_j}}\left({\rho }_m{u_{i}u}_j\right)=-{\displaystyle \frac{\partial p}{\partial x_i}}+{\displaystyle \frac{\partial }{\partial x_j}}\left[\left({\mu }_m+{\mu }_t\right)\left({\displaystyle \frac{\partial u_i}{\partial x_j}}+{\displaystyle \frac{\partial u_j}{\partial x_i}}\right)\right]+{\rho }_m\overrightarrow{g}$$

(3)

Here，p is pressure, $\overrightarrow{g}$ is gravitational acceleration. μ_m and μ_t is the laminar and turbulent viscosity of the two-phase mixture, respectively.

3) Energy equation

$${\displaystyle \frac{\partial }{\partial \mathrm{t}}}\left({\rho }_m{h}_m\right)+{\displaystyle \frac{\partial }{\partial x_j}}\left[{(\rho }_m{h_m+p)u}_j\right]={\displaystyle \frac{\partial }{\partial x_j}}\left[\left({\lambda }_m+{\lambda }_t\right){\displaystyle \frac{\partial }{\partial x_j}}\left(T\right)\right]+S_E$$

(4)

Here, h_m is the enthalpy of the mixture, λ_m, and λ_t is the laminar and turbulent thermal conductivity, respectively. T is temperature and S_E is the volume heat source due to phase change:

S_E=-ht·R

(5)

Here, h_t is the latent heat of phase change, and R is the net mass transfer rate between the phases.

4) Vapor volume transport equation

$${\displaystyle \frac{\partial }{\partial \mathrm{t}}}\left({\alpha }_v{\rho }_v\right)+{\displaystyle \frac{\partial }{\partial x_j}}\left({{\alpha }_v{\rho }_{v}u}_j\right)=R_E-R_C$$

(6)

Here, R_E and R_C are the mass sources due to evaporation and condensation, respectively.

In addition, the ideal gas equation is used to calculate the vapor density. The other physical properties of fluids in the simulation, such as density ρ_l, dynamic viscosity μ, specific heat capacity c_p, thermal conductivity λ_m, and saturation pressure p_v are all functions of temperature, obtained from REFPROP 8.0 [5].

5) Turbulent closure

Due to the Reynolds stress in the Realizable k-ε model being more in line with mathematical constraints, consistent with turbulence mechanisms, and providing better descriptions of streamline and vortex characteristics, the turbulent model has been widely used to simulate cryogenic cavitating flow. The equations for turbulent kinetic energy k and turbulent kinetic energy dissipation rate ε are as follows.

$${\displaystyle \frac{\partial }{\partial \mathrm{t}}}\left({\rho }_{m}k\right)+{\displaystyle \frac{\partial }{\partial x_j}}\left({\rho }_m{ku}_j\right)={\displaystyle \frac{\partial }{\partial x_j}}\left[\left({\mu }_m+{\displaystyle \frac{{\mu }_t}{{\sigma }_k}}\right){\displaystyle \frac{\partial k}{\partial x_j}}\right]+G_k-{\rho }_m\epsilon +G_b-Y_M$$

(7)

$${\displaystyle \frac{\partial }{\partial \mathrm{t}}}\left({\rho }_m\epsilon \right)+{\displaystyle \frac{\partial }{\partial x_j}}\left({\rho }_m{\epsilon u}_j\right)={\displaystyle \frac{\partial }{\partial x_j}}\left[\left({\mu }_m+{\displaystyle \frac{{\mu }_t}{{\sigma }_{\epsilon }}}\right){\displaystyle \frac{\partial \epsilon }{\partial x_j}}\right]+{\rho }_m{C}_{1}S\epsilon -{\rho }_m{C}_2{\displaystyle \frac{{\epsilon }^2}{k+\sqrt{\nu \epsilon }}}+C_{1\epsilon }{\displaystyle \frac{\epsilon }{k}}{C}_{3\epsilon }{G}_b$$

(8)

Where, $\eta =S{\displaystyle \frac{k}{\epsilon }}$，$S=\sqrt{2S_{ij}{S}_{ij}}$，$S_{ij}={\displaystyle \frac{1}{2}}({\displaystyle \frac{\partial u_i}{\partial x_j}}+{\displaystyle \frac{\partial u_j}{\partial x_i}})$，G_k is the turbulent kinetic energy generated by the average velocity gradient, G_b is the turbulent energy generated by buoyancy, and Y_M is the total turbulent energy dissipation rate due to the compressible effect. The constant terms in the model are: σ_k=1.0, σ_ε=1.2, C_μ=0.09, C₁=max[0.43, η/(η+5)], C₂=1.9, C_1ε=1.44.

6) Cavitation model

The modeling of cryogenic cavitation needs to consider the thermodynamic effects. The Schnerr-Sauer model considers the thermodynamic effects and has proved that when the bubble number density is n_b=10⁸, the model can better predict the cavitating flow of cryogenic fluids [14]. The final form equation is [15]:

$${p\le p}_v\left(T\right),R_E=3{\displaystyle \frac{{\rho }_v{\rho }_l}{{\rho }_m}}{\alpha }_v\left(1-{\alpha }_v\right){\left({\displaystyle \frac{{\alpha }_v}{1-{\alpha }_v}}{\displaystyle \frac{3}{4\pi n_b}}\right)}^{-{\displaystyle \frac{1}{3}}}\sqrt{{\displaystyle \frac{2}{3}}{\displaystyle \frac{p_v\left(T\right)-p}{{\rho }_l}}}$$

(9)

$${p>p}_v\left(T\right),R_C=3{\displaystyle \frac{{\rho }_v{\rho }_l}{{\rho }_m}}{\alpha }_v\left(1-{\alpha }_v\right){\left({\displaystyle \frac{{\alpha }_v}{1-{\alpha }_v}}{\displaystyle \frac{3}{4\pi n_b}}\right)}^{-{\displaystyle \frac{1}{3}}}\sqrt{{\displaystyle \frac{2}{3}}{\displaystyle \frac{p-p_v\left(T\right)}{{\rho }_l}}}$$

(10)

3. Model Validation

The modeled geometry consists of a perforated plate and the upstream and downstream pipelines with an inner diameter of D=50 mm. Due to the symmetrical structure of the plate, to ensure both high efficiency and accuracy of the simulation, the model is simplified to a quarter of the total. Figure 1 shows the three-dimensional structure and cross-sectional view of the computational domain. To ensure the full development flow upstream and complete recovery of static pressure downstream, the lengths of the straight pipes upstream and downstream of the plate are set to 10D and 15D, respectively.

The equivalent diameter ratio β=0.1-0.75 is recommended, β=D_e/D, D_e is the equivalent diameter of the total opening area of the perforated plate. Correspondingly, the single hole size and distribution of perforated plates are preliminarily designed, as shown in Figure 2, which consists of a central hole and a circle of uniformly distributed small holes around it. Among them, D₁=D₂=10 mm and D_r=30.86 mm are used to represent the diameter of the central hole, surrounding holes, and distributed circle, respectively. t= 3 mm is the plate thickness.

The sub-zone grid scheme is adopted, that is, the grid inside the holes and the surroundings is specially densified, while the grid of the pipelines upstream and downstream are sparsed to save computing resources, as shown in Figure 3.

The inlet and outlet of the pipeline are set as velocity inlet and pressure outlet boundary conditions respectively. During the calculations, through fixed outlet pressure and adjustable inlet velocity, different inlet Re numbers.

The coupling between pressure and velocity adopts the "Coupled" algorithm. The discretization of the pressure and vapor phase volume fraction terms respectively adopts the "PRESSO!" and "QUICK" algorithms; The discretization of momentum, energy, and turbulence terms all adopts the second-order upwind algorithm. The convergence satisfies the following criteria simultaneously: (1) the residual values of the continuity and momentum equation are less than 10^-3, and that of other indicators is less than 10^-6; (2) The error of inlet and outlet volume flowrate is less than 0.1%.

To verify grid independence, grid numbers 1007340, 1519930, 2032520, and 2545110 are checked for calculation. The main difference between them is the number of nodes inside the holes of the plate. LN₂ was used as the working fluid with inlet velocities of 2.4m/s, 3.2m/s, 4.0m/s, 4.8m/s, and 5.6m/s, which includes non-cavitation and cavitation conditions. The flow coefficient C and static pressure loss coefficient ζ are obtained, as shown in Figure 4, where C is defined as the ratio of measured flow rate to theoretical value, and ζ =ΔP/(0.5ρu²), ΔP is the pressure difference between the tap. The changes in C and ζ for grid 2032520 and 2545110 are less than 0.3% and 0.5%, respectively. Therefore, to improve computational speed, the following calculations adopt a mesh scheme with a grid number of 2032520.

According to the geometric parameters in Figure 2 and concerning the international standard ISO 5167-2 [13], we designed and manufactured a DN50 balanced flowmeter, as shown in Figure 5. the parameters are listed in Lable 2.

Table 2. Geometrical parameters for the test perforated plate.

Type	Inner diameter D (mm)	β	Thickness t (mm)	Number of surrounding holes N	D1/D2	Dr (mm)
Orifice	50	0.6	3	0	\	\
Perforated				8	1	32.26

Table 2. Geometrical parameters for the test perforated plate.

Type	Inner diameter D (mm)	β	Thickness t (mm)	Number of surrounding holes N	D1/D2	Dr (mm)
Orifice	50	0.6	3	0	\	\
Perforated				8	1	32.26

Due to the scarcity of experimental data for LH₂ balanced flowmeters, we constructed a DN50 LN₂ experimental rig for flowrate measurement based on the standard flowmeter, as shown in Figure 6. The main components include: an injection Dewar with a volume of 300L and rated pressure of 0.6MPa, a standard LN₂ flowmeter (Hoffer-hfc2000, Accuracy: ± 0.2% FS), a vacuum chamber with the test flowmeter inside, a cryogenic control valve, and a collection Dewar. The LN₂ is pressurized out of the Dewar by the high-pressure nitrogen gas, and due to the vertical flow of LN₂ from bottom to top, the single-phase state and measurement accuracy of the standard cryogenic flowmeter are ensured. The temperature is measured by a PT100 thermometer with an accuracy of ± 0.05k.

Figure 7 shows the numerical and experimental results of the flow coefficient C and static pressure loss coefficient ζ of the orifice and perforated plate with different Re. Compared with orifice plate. It is found that the average flow coefficient C of the perforated plate increased by 16.1%, and the average static pressure loss coefficient ζ decreased by 9.8%. Thus, the preferred porous plate can be better used for the measurement of LN₂ flow. The quantitative deviation of numerical results compared to the experimental one is about 5%, which may be caused by the fact that the working conditions in the actual measurement are difficult to be completely consistent with those in the numerical simulation. Overall, the simulation results are consistent with the experimental results qualitatively, and the accuracy of the numerical model is verified.

4. Modeling and Analysis of LH₂ Balanced Flowmeter

Taking the above developed perforated plate as the research object, water, LN_2, and LH₂ are selected and modeled for comparison with the inlet temperature is 300K, 77.36K, and 20.37K, respectively, and the out pressure is all set as 0.2MPa. The results are shown in Figure 8 and Figure 9. It can be seen that compared with water, the variation range of flow coefficient C and static pressure loss coefficient ζ of LN₂ before the onset of cavitation is similar and smaller, the average C value decreases by 1.0%, the average ζ value increases by 1.9%, while, the upper limit of Re number Re_U for LN₂ fluid increases from 2.99 × 10⁵ to 10.65 × 10⁵. Compared with LN₂, the variation of the C and ζ of LH₂ under non-cavitating conditions are relatively consistent, the average C and ζ are relatively close, while the Re_U increases show a significant improvement from 10.65 × 10⁵ to 37.82 × 10⁵. This means that the Re range of LH₂ is greatly improved compared to LN₂. In addition, when Re<10.65×10⁵, due to the almost identical C and ζ values of LH₂ and LN₂, the flowmeter can be directly used for the measurement of two fluids. The C value calibrated with LN₂ has an error of 0.4% when directly used for LH₂, which can be eliminated by correction.

Figure 10 shows the axial velocity and turbulence intensity contours of the three fluids under non-coitation conditions. It can be seen that compared with water, at the same inlet velocity, the velocity of LN₂ and LH₂ in the hole is slightly higher, the area of the vortex zone at the hole inlet and downstream of the perforated plate is slightly larger, and the turbulence intensity at the two places is also slightly higher. In addition, Figure 11 shows the cavitation and turbulence intensity contours of the three fluids under the cavitation condition. It indicates that the LN₂ is more sensitive to cavitation than the other two fluids, and the turbulence intensity downstream of the perforated plate is also high. However, as the inlet velocity continues to increase, the pressure drop of all fluids through the plate, as well as the cavitation intensity, also continues to increase. Figure 12 shows the axial cavitation cloud and temperature contours of LH₂, which show that the cavitation process of LN₂ and LH₂ is accompanied by a decrease in temperature. Therefore, compared with water, the LN₂ and LH₂ flow slightly increase the losses at the entrance of the hole and downstream of the perforated plate. In the non-cavitation condition, the flow field characteristics of LN₂ and LH₂ are relatively consistent, but in the cavitation condition, the onset of cavitation of LN₂ occurs earlier. The cavitation vapor in the orifice reduces the effective flow area and, therefore, leads to higher static pressure loss.

5. Conclusions

The numerical study of a cryogenic balanced flowmeter with LN₂ and LH₂ as working fluid was carried out. The computational framework was based on the Mixture "single-phase flow" model, the Realizable turbulence model, and the Schnerr Sauer cavitation model. The model was verified by comparing it with the experimental data of the LN₂ flow rate. Then, the flow field characteristics of the perforated plate under non-cavitation and cavitation conditions are simulated and compared between the results of water, LN_2, and LH₂. Also, the measurement performance of the orifice plate and perforated plate was compared. The following conclusions were obtained:

1. In non-cavitation conditions, the pressure drop coefficient hardly changes with Re; Under cavitation conditions, cavitation first occurs at the entrance of the hole, which leads to a decrease in the effective flow area of the liquid and an increase in static pressure loss. Compared with the orifice plate, the perforated plate significantly increases the Re upper limit of the non-cavitation zone, thus significantly expanding the measurement range.

2. Different working fluids affect measurement stability and Re upper limit of the balanced flowmeter. Compared with water, LN₂ can obtain a relatively stable flow coefficient C, but cavitation occurs earlier and the measurement upper limit is lower compared to LH₂. The lower density of LH₂ results in a lower pressure drop and higher measurement limit.