Impact of modelling assumptions in cavitating flow of simplified injector

Here, three-dimensional Reynolds-averaged Navier–Stokes (RANS) simulations are employed for in-nozzle flow modelling. The aim is to evaluate the capabilities and sensitivities of simulating in-nozzle flow phenomena by using homogeneous phase change models. Two commonly used homogeneous models are employed – (i) Homogeneous Equilibrium Model (HEM) and (ii) Homogeneous Relaxation Model (HRM) A simplified nozzle with reference experimental data is modelled in the developing, super cavitation, and hydraulic flip regimes using OpenFOAM. The influence of the inlet boundary condition and turbulence models are investigated. With HEM, three barotropic compressibility models are tested. The research seeks to understand the strengths and limitations of each approach by comparing different sub-models and their interactions with the flow fields and phase change phenomena. Ultimately, this study contributes to a better understanding of how different modelling choices can influence the accuracy and reliability of cavitation simulations using homogeneous phase change models. While the simulations demonstrated reliable predictions for low-order velocity statistics, substantial discrepancies were identified in the prediction of the cavitating region by both homogeneous models. Additionally, the relaxation model predicted only minor cavitation. In contrast, acceptable flow predictions were achieved using different barotropic compressibility models and inlet boundary condition types, as well as most RANS turbulence models, where the RNG 𝑘 - 𝜀 model deviated more from the other turbulence models. Therefore, the added value of the present study is summarised as: (a) Four different cavitation regimes are simulated with HEM and HRM in a single study. (b) The more comprehensive HRM approach does not produce better velocity or cavitation predictions for a low-speed flow using the default relaxation time. (c) The choice of turbulence model can radically affect the cavitation prediction, in addition to velocity fields. (d) The HEM flow simulation has little sensitivity to the barotropic compressibility model.


Introduction
Cavitation is a phenomenon of local phase change caused by decrease below the vapour pressure (Andriotis et al., 2008).This phenomenon occurs in various engineering applications such as propeller blade tips, pumps, hydrofoils, and injector nozzles.Based on previous literature, the primary mechanisms for cavitation formation in turbulent flows have been proposed to be as follows: (a) shear-driven cavitation, induced by sharp edges creating low-pressure regions that may promote vapour formation at nucleation sites (Franc and Michel, 2006) and (b) cavitating vortices or string-type cavitation, formed due to pressure drop in the vortex core (Schmidt and Corradini, 2001).
Moreover, turbulent cavitating flow is often a mutual interaction or reciprocal process, especially in the cavitation inception case.E.g. cavitation may influence turbulence by enforcing eddy break-up and, resolving the flow field to the level of the inertial sub-range (Koukouvinis et al., 2017).Hence, understanding and managing cavitation effects are crucial for optimising both atomisation efficiency and the durability of high-pressure liquid injectors in various applications.
The collapse of cavitation structures could cause a major pressure and temperature increase -at a scale of 300 MPa and 1000 K (Kyriazis et al., 2017), which makes surrounding structures susceptible to erosion (Franc and Michel, 2006).It is reported that cavitating vortices can be more erosive than spherical bubbles when they collapse (Ji et al., 2013).
Cavitating vortices have been observed in diesel injectors (Schmidt and Corradini, 2001) making the nozzles susceptible to erosion.The erosion could alter the intended spray pattern, impairing injector performance (Greif and Srinivasan, 2011).Furthermore, cavitation structures exiting the nozzle are one of the main causes for spray hole-tohole variations (Payri et al., 2004).Nevertheless, cavitation in diesel injector nozzles could also be advantageous.For instance, cavitation strings when released from an injector nozzle can improve atomisation by breaking up the spray into smaller ligaments and widening the spray angle (Hiroyasu, 1991;Afzal et al., 1999).
In modern internal combustion engines, injection pressures have increased to thousands of bars which has increased the local velocities in the nozzle to the upper limit of the subsonic flow (Duan et al., 2016).Furthermore, it has been reported that shock waves may exist at the spray edge (Hillamo et al., 2010).Obtaining insight into the internal flow of specific injectors often necessitates measurements, which can be challenging in the confined high-pressure spaces within injectors.Another approach to acquire this knowledge is through simulations.Nevertheless, simulating these flow conditions requires the use of specialised numerical techniques, such as compressible flow solvers and phase change models, along with turbulence models.Simulating cavitation requires the incorporation of appropriate models that capture the phase change dynamics and the effects of vapour formation and collapse on the flow field.The quantity of dissolved gas should also be taken into account since it can influence amount of void and its distribution in the domain (Battistoni et al., 2015).
There are two main approaches for modelling multiphase flows: (a) Eulerian-Lagrangian and (b) Eulerian-Eulerian methods.The present study uses the latter approach.In the scope of phase-change modelling, commonly used Eulerian-Eulerian methods are (in decreasing order of complexity): 1. Volume-of-Fluid (VoF), 2. Mixture models and/or bubble transport models, such as Schnerr and Sauer (2001) model or a thermodynamic closure model by Kyriazis et al. (2018), 3. Homogeneous or pseudo-fluid models, such as Homogeneous Equilibrium Model (HEM) by Clerc (2000) and Homogeneous Relaxation Model (HRM) by Downar-Zapolski et al. (1996).
The aforementioned models can be combined and extended to simulate phase-change phenomena involving multiple fluids using Eulerian-Lagrangian or Eulerian-Eulerian methodologies, employing diffused (i.e., mixture models) or captured interphase approaches (i.e., VoF).
In the context of injector simulations, the choice of modelling approach depends on the specific objectives of the study.If the primary focus is on capturing the intricate details of fuel-air mixing, breakup, and atomisation, the VoF method or mixture model may be more appropriate (Tretola et al., 2022).Yet the resolution demands of VoF are extraordinary due to the multi-scale nature of the bubbly flow.However, in flow situations where the overall flow behaviour, such as macroscopic characteristics of internal injector flow or global flow patterns, the mass transfer or homogeneous models may provide reasonable approximations with lower computational cost.Conversely, the mass-transfer models have notable limitations due to their reliance on empirical constants, which can introduce significant variability as previously highlighted by Ducoin et al. (2012) and Capurso et al. (2017).As pointed out by Capurso et al. (2017), mass-transfer models are commonly employed in scenarios where baroclinic torque plays a significant role, such as in hydrofoils and propellers.Hence, with relevance to the present investigations, homogeneous models may be preferred in the in-nozzle flow context, not only because they are free of empirical coefficients, but also because baroclinic torque is typically not a dominant factor in such applications (Egerer et al., 2014).
The present study is focused on homogeneous models within category 3 above, sometimes called pseudo-fluid models (Schmidt et al., 2010).This choice is influenced by two key factors: a. Demands on mesh resolution.To achieve meaningful and accurate interface behaviour using VoF method, the mesh resolution requirements can be substantial.This is underscored by the work of Arienti and Sussman (2014), emphasising the importance of a finely resolved mesh for capturing interface dynamics effectively.b.High Weber numbers.The flow regimes explored in this study are characterised by high Weber numbers.In such scenarios, the surface tension impact is typically low.Consequently, detailed interface resolution becomes less crucial.
These considerations inform the selection of homogeneous models in the current scope.More comprehensive descriptions of VoF and mixture models can be found in dedicated reference materials (Schnerr and Sauer, 2001) and handbooks (Franc and Michel, 2006;Ghiaasiaan, 2017).Homogeneous models require a continuous speed of sound relation between liquid and vapour phases, often referred to as a barotropic compressibility model.The most commonly known relations include the linear, Wallis (1969), and Chung (Chung et al., 2004) models.The Wallis model, originally employed for HEM by Schmidt et al. (1999), has been utilised by researchers such as Koukouvinis et al. (2017), Kyriazis et al. (2018) and Kyriazis (2018), among others.Nevertheless, for system simplicity, simulation convergence and stability, authors more often opt for the linear model (Kärrholm et al., 2007;Salvador et al., 2013Salvador et al., , 2018;;McGinn et al., 2021).Erney (2009) compared the aforementioned relations using OpenFOAM for a cavitating NACA hydrofoil, reporting no significant differences.Generally, the Chung model is relatively less frequently used.However, there is a lack of comparison of different compressibility models specifically regarding in-nozzle flow scenarios.
With Eulerian-Eulerian phase change models, turbulence effect on cavitation (and vice versa) has been studied numerically alongside various experimental investigations.For instance, studies employing particle image velocimetry (PIV), e.g.Iyer and Ceccio (2002), or direct numerical simulation (DNS) like (Okabayashi and Kajishima, 2009) have examined how cavitation growth and collapse influence shear flows.These investigations have revealed that the presence of cavitation bubbles can alter the shear layer's growth rate, modify downstream flow characteristics, and reduce cross-stream velocity fluctuations and Reynolds stresses.These changes indicate a decreased coupling between stream-wise and cross-stream velocity fluctuations affected by cavitation.Sou et al. (2014) conducted a well-documented experiment on shear-driven cavitation using a transparent square-shaped nozzle and Laser Doppler Velocimetry (LDV).Their study not only provided visualisations of cavitation regimes but also time-averaged velocity profiles at different axial positions from the nozzle edge.Additionally, CFD simulations were performed to compare to experimental data.The results indicated that conventional turbulence models struggled to accurately predict the inception of shear-driven cavitation.Koukouvinis et al. (2017) replicated the experimental setup using different cavitation models and various turbulence modelling approaches.Their findings suggested that properly executed Large Eddy Simulation (LES) outperformed Reynolds-averaged Navier-Stokes (RANS) models, particularly when capturing small eddies and nucleation sites with low-pressure drops was crucial.Edelbauer (2017) conducted simulations of all cavitation regimes in the same experimental setup using LES, validating a combination of Eulerian-Eulerian (internal flow) and Volume of Fluid (VoF) (external spray flow) approaches.Biçer and Sou (2016) successfully captured the transient motion of cavitation in the recirculation zone by employing the Schnerr-Sauer mixture model (Schnerr and Sauer, 2001), coupled with the Rayleigh-Plesset equation based on critical pressure (rather than vapour saturation pressure) proposed by Franc and Michel (2006).Their approach combined VoF with RANS turbulence modelling, specifically the Re-Normalisation Group (RNG) - model.
Despite using different cavitation models, the results in the above studies raise questions about the extent to which RANS can be reliably used for modelling cavitation.The present paper explores the choices for homogeneous phase-change models along with their sub-models and analyses the performance of each turbulence model described in the previous studies.Contrary to LES, RANS offers a computationally moderate solution that does not necessitate complex inlet boundary conditions, making it a valuable tool for specific cavitation simulations.Generally, based on the literature survey several research gaps emerge: (a) the need to systematically compare the performance of different cavitation regimes simulated with both HEM and HRM within a single study.(b) The recognition that the choice of turbulence model can significantly impact cavitation predictions, as well as velocity fields, underscores the need for further investigation into this area.(c) The existing literature lacks comprehensive comparative studies specifically evaluating the performance of different barotropic compressibility models in simulating in-nozzle flow scenarios using the HEM.
This study focuses on the rectangular nozzle in line with the experiments by Sou et al. (2014).Thirteen cases were simulated using both the homogeneous equilibrium and relaxation models representing the developing, super cavitation, and hydraulic flip regimes with fixed mass-flow rate, pressure inlet boundary conditions and different turbulence models.The main objectives of the study are: 1. determine the extent to which the HEM and HRM models can accurately predict the cavitation phenomena in RANS context 2. assess the effect of various turbulence models on cavitation 3. study the effect of different barotropic compressibility models 4. quantify the effect of inlet boundary condition type on flow characteristics

Numerical approach
Here, a brief description of the models are provided.The description covers cavitation and turbulence models, all of which have been employed in OpenFOAM (Weller et al., 1998).

Governing equations
As described in Section 1, several different types of cavitation models exist.The present study employs homogeneous models whose common assumption is that the liquid and vapour phases are considered to be in mechanical equilibrium (Clerc, 2000):   =   =  ;   =   = .This approach is also called the pseudo-fluid model by Moody (1965), and it utilises a typical set of governing equations for incompressible single-phase flow, i.e. mass conservation (1) and momentum conservation (2) accompanied with an equation of state (2) Despite the use of a similar system of conservation laws (1)-(2) as for single phase flow, both models still require closure through an equation of state.The key distinction lies in the specific equation of state employed by each model.

Homogeneous equilibrium model
The HEM assumes a uniform distribution of vapour and liquid phases within the system, assuming not only mechanical equilibrium but thermal equilibrium within phases (Clerc, 2000).This model simplifies the complex interactions between phases and neglects spatial variations (Dumont et al., 2001).While it may not capture the detailed dynamics and complexities of real-world systems, it provides a useful approximation for certain engineering applications and serves as a foundation for more advanced models, such as the HRM (Koukouvinis et al., 2017).
In general, for Eqs. ( 1)-( 2), HEM uses an equilibrium approach where the density depends on pressure and enthalpy, i.e.  , =  ( , , ℎ , ) by Bilicki and Kestin (1990).For most homogeneous multiphase flows, a barotropic fluid equation of state  , =  ( , ) associates the speed of sound with density as in (see Brennen (2005, p.220)) One of the main challenges of the barotropic fluid relation ( 3) is to assess the compressibility  or speed of sound of the two-phase mixture (0 <  < 1).At low void fractions, the speed of sound is primarily governed by the liquid phase, while at higher void fractions, the presence of gas phase and phase change phenomena contribute to the observed decrease and subsequent increase in the speed of sound, respectively (Wilson and Roy, 2008).Therefore, the speed of sound curve in homogeneous two-phase flow against void fraction exhibits a U-shaped form.The speed of sound or compressibility in a two-phase mixture can be predicted by the following models: • Linear model.Naive and the simplest model.Based on weighting the compressibility of each phase by the void fraction.
• Wallis model.Takes into account the equilibrium state of the phases and provides a more accurate prediction of compressibility than the linear model (Wallis, 1969).It uses a homogeneity assumption of the two-phase flow and is derived from the Wood (1930) equation for the speed of sound, which considers the saturation densities of the vapour and liquid phases and their respective compressibilities to calculate the overall compressibility of the mixture (5) • Chung model.Based on the hyperbolic two-fluid model for nonequilibrium and non-homogeneous flow.The hyperbolic two-fluid model uses the interfacial pressure jump terms in the momentum equations derived from the surface tension concept by Chung et al. (2004).The model promises better agreement with experimental results than the Wallis model for lower void fractions, i.e. for  ≪ 0.3. where The speed of sound for each barotropic compressibility model is illustrated in Fig. 1 for a water vapour-liquid mixture.After defining compressibility , the density  is defined, which can be coupled with the vapour volume fraction .The equilibrium quality or dryness is defined by Ghiaasiaan (2017, 103) as: which the same either for specific enthalpy ℎ, specific entropy , specific volume , or specific internal energy .ℎ   − ℎ   of Eq. ( 7) is equal to the latent heat of evaporation.
As stated above, for a single-component (pure vapour-liquid) mixture, thermodynamic equilibrium requires that   =  (Ghiaasiaan, 2017, p.104).Therefore, one can use the density to specific volume, denoted as  =  −1 .By substituting this relation into Eq.( 7) and performing elementary transformations, one can establish the relation for vapour volume fraction and one can obtain an equation of state to close the system of Eqs. ( 1)-( 2): where The equation of state ( 9) for each barotropic compressibility model is illustrated in Fig. 1.
The limitations of HEM are widely known.The model is unable to accurately capture strong kinetic or thermodynamic non-equilibrium effects that are present in certain flow conditions, such as annular flows or gas flows with droplets (Schmidt et al., 1999).In cases where non-equilibrium effects are relatively small, they can be addressed by incorporating correction terms, such as the drift flux velocity.However, when non-equilibrium effects become more significant, additional equations are required to accurately predict the system behaviour (Clerc, 2000), e.g. by introduction of momentum balance equations.In general, while the HEM has its advantages in certain engineering applications, it falls short in accurately capturing nonequilibrium effects which are important in many applications.To address these limitations, more advanced models like the HRM have been developed, providing a more realistic representation of phase change dynamics in cavitation and related phenomena.

Homogeneous relaxation model
Built upon the HEM, the HRM also assumes a homogeneous distribution of phases (Bilicki and Kestin, 1990).However, the HRM goes a step further by considering the dynamic evolution of phase change phenomena (Downar-Zapolski et al., 1996).It acknowledges that the transition between phases is not an instantaneous process but involves relaxation processes that gradually bring the system to an equilibrium state (Schmidt et al., 1999), which is approached by augmenting the system of Eqs. ( 1), ( 2) with conservation of energy.
However, the system of Eqs. ( 1), ( 2) and ( 10) is still not closed since in cases where non-equilibrium heat transfer dominates flow dynamics, conventional equations of state may not be suitable (Brennen, 2005).This is especially true for two-phase mixtures represented by the pseudo-fluid assumption, as they are not in thermodynamic equilibrium, and therefore, alternative approaches are required to account for the complex interactions between phases (Schmidt et al., 2010).
The HRM uses a relaxation approach to reach thermal equilibrium and therefore to close the system.It was originally introduced by Bilicki and Kestin (1990) for modelling one-dimensional flashing flow and reproduces the mass and heat exchange process between two phases.The relaxation itself is determined as which describes deviation of the quality  from equilibrium quality   (see Eq. ( 7)) over some timescale .The quality or vapour mass fraction  can be defined from the vapour volume fraction  (see Eq. ( 8)).
The timescale  in general has a form proposed by Downar-Zapolski et al. (1996), and the current study uses formulations and values suggested by Downar-Zapolski et al. (1996) for flash-boiling flows with pressures below 1 MPa Here  is a dimensionless pressure difference between the local static pressure and the vapour pressure, and  0 ,  and  are empirical coefficients 6.51 ⋅ 10 −7 [s], −0.257 and −2.24 accordingly.The relaxation Eq. ( 11) is used to obtain the pressure while mixture density and velocity are found respectively by continuity (1) and momentum (2) equation.The HRM solves the pressure to satisfy the chain rule and employs the continuity equation indirectly.Thus, the model does not couple the pressure using the continuity equation as is typically done for incompressible Navier-Stokes solvers.Through the chain rule, the pressure responds to both compressibility and density change due to phase change by Schmidt et al. (2010).In general, for thermodynamic non-equilibrium flows the density is the function of pressure, quality, and enthalpy by (Bilicki and Kestin, 1990): By subtracting the above Eq.( 13) from the conservation of mass (1) one can obtain a velocity divergence equation where the last term can be neglected under the isenthalpic process assumption, as made in the present study, where no heat exchange between the gas flow and the surrounding walls is considered.The Eq. ( 14) is furthermore coupled to the pressure and eventually reflected to the velocity field as A more detailed exposition of the HRM's governing equations and numerical implementation can be found in the original study by Schmidt et al. (2010).
Overall, the inclusion of relaxation processes and the consideration of the system's dynamic evolution make the HRM particularly wellsuited for simulating transient behaviour and non-equilibrium phenomena.Cavitation is inherently transient and involves rapid changes in pressure, vaporisation, and collapse of vapour cavities.The ability of HRM to simulate these time-dependent processes should provide a more realistic representation of phase change dynamics in cavitation and related phenomena (Schmidt et al., 2010).

Turbulence models
Cavitation involves the formation, growth, and collapse of vapourfilled cavities within a fluid flow, leading to significant changes in pressure, flow patterns, and phase transitions (Schmidt and Corradini, 2001).Turbulence plays a key role in these processes, and selecting an appropriate turbulence treatment is crucial for reliable predictions.The most fundamental approach to simulate turbulence is Direct Numerical Simulation (DNS).DNS is computationally heavy as it accurately resolves all turbulent scales, making it impractical for most industrial problems (Moin and Mahesh, 1998).The solution to this issue is to employ turbulence modelling at either sub-filter scales (LES) or over all turbulent scales (RANS).LES provides more detailed representation of the energy dissipation, taking into account larger eddies (Speziale, 1991;Sagaut, 2005).The current study uses the unsteady RANS (URANS) approach, particularly the following linear two-equation eddy viscosity models.Generally the - model kinetic energy equation and dissipation equation have the following form accordingly where Three models are utilised in the study: • the Re-Normalisation Group (RNG) -.The model is an enhanced version of the standard - model (where the turbulent kinetic energy and the turbulence dissipation,  and  accordingly, are modelled using transport equations with diffusion and source terms by Koukouvinis et al. (2017)).This model employs the Renormalisation Group (RNG) approach to improve the behaviour of the model near solid boundaries and to better predict flow separation and reattachment (Han and Reitz, 1995).This approach addresses the non-linearities in the model equations, which can help to improve the model's performance in flows with strong non-linearity according to Wilcox (2006).The model's source terms and coefficients are where  1 and  2 = 1.9,  3 = 1.62 and further closure coefficients can be found in the original study Yakhot and Smith (1992).• realisable -.The model addresses certain limitations of the standard - model, i.e. it reduces unphysically high production of turbulence (Shih et al., 1995).The model allows for more accurate predictions in flows with strong pressure gradients by separating the production and destruction of turbulence kinetic energy and dissipation rate (Pope, 2000).The model's source terms and coefficients are where ) and  2 = 1.9 and further closure coefficients can be found in the original study by Shih et al. (1995).
• - Shear Stress Transport (SST).The model is a hybrid model that combines aspects of the - and - models.The core idea is to use better predictions in the near-wall region from the former model and in the far-field for the latter by using blending functions (Menter, 2001;Pope, 2000).There are some changes from the ''parent'' models: turbulent viscosity is changed to incorporate turbulent shear stress and the empirical constants are tuned.Nevertheless, the model underpredicts normal turbulent stresses for low turbulence flows (Zhang et al., 2007).Instead of the dissipation equation, the model utilises a turbulence specific dissipation rate or turbulent frequency equation.

Discretisation methods
The simulations were performed using the open-source CFD software OpenFOAM (Weller et al., 1998).The effect of turbulence on the resolved flow field was modelled using an unsteady RANS (URANS) turbulence modelling approach.The modelling of phase change caused by cavitation was modelled using the homogeneous equilibrium model (HEM) by Clerc (2000) and homogeneous relaxation model (HRM) by Bilicki and Kestin (1990) and Schmidt (2023) with OpenFOAM 10 and v2206, respectively.
The temporal discretisation in the simulation was achieved using a second-order implicit scheme, while the diffusion terms were discretised using a second-order central scheme.The locally dissipative non-linear flux-limiting scheme developed by Jasak et al. (1999) was employed in the simulation for the advective terms.A control parameter value of 0.3 was used to discretise the momentum convection, while a value of unity was chosen for scalar convection to ensure a bounded Total Variation Diminishing (TVD) solution.A variable time-step algorithm was employed with a maximum Courant-Friedrichs-Lewy (CFL) number of 0.1 and 0.65 for HEM and HRM accordingly.In the HEM simulations, an additional limitation was imposed by restricting the acoustic CFL number to 8. Block-structured meshing was used (see in Fig. 2(b), (c)) with a resolution of 1.675 million cells and the smallest cell size of 44 μm.Detailed insights into the impact of different mesh resolutions on the simulation results can be found in Appendix A.1 where comprehensive mesh sensitivity tests are presented.Run times for developing cavitation cases were 466.9 and 131.4 h using 80 cores for the HEM and HRM accordingly.

Geometry specification and operating conditions
The simulations represent an experimental configuration by Sou et al. (2014).It is a scaled-up rectangular nozzle with a sharp edge on one side -the nozzle length and width is 8 and 1.94 mm accordingly (see in Fig. 2(a)).The orifice walls are parallel to each other (no profile divergence, i.e. -factor ≡ 0).The injected fluid is water at a temperature of 20 • C -fluid properties for the given temperature are represented in the Table 1.Pure, non-degassed water was assumed in the simulations due to the face that Sou et al. (2014) did not provide the non-condensable gas quantity dissolved in the ''filtered tap water'' used in their study.The properties in the Table 1 were used for the HEM simulations.For the HRM simulations tabulated, fluid properties were obtained from NIST REFPROP by Lemmon et al. (2018).
Three different cavitation regimes were studied.The regimes are generally defined using cavitation number defined as where  ,  and    mean liquid velocity (in the nozzle), local and vapour saturation pressure, respectively.The cavitation regimes are generally characterised as (see in Fig. 3(a)): i. Developing Cavitation (0.85 ≤  ≤ 1.2).In this regime, cavitation bubbles or strings predominantly appear in the upper half of the nozzle (Sou et al., 2007).ii.Super Cavitation (0.75 ≤  ≤ 0.85).The cavitation zone extends from the inlet to just above the nozzle exit (Sou et al., 2007).This regime is also often characterised by the presence of a cavitation film, which acts as a layer isolating the surface from the surrounding flow (Brennen, 2005).iii.Hydraulic Flip ( ≪ 0.75).The liquid flow separates at one inlet edge and does not reattach to the side wall (Sou et al., 2007).
More comprehensive cavitation regimes description can be found in handbooks (Franc and Michel, 2006;Brennen, 2005).Each mentioned regime was studied experimentally and reported by Sou et al. (2014).
In the current study, an investigation is conducted for each of these regimes, as depicted in Fig. 3(b).Three turbulence models were used which were tested using both equilibrium models and the homogeneous relaxation model.Moreover, three different types of inlet boundary conditions were investigated: (1) fixed volumetric flow rate (uniform velocity), ( 2) fixed uniform pressure, and (3) mapped (non-uniform) boundary conditions from developed channel flow profiles.The total simulation matrix is represented in Table 3.

Results and discussion
In the first subsection, the validation of a baseline setup is provided.The second subsection focuses on the effect of barotropic compressibility models for HEM simulations.In the third subsection, the effect of boundary conditions is discussed.The fourth subsection focuses on the analysis of homogeneous models and their impact.Finally, the influence of the selected turbulence model is assessed.Mesh sensitivity analysis of the baseline setup is provided in the Appendix.

Flow validation
The numerical methods are validated against the baseline condition (see Table 2) for both homogeneous models.Validation is performed by comparing three downstream profiles of mean velocity and root mean square (RMS) velocity fluctuation against both experimental (LDV by Sou et al.) and numerical results (RANS and LES by Koukouvinis et al. (2017) and Edelbauer (2017)).However, it is essential to acknowledge that there are significant differences between different numerical studies simulating the problem.It is observed that both homogeneous models provide satisfactory predictions for mean velocity in the upper part of the nozzle when compared to measurements but the results exhibit larger deviations nearer the nozzle exit.This study also observes discrepancies in the RMS prediction near the nozzle exit, similar to the findings reported by Koukouvinis et al. (2017) where they attributed those issues to the assumption that the -based RMS assumes isotropic velocity fluctuations ( ′  ′ ,  ′  ′ ,  ′  ′ ).The mean velocity and RMS of turbulent velocity fluctuations are shown in Fig. 4.
The cavitation presence may lead to transient flow features characterised by large-scale, periodic vortices (Arndt, 2002;Stanley et al.,   2014; Li et al., 2023).In essence, periodic cavitation may induce transient flow features that are resolvable on the grid.The flow velocity fluctuation phenomena are represented in the literature across various scenarios.The most well-known example is the Von Kármán (2004) vortex street.Notably, the simulation of Kármán vortices, characterised by a strongly periodic flow, does not require turbulence model utilisation.Nevertheless, the case's inherent flow dynamics can be accurately captured despite the absence of turbulence model (Beaudan, 1995;Murman, 2000), underscoring that not all flow fluctuations are inherently turbulent.In engineering contexts, similar phenomena may occur in the flow through hydroturbine blades (Zhang et al., 2018) and around hydrofoils (see e.g.Huang et al. (2014) and Ji et al. (2017)).
Our hypothesis regarding cavitation-induced fluctuations underscores the fact that velocity fluctuations in the non-cavitating regime are negligible compared to cavitating regimes.To properly assess those effects, when assessing the root mean square (RMS) fluctuations, one needs to consider not only the turbulent kinetic energy () but also the velocity fluctuations ( ′  ′ ), i.e.The RMS profiles exhibit more notable deviations, in contrast to the mean velocity profiles, from both the measurements and the simulations conducted by Sou et al. (2014) and Koukouvinis et al. (2017).The HEM simulations were able to improve the observed variations reproduced in the measurements by including velocity fluctuations (prime values) in the RMS calculations.It is important to note that the turbulent distribution obtained from the HRM model exhibits the correct pattern and matches better the Koukouvinis et al. (2017) RANS results.Despite the slight underestimation against the lowest measurement, the HRM model captures the general trend and characteristics of the turbulent fluctuations at all locations.Also for cross-comparison with other RANS results it is important to note that these analyses likely used the isotropic assumption.Therefore, taking anisotropy into account could potentially enhance predictive accuracy for those studies.

Effect of barotropic compressibility models
As described in Section 4.2, to adequately capture the liquid-vapour phase transition and intermediate states under the homogeneous phase distribution assumption, the pseudo-fluid models require continuous speed of sound relations.However, many researchers (e.g.Kärrholm et al. (2007), Salvador et al. (2013, 2018) and McGinn et al. (2021), etc.), for simulation convergence and stability, satisfy this requirement by employing the physically less realistic linear model.
Nonetheless, in this comparison, the linear, Wallis and Chung models are evaluated at the developing cavitation condition.Fig. 5(a,b,c) shows a time-averaged cavitation cloud comparison for the models.The results indicate that there are no significant differences among the three models in terms of predicting the formation and behaviour of the cavitation cloud and its penetration depth.
However, differences emerge in the lower void fractions areas ( ≪ 0.3), particularly evident in the lowest downstream profile, as depicted in Fig. 6.This aligns with Chung et al. (2004), who observed improved behaviour of their model for smaller void fractions.Notably, unlike other models, which predicted a void fraction of about 5 ⋅ 10 −7 , the Chung model predicted absolutely no vapour further from the nozzle centre-line.This qualitative alignment may correspond with the Sou et al. ( 2014) findings, who reported the bubble absence in that region.2).The profiles are selected based on the most substantial differences.Note that the Wallis model was utilised in the barotropic case of Koukouvinis et al. (2017).
Conversely, for RMS of turbulent fluctuations, the only differences are in the opposite region, which may also be connected to the fact that left from the centre-line, void fractions do not match either by maximum value or its location.However, again, the Chung model remains the most trustworthy since it is a low-vapour region.
As mentioned in Section 1, Erney (2009) conducted pressure coefficient comparison of a cavitating hydrofoil using the same barotropic models set.While minimal differences were observed, the study lacked low-fraction prediction in-depth analyses.In contrast to both our results and Erney's, where minimal differences between the models were observed, we employ the Chung model for further simulations within the study.This inclusion aims to address potential discrepancies for lower void fractions and enhance simulations robustness.

Effect of inlet boundary condition
The numerical solution sensitivity to the inlet boundary condition was studied at the developing cavitation regime.The variations included a uniform fixed volume flow rate, a uniform fixed total pressure, and mapped (non-uniform) profiles obtained from a precursor simulation.
Table 4 presents macro flow quantities, i.e. predicted pressure drop and flow rate through the domain.In the table, the error is defined against the reference values presented by Sou et al. (2014): 48 ml/s and 0.12 MPa accordingly -see Table 2.For fixed velocity cases (de.eqn.Ch.SST.U, de.eqn.Ch.SST.map) the pressure loss was overestimated, and for the fixed pressure case it was opposite of the  mass flow rate prediction.Thus, generally, the losses were slightly underestimated regardless of the selected boundary condition type.
Both downstream and stream-wise profiles (Figs. 8, 7 accordingly) showed a remarkable similarity with minor deviations observed for the simulation with pressure inlet boundary conditions.Notably, timeaveraged velocity and root-mean-square profiles (z = 3.6 mm in Fig. 8) exhibited increased deviations near the nozzle outlet.A similar trend was observed in stream-wise profiles in Fig. 7, with the velocity inflection point shifting downstream and turbulent kinetic energy presenting a broader, albeit lower, peak.The above observations indicate that the flow behaviour was moderately independent of the specific type of inlet boundary condition.
Overall, the findings suggest that the choice of the exact inlet boundary condition type has only a modest influence on the flow characteristics.This implies that the simulations are reasonably robust and the flow is relatively insensitive to the specific boundary condition details, allowing for a certain flexibility level in selecting the appropriate boundary condition for further simulations.Consequently, all subsequent simulations will employ either fixed pressure or fixed inlet boundary conditions.For instance, the former will be used for HEM simulations, while the latter will be applied for HRM simulations.

Effect of homogeneous models
Fig. 9 shows instantaneous qualitative high-speed imagery comparison from Biçer and Sou (2016) with a line-of-sight representation of CFD results.The representation is computed using the approach described by Rachakonda et al. (2018), which employs the Beer-Lambert intensity law with Mie regime assumption.The figure does not depict the no cavitation regime since both HEM and HRM (no.eqn.Ch.SST.U and no.relax.SST.pcase indices accordingly) successfully predicted the vapour absence for this particular regime.For the developing cavitation regime (see in Table 2), the HEM predicts the cavitation phenomenon accurately.However, the HRM produces hardly any vapour for this regime, as the minimum pressure fluctuates around the saturation pressure throughout the entire simulation.For higher cavitation regimes, both models tend to underestimate either the intensity or the penetration length of cavitation.Unlike the simulation by Biçer and Sou (2016), there is no attachment of the cavitation cloud to the walls for developing cavitation regime in our study, which aligns well with the experimental findings.As described in Section 3.2, non-degassed water was employed since the reference study by Sou et al. (2014) did not provide the non-condensable gas amount in their ''filtered tap water'', which could contain noncondensable gases up to a mass fraction of 3 ⋅ 10 −5 (see e.g.Battino et al. (1984)).Considering that by Battistoni et al. (2015) the noncondensable gas quantity can significantly affect the presence of voids, this discrepancy in dissolved gas could be a potential explanation of the observed cavitation underestimation in the simulations using both HEM and HRM.On the contrary, an alternative perspective from experimental observations could suggest that even a small amount of vapour at a specific point along the line of sight could effectively block a significant portion of light, potentially rendering other points along the same line of sight irrelevant.This threshold-based approach, while enhancing the both models representation, still under-predicted the cavitation cloud.However, due to the absence of a specific reference for determining the light blocking threshold, we have opted to employ the Beer-Lambert intensity law as a pragmatic solution.
In addition to under-predicting cavitation, the HRM also tends to underestimate the RMS (see Fig. 4), particularly in the recirculation zone near the nozzle outlet.While the HRM is a more complex model it encounters challenges in predicting cavitation inception for this particular geometry.It is important to note that the HRM was initially designed to simulate flash-boiling and incorporates an empirical timescale correlation (Eq.( 12)) to address thermal non-equilibrium effects between phases.Thermal non-equilibrium effects are significantly more pronounced in flash boiling flows than in cavitating flows, which may lead to the observed under-predictions for the cavitation inception in room temperature liquid.Nevertheless, the HRM is an established model which appears in the most common CFD codes (while it may not be as extensively utilised with Fluent and StarCCM+ compared to CONVERGE).Despite the lower robustness of the HEM for this particular geometry and regime, this model still has higher scientific interest since it includes an energy equation.This makes it possible to estimate heat exchange and non-equilibrium effects.

Effect of turbulence models
In this section, the realisable -, the RNG -, and the - SST models are compared for the developing cavitation regime.These specific models were chosen to combine the results obtained from the RNG - model, as conducted by Biçer and Sou (2016), with those from  the realisable - and - SST models, as conducted by Koukouvinis et al. (2017), within a single study using the same cavitation model.The HRM model is excluded from the comparison due to its underprediction of the cavitation inception (the only regime with reference experimental data).Therefore, the HEM model is employed here for the turbulence model comparison.
Each turbulence model matches well with the mean velocity reference data, with the RNG - model slightly outperforming others in the recirculation zone (see z = 6 mm in Fig. 12).The RMS results from the current study exhibit an improved match with the experimental pattern compared to the previous RANS study by Koukouvinis et al. (2017).
The more extensive prediction of cavitation from the realisable - model can be traced back to a specific feature within the flow field, notably the larger simulated recirculation region at z = 3 mm.This is clearly depicted in the mean velocity profiles presented in Figs.11 and 12.The expanded recirculation zone corresponds to an increased flow acceleration into the nozzle, causing a reduction in pressure and, consequently, an increase in cavitation (Kumar et al., 2020).This observation hints at the superior performance of the realisable model, especially within the z = 3 mm region, where it closely matches the experimental mean velocity profiles (Fig. 12).
Furthermore, one can observe a distinct pattern for the RNG - model in Fig. 10(I).The vapour generated by this model disperses immediately after the nozzle edge (see in Fig. 5,e), without forming larger cavitation bubbles downstream, contrary to the behaviour exhibited by the other two models and observed in the experiment (see in Fig. 9).The distinction in behaviour between the RNG and realisable models can be attributed to the treatment of the turbulent variable  which has a direct impact on turbulent viscosity   .The realisable model incorporates  as a variable while the RNG model uses it as a constant.A similar analysis by Kumar et al. (2020) has shown that the turbulence viscosity over-prediction may diffuse pressure gradients and constrict the recirculation region by decreasing pressure gradients which, therefore, reduce the extent of the cavitation prediction.
Despite the fact that the - SST model employs a different formulation for turbulent viscosity, it predicts a comparable level of turbulence with the realisable - model (see   in Fig. 11).Furthermore, the mean void fraction exhibits a similarity between these two models (see (c,d) in Fig. 5 and  in Fig. 11).This observation reinforces the trend: a decrease in turbulent viscosity corresponds to an increase in cavitation.
Therefore, the observations in the present paper support the observations by Kumar et al. (2020) that an over-prediction of turbulence viscosity, as seen in models such as RNG -, can have adverse effects on the cavitation prediction accuracy.In such cases, where excessive turbulent viscosity dampens cavitation, the choice of a turbulence model becomes critical.That is, models that tend to overestimate turbulence viscosity may not be the most suitable option for accurately capturing cavitation.Overall, each turbulence model demonstrated improvement in predicting the cavitation presence compared to the challenges reported in the RANS part of the study by Koukouvinis et al. (2017) (see Table 5).
In considering the developing cavitation region several observations emerge: 1.The - SST and the realisable - models predict at least one major cavitation string, as in Fig. 10(I), which is also observed on the line-of-sight views (Fig. 9).2. The - models demonstrate the attached string bubble formation at the nozzle edge, see Fig. 10(III, IV).In contrast, the SST model does not exhibit this cavitating structure, see Fig. 10(a).3.For the SST model, a gradual reduction in the cavitating string size is observed whereas for the realisable model, the cavitation string breaks up into smaller bubbles with a complex shape, see Fig. 10(II).4. The RNG model does not predict a major cavitation string, e.g. as depicted in Fig. 10(I).Instead, it exhibits a more dispersed bubble region -see in Fig. 10(IV) which is particularly noticeable in time-averaged profiles (Fig. 5(e)) and on the line-of-sight views (Fig. 9).

Conclusions
A rectangular nozzle was simulated reproducing the developing, super cavitation, and hydraulic flip regimes (Sou et al., 2014) using an Eulerian-Eulerian approach -homogeneous models (HEM and HRM), particularly.For the HEM, three different barotropic models were tested: linear, Wallis, and Chung barotropic compressibility models.Turbulence was modelled using the RANS approach and three twoequation models were used: the RNG -, the realisable -, and the - SST.The current results were compared not only to the experimental data but with other simulations, including simulations with more complex turbulence modelling approaches, i.e.LES.The comparisons were made by examining macro flow quantities, the averaged (along-, downstream of velocity and RMS profiles) and instantaneous parameters, such as qualitative line-of-sight views and iso-surfaces.Velocity fluctuations were incorporated in the result processing of averaged profiles.
The main findings of the comparisons are summarised according to the objectives list (see Section 1) as follows: 1.The more comprehensive HRM model did not yield better velocity or cavitation predictions in comparison to HEM in the present rectangular nozzle using RANS modelling.2. The choice of turbulence model affected cavitation behaviour.
This can be attributed to variations in formulations of the turbulent variable   which altered turbulence viscosity values.Specifically, the RNG - model yielded lower cavitation predictions.Such an observation aligns with the findings by Kumar et al. (2020), who also reported a similar dependency of turbulent viscosity on cavitation predictions.3. HEM shows satisfactory predictions when utilising any barotropic compressibility model, including even the simplest linear one.4. The inlet boundary condition type did not have a substantial impact on the predictions, although for each simulation the losses were slightly underpredicted, i.e. the domain pressure drop or volume flow rate difference was within 10% from the  reference.Therefore, the choice of the inlet boundary condition may not be a critical factor in predicting cavitation inception in the studied scenarios.
Overall, reasonable flow predictions were obtained with various RANS models and cavitation models, as well as the different barotropic compressibility models and inlet boundary condition types.The most significant disparities were observed in the cavitating region prediction.The present study demonstrated that various RANS models generally provided reasonable predictions of flow behaviour, with the exception of the RNG - model.The HRM faced challenges in reproducing experimental cavitation patterns in this low-speed test.However, in the higher velocity conditions typical of engine injector simulations, in scenarios where an equilibrium state, and particularly, where heat contributes to cavitation (i.e., flash boiling), HRM is expected to perform better.To date, most of the validation tests for HRM under cavitating conditions have been with much higher upstream pressure (see e.g.Duke et al. (2014Duke et al. ( , 2016Duke et al. ( , 2013) ) and Mitra et al. (2019)).It is possible that the low-speed flow of the present experiment is beyond the range of applicability of the HRM model.
For the present study, the HEM model was found to be preferable, except for the higher computational cost.This disadvantage does not extend to higher injection pressure situations where the HEM's compressible formulation becomes more efficient due to the system of equations becoming less stiff, leading to a less restrictive acoustic CFL number.Additionally, the present study found that the HEM's barotropic compressibility model and the inlet boundary conditions had a relatively minor impact on the results.The above observation potentially offers more flexibility and simplifies the modelling process for more complex realistic injectors.
However, some limitations arise from the simplified enlarged injector flow setup in the present study.Specifically, our findings may not directly extrapolate to scenarios involving low upstream pressure typical of many real-world engineering applications, such as hydrofoils and propellers.This limitation is twofold: (a) as described in Section 1, simulating such flows requires capturing baroclinic torque (see Capurso et al. (2017)), which is not feasible using the employed homogeneity assumption and (b) in these scenarios, numerous accompanying phenomena influence cavitation and turbulence dynamics, further complicating the predictive accuracy of the modelling assumptions applied in the present study.

Fig. 1 .
Fig. 1.Barotropic HEM equation of state and speed of sound variation for water at 20 • C using linear (olive), Wallis (purple) and Chung (pink) compressibility models.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 2 .
Fig. 2. Scheme of nozzle nominal geometry (a), finite-volume mesh of the domain overall (b) and the nozzle region (c).

Fig. 3 .
Fig. 3. Schematic representation of void fraction for developing cavitation (i), super cavitation (ii) and hydraulic flip (iii) regimes (a) based onSou et al. (2007).Approximate distribution of these regimes with reference to experimental(Sou et al., 2014), numerical(Koukouvinis et al., 2017) and present studies of the HRM with - SST RANS model and the HEM with - SST, RNG - and realisable - (b).

Fig. 4 .
Fig. 4. Time-averaged velocity and RMS of turbulent fluctuations profiles at developing cavitation regime using homogeneous equilibrium and relaxation models and - SST turbulence model.whereonly the lateral component is included, as the experimental measurements provided bySou et al. (2007) andSou et al. (2014) are provided for this specific component.Given that these fluctuations are distinct from turbulence, they are not incorporated within the turbulent kinetic energy () field.Therefore, including them into the RMS calculation does not lead to a duplication of the turbulent fluctuations.Moreover, the RMS fluctuations we are comparing against are experimental, disregarding whether resolved fluctuations are present or not.The RMS profiles exhibit more notable deviations, in contrast to the mean velocity profiles, from both the measurements and the simulations conducted bySou et al. (2014) andKoukouvinis et al. (2017).The HEM simulations were able to improve the observed variations reproduced in the measurements by including velocity fluctuations (prime values) in the RMS calculations.It is important to note that the turbulent distribution obtained from the HRM model exhibits the correct pattern and matches better theKoukouvinis et al. (2017) RANS results.Despite the slight underestimation against the lowest measurement, the HRM model captures the general trend and characteristics of the turbulent fluctuations at all locations.Also for cross-comparison with other RANS results it is important to note that these analyses likely used the isotropic assumption.Therefore, taking anisotropy into account could potentially enhance predictive accuracy for those studies.

Fig. 6 .
Fig. 6.Time-averaged vapour volume fraction and RMS of turbulent fluctuations profiles for linear, Wallis and Chung barotropic compressibility models at the developing cavitation regime (see in Table2).The profiles are selected based on the most substantial differences.Note that the Wallis model was utilised in the barotropic case ofKoukouvinis et al. (2017).

Fig. 7 .
Fig. 7. Time-averaged velocity and turbulent kinetic energy profiles along the nozzle at 22 μm from the nozzle edge in the developing cavitation regime using the homogeneous equilibrium model with fixed velocity, mapped profiles, and fixed pressure inlet boundary conditions.

Fig. 8 .
Fig. 8. Time-averaged velocity and RMS of turbulent fluctuations profiles at developing cavitation regime using homogeneous equilibrium model with fixed velocity, mapped profiles and fixed pressure inlet boundary conditions.

Fig. 11 .
Fig. 11.Time-averaged profiles along the nozzle at 47 μm from the nozzle edge at developing cavitation regime using the homogeneous equilibrium model with Chung barotropic compressibility models and - SST, RNG -, realisable - turbulence models.

Fig. 12 .
Fig. 12. Time-averaged velocity and RMS of turbulent fluctuations profiles in the developing cavitation regime using the homogeneous equilibrium model with - SST, RNG - and realisable - turbulence models.

Fig. A. 2 .
Fig. A.2. Time-averaged velocity and RMS of turbulent fluctuations profiles in the developing cavitation regime using the homogeneous equilibrium model with the - SST turbulence model for 59, 44 and 30 μm nominal cell size.

Fig
Fig. A.3.Time-averaged velocity and turbulent kinetic energy profiles along the nozzle at 0.97 and 0.72 mm from the nozzle edge in the no cavitation regime using the HEM and HRM.The profiles are selected based on the most substantial differences observed between the two codes.

Table 3
Simulation matrix (mesh sensitivity cases are not included).

Table 4
Estimated and predicted volumetric flow rate and pressure difference.

Table 5
Minimum domain pressure relative to saturation pressure (calculated as follows: []  = min () −   ) at developing cavitation regime.
Fig. 9. Instantaneous qualitative line-of-sight vapour volume fraction views for experiment (left column) and present simulation results in developing, super cavitation and hydraulic flip regimes (represented as rows according to Table2).Views have been selected to illustrate typical flow patterns and may not correspond to the same time-step.