A method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework

Abstract

This invention provides a method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework. Step S1 involves constructing a Schnerr-Sauer three-phase cavitation model based on a multiphase flow interface capture method and a two-phase Schnerr-Sauer cavitation model. Step S2 involves constructing a base class for the Schnerr-Sauer three-phase cavitation model. Step S3 involves constructing the three-phase phase volume fraction transport equations as the main solution program based on the two-phase Schnerr-Sauer cavitation solver. This method can construct a solver suitable for three-phase cavitation problems, filling the gap in existing technology for handling three-phase cavitation problems involving water, vapor, and non-condensable air phases, and providing a convenient and efficient numerical solution method for related engineering fields.

Description

Technical Field

[0001] This invention belongs to the technical field of fluid mechanics, specifically relating to a method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework. Background Technology

[0002] Cavitation is a phenomenon that occurs when the static pressure in a localized area of ​​a liquid drops below a critical pressure threshold, causing an evaporative phase change from liquid to gas. This phenomenon is widespread in various engineering applications and can negatively impact structural or machine performance, including causing noise, vibration, and surface erosion of materials.

[0003] The mainstream approach to solving cavitation problems in engineering currently relies on computational fluid dynamics (CFD) technology to construct corresponding mathematical models. These models are then developed using commercial software (such as Fluent) or open-source software (such as OpenFOAM), combined with cavitation models (such as the Schnerr-Sauer model) to construct mass transfer terms, thereby achieving accurate simulation of the cavitation process. However, current commercial and open-source software has limitations: it cannot add a three-phase cavitation solver to the existing aqueous and vapor phases (non-condensable air phase). Therefore, it cannot provide a convenient numerical solution for three-phase cavitation problems, including those in spillways and hydrofoils.

[0004] Although existing technologies have proposed several methods for constructing Schnerr-Sauer three-phase cavitation solvers, none of these methods specifically design and optimize the flux-limited calculation method for the unique condition of three-phase systems. This makes it difficult to maintain high accuracy of the discretization scheme while ensuring that the calculation results do not exceed limits (i.e., maintain physical rationality) during the simulation process. Therefore, in practical applications, most can only use a first-order precision discretization scheme to guarantee the physical meaning of the solution, but this also means sacrificing computational accuracy, leading to higher numerical dissipation problems and affecting the accuracy and reliability of the simulation results.

[0005] Therefore, there is a need to provide an improved technical solution that addresses the shortcomings of the existing technology. Summary of the Invention

[0006] The purpose of this invention is to address the lack of three-phase cavitation solvers currently available for engineering problems by constructing a highly customizable three-phase cavitation solver with high degrees of freedom. This allows users to add appropriate cavitation models according to their own requirements, thereby enabling accurate prediction of three-phase cavitation phenomena with relatively low computational cost and high precision.

[0007] To achieve the above objectives, the present invention provides the following technical solution:

[0008] A method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework includes the following steps:

[0009] Step S1: Construct a Schnerr-Sauer three-phase cavitation model based on the multiphase flow interface capture method and the two-phase Schnerr-Sauer cavitation model;

[0010] Step S2: Construct a base class for the Schnerr-Sauer three-phase cavitation model;

[0011] Step S3: Construct the three-phase phase volume fraction transport equation as the main solution program based on the two-phase Schnerr-Sauer cavitation solver.

[0012] The process of constructing the volume fraction transport equations for the three phases includes:

[0013] Add a mixture model to the header file used by the two-phase cavitation solver to initialize the field, and use the mixture model to construct the three-phase volume fraction field and density constant;

[0014] Obtain the parameters of the mixture model from the header file of the two-phase cavitation solver used to calculate the phase volume fraction transport equation;

[0015] The three-phase fluxes of the phase volume fraction transport equation are constructed. After the three-phase fluxes are constructed, flux limits are imposed on each phase flux based on global and local extrema. The maximum absolute value of the allowable inflow and outflow anti-diffusion fluxes of each unit is calculated.

[0016] The reduction factor λ of each unit limiter is calculated based on the calculated maximum anti-diffusion flux. The outflow and inflow flux limiters of each unit limiter are calculated by iterative method to obtain the iterative formula. The reduction factor λ of each unit limiter is discretized to each interface, and then the flux of each interface is limited so that the sum of the corrected fluxes of each phase at each interface is 0 and the sum of the volume fractions inside the unit limiter is 1.

[0017] The base flux is added to the modified flux to form the complete flux, which is then solved explicitly using the explicitSolve function of the MULES algorithm, thereby constructing the three-phase volume fraction transport equation.

[0018] Step S4: Construct a Schnerr-Sauer three-phase cavitation derived class model based on the two-phase Schnerr-Sauer cavitation model.

[0019] Step S5: Compile the source files of the Schnerr-Sauer three-phase cavitation derived class model and the corresponding base class into a dynamic library, compile the volume fraction transport equations of the three phases and link them with the dynamic library to build a new Schnerr-Sauer three-phase cavitation solver.

[0020] The usability of the new Schnerr-Sauer three-phase cavitation solver and the effectiveness of model corrections were tested based on native two-dimensional example meshes.

[0021] Preferably, step S1 includes:

[0022] Step S11: Construct the governing equations for the three-phase cavitation of water, steam, and air;

[0023] Step S12: Construct a three-phase cavitation model based on the two-phase Schnerr-Sauer model to provide a method for calculating the mass transfer term;

[0024] S13. Select the discretization method and construct the pressure-velocity decoupling algorithm according to the problem to be solved.

[0025] Preferably, step S2 includes:

[0026] Step S21: Extend the base class of the two-phase cavitation mixture to the base class of the three-phase cavitation mixture;

[0027] Step S22: Make the base class of the three-phase cavitation mixture inherit the interface feature module and the basic mixture model of OpenFOAM;

[0028] The constructor was modified synchronously, and the original correct() function was divided into viscosity coefficient correction and interfacial tension correction.

[0029] Add the source term function vDotAlphav for the vapor phase volume fraction transport equation and perform the calculation.

[0030] Preferably, in step S3,

[0031] When limiting flux, first calculate the maximum absolute value of the allowed inflow and outflow of the back diffusion flux for each unit limiter.

[0032] Preferably, the iterative formula for the reduction factor λ of the unit limiter is as follows, where λ 0 =1:

[0033]

[0034] By discretizing the reduction coefficient λ of the unit limiter onto each interface, the flux at each interface is limited.

[0035] Preferably, the sum of the corrected fluxes on each facet is calculated:

[0036] sum=ΣΦ i

[0037] Calculate its residual R:

[0038] R = sum - 0

[0039] Based on the sign and magnitude of the residuals, the inflow correction flux or outflow correction flux at the interface is reduced.

[0040] If R is positive, the outflow correction flux is reduced; if R is negative, the inflow correction flux is reduced.

[0041] This results in the sum of the modified fluxes of each phase at each interface being 0, and the sum of the volume fractions inside the cell confinement being 1.

[0042] Preferably, in step S4, the new Schnerr-Sauer three-phase cavitation derived class model is constructed using the following method:

[0043] First, the rRb function in the Schnerr-Sauer source code was rewritten accordingly;

[0044] Then the mDotAlphal function used for the source terms of the phase volume transport equation and the mDotP function used for the source terms of the pressure equation are rewritten respectively.

[0045] This completes the construction of the cavitation model.

[0046] Preferably, in step S11, the continuity equation in the governing equations for the three-phase cavitation of water, steam, and air is:

[0047]

[0048] According to the assumptions of the VoF model, ρ m The average density is obtained by weighting the volume fractions of each phase within a control body.

[0049] In step S12, since the surface tension between air and water vapor is too small, only the surface tension between water and the other two phases is considered.

[0050] Preferably, in step S13, the equations are closed using the Schnerr-Sauer three-phase cavitation model;

[0051] The phase transition coefficient is adjusted based on the turbulence timescale and corrected using strain rate according to the turbulence model assumptions.

[0052] After the model is built, the finite volume method (FVM) is selected to discretize the equations, and a co-located grid is selected to store the physical quantities.

[0053] Beneficial effects: This invention can construct a solver suitable for three-phase cavitation problems, filling the gap in the existing technology for handling three-phase cavitation problems including water phase, vapor phase and non-condensable air phase, and providing a convenient and efficient numerical solution method for related engineering fields.

[0054] By specifically constructing and limiting the three-phase flux, this invention can maintain the high precision of the discrete scheme while ensuring that the calculation results do not exceed the limits. This avoids the numerical dissipation problem caused by the use of first-order precision discrete schemes in the prior art, and improves the accuracy and reliability of the simulation results. Attached Figure Description

[0055] Figure 1 Flowchart of a method for constructing a three-phase cavitation solver based on the OpenFOAM framework;

[0056] Figure 2 The flowchart shows the calculation process of the three-phase cavitation model using the PIMPLE algorithm.

[0057] Figure 3 This is a UML framework diagram of the cavitation model framework in this invention;

[0058] Figure 4 Here is the mesh diagram for the example;

[0059] Figure 5 This is a diagram showing the pressure distribution along the central axis of the throttle valve channel after the two-phase air conditioning condition has stabilized.

[0060] Figure 6 The cavitation number distribution along the upper wall of the throttle valve channel after stabilization under two-phase non-cavitation conditions;

[0061] Figure 7 This is a diagram showing the distribution of cavitation cavities after stabilization under three-phase cavitation conditions.

[0062] Figure 8 This represents the pressure distribution along the central axis of the first half of the channel after stabilization under three-phase air conditioning conditions. Detailed Implementation

[0063] The present invention will now be described in detail with reference to the accompanying drawings and embodiments. It should be noted that, unless otherwise specified, the embodiments and features described herein can be combined with each other.

[0064] OpenFOAM is an open-source computational fluid dynamics (CFD) software with completely open source code, allowing users to freely view, modify, and extend it. This makes it a flexible and powerful tool suitable for various complex simulation needs, and it has wide applications in multiple simulation analysis fields such as aircraft aerodynamic design, automotive aerodynamics, and climate change simulation. This application presents a method for constructing a three-phase cavitation solver based on the Schnerr-Sauer cavitation model within the OpenFOAM framework, used for high-precision numerical simulation of cavitation phenomena in liquid-vapor-solid / gas multiphase systems. This addresses the problems of existing technologies being unable to effectively handle three-phase cavitation problems and having low computational accuracy due to flux limitations.

[0065] Specifically, such as Figure 1 The flowchart of the method for constructing a three-phase cavitation solver is shown, mainly including five steps. Step S1: This application constructs a Schnerr-Sauer three-phase cavitation model based on the multiphase flow interface capture method and the two-phase Schnerr-Sauer cavitation model. The multiphase flow interface capture method is a numerical simulation technique based on volume fraction tracking of multiphase flow interfaces, widely used in the interface evolution analysis of gas-liquid and liquid-liquid multiphase flow systems. The two-phase Schnerr-Sauer cavitation model is a high-order solver specifically designed to handle two-phase flows accompanied by phase changes. By coupling the VOF interface capture method with the phase change physics model, it becomes an industrial-grade tool for simulating cavitation, boiling, and condensation phenomena. Following the construction of the three-phase cavitation model, it is also necessary to select the discretization method and construct a pressure-velocity decoupling algorithm based on the problem being solved.

[0066] Step S2 involves constructing a base class for the three-phase cavitation model based on library files and three-phase flow and interface processing classes. This base class is responsible for calculating the physical properties, phase transition processes, and mixture properties of the three-phase system (water, steam, and air). The library files contain definitions of containers, fields, meshes, and mesh features used in various computational operations. These libraries provide OpenFOAM with rich functionality, such as finite volume discretization, boundary condition processing, turbulence models, and transport models, enabling users to flexibly construct and solve complex continuum mechanics problems. The three-phase flow and interface processing classes primarily focus on mesh management, field operations, matrix solving, and boundary condition processing. In this step, the constructed Schnerr-Sauer three-phase cavitation model forms the three-phase cavitation base class through inheritance and other methods.

[0067] Step S3: Construct the three-phase phase volume fraction transport equation as the main solution program based on the two-phase Schnerr-Sauer cavitation solver.

[0068] The construction of the three-phase volume fraction transport equations includes:

[0069] Add a mixture model to the header file used for initializing the field in the two-phase cavitation solver, and use the New function to point to a subclass object using a dynamic selection mechanism. During later runtime, the SchnerrSauer class can be selected using TransportProperties, and the mixture model can be used to construct the three-phase volume fraction fields and density constants. Obtain the parameters of the mixture model from the header file used for calculating the phase volume fraction transport equations in the two-phase cavitation solver, and then construct the fluxes for each phase volume fraction transport equation. Since three-phase calculations are involved, the three-phase fluxes for the phase volume fraction transport equations are constructed first, including the construction of flux and artificial compression flux, where artificial compression flux is only added between water-vapor and water-air.

[0070] After constructing the three-phase fluxes, flux limits for each phase are imposed based on global and local extrema. The maximum absolute value of the allowable inflow and outflow anti-diffusion fluxes for each unit is calculated. Based on the calculated maximum anti-diffusion flux, the reduction factor λ of each unit's limiter is calculated. The flux limiters for outflow and inflow of each unit's limiter are calculated using an iterative method to obtain the iterative formula. The reduction factor λ of each unit's limiter is discretized onto each interface, thereby limiting the flux at each interface so that the sum of the corrected fluxes of each phase at each interface is 0 and the sum of the volume fractions within the unit limiter is 1. Finally, the base flux is added to the corrected flux to form the complete flux, which is then explicitly solved using the explicitSolve function of the MULES algorithm to construct the three-phase volume fraction transport equation. To ensure that each phase does not exceed the limits while the sum of the three-phase volume fractions within each unit remains 1, further adjustments to the corrected fluxes are required. For each interface, the sum of the corrected fluxes of each phase should be 0 in order to satisfy the physical definition that the sum of the volume fraction fluxes of each phase on this surface is 1, and at the same time, it will also ensure that the sum of the volume fractions within the cell is 1.

[0071] This allows us to avoid obtaining non-physical results (i.e., results greater than 0 or less than 1) while selecting a high-precision discretization scheme during the solution of the phase volume fraction transport equation.

[0072] Step S4: Construct the Schnerr-Sauer three-phase cavitation derived class model based on the Schnerr-Sauer class. Specifically, this involves reconstructing the functions for the source terms of the phase volume transport equation and the functions for the source terms of the pressure equation, thus completing the construction of the cavitation model.

[0073] Step S5: Compile the source files of the Schnerr-Sauer three-phase cavitation derived class model and its corresponding base class into a dynamic library. Compile the volume fraction transport equations of the three phases and link them with the dynamic library to construct a new Schnerr-Sauer three-phase cavitation solver. To test the solver's performance, based on the original two-dimensional example mesh, the mesh was refined in three layers on both sides of the throttling channel to ensure the accuracy of physical quantity calculations in the high gradient region. The calculation results were compared with those of the new Schnerr-Sauer three-phase cavitation solver under two-phase cavitation conditions and no cavitation conditions. Furthermore, under the three-phase cavitation condition, the effectiveness of the phase transition coefficient correction was tested to test the usability of the new Schnerr-Sauer three-phase cavitation solver and the effectiveness of the model correction.

[0074] like Figure 2 The calculation flow of the three-phase cavitation model shown is based on the PIMPLE algorithm, a pressure-velocity coupled solution algorithm in OpenFOAM used for transient simulations. The specific calculation flow includes: 1. At the beginning of each time step, read or initialize the velocity field (U), pressure field (p), and volume fraction field (alpha) of each phase, then initialize the velocity field, pressure field, and volume fraction field within the three-phase cavitation model. 2. Calculate the surface flux (phi) using the current velocity field, then update the flux to ensure that the velocity and pressure fields satisfy the continuity equation. 3. Solve the phase volume fraction equation to update the distribution of each phase in the mesh, obtaining a new phase fraction field, which in turn affects the fluid properties (such as density and viscosity). 4. Calculate the density and viscosity fields of the mixture based on the updated phase volume fractions. 5. Solve the momentum equation (ignoring the pressure gradient term) to obtain the predicted velocity field for momentum prediction. 6. Solve the pressure equation and execute the PISO loop to obtain more accurate pressure and velocity on a non-orthogonal mesh. 7. After pressure correction, solve the turbulent transport equations (such as the turbulent kinetic energy k and dissipation rate ε equations) and update the turbulent viscosity. Determine whether to start the PIMPLE external loop based on the solution of the turbulent transport equations. 8. After completing all calculations in the current time step, advance the time step to the next time point and repeat the above process until the calculation is complete.

[0075] like Figure 3The UML framework diagram of the cavitation model shown includes the basic mixture model (IncompressibleThreePhaseMixture), which provides the foundation for calculating the physical properties of a three-phase coexisting system; the interface properties module (threePhaseInterfaceProperties), the cavitation model base class (threePhaseCavitationMixtures), the Schnerr-Sauer concrete implementation (Schnerr-Sauer Model), and other cavitation models (Other Cavitation Model). This framework first achieves strict mass conservation of the three-phase system at the basic physical property layer, and ensures the physical continuity of physical parameters at the phase interface through intelligent coupling calculation of multiphase viscosity and density; at the same time, it introduces an adaptive compression mechanism and contact angle correction algorithm at the interface dynamics layer, which effectively suppresses the interface ambiguity phenomenon caused by numerical diffusion.

[0076] The following section provides a detailed description of the solver components based on OpenFOAM. In step S1, a three-phase cavitation model is constructed using the VoF multiphase flow interface capture method and the two-phase Schnerr-Sauer method. VoF is a numerical technique in computational fluid dynamics (CFD) used to track free surfaces or multiphase flow interfaces. In OpenFOAM, the VOF method is typically used to simulate the flow of two or more immiscible fluids (such as water and air), where the phase fraction (usually alpha) represents the volume fraction of one fluid in each cell.

[0077] The specific steps include: Step S11, constructing the control equations for the three-phase cavitation of water, steam and air;

[0078] In the governing equations for three-phase cavitation of water, steam, and air, the continuity equation is:

[0079]

[0080] According to the assumptions of the VoF model, ρ m The average density is obtained by weighting the volume fractions of each phase within a control volume, where U is the velocity vector and t is time.

[0081] The momentum equation is:

[0082]

[0083] Where p is pressure, τ is viscous force, g is gravitational acceleration, and f σ Let be the surface tension between the three phases, and its expression is as follows:

[0084]

[0085] Where σm is the surface tension coefficient and κ is the curvature. The unit normal vector of the interface.

[0086] (x-xs) represents the distance to the interface.

[0087] Step S12: Construct a three-phase cavitation model based on the two-phase Schnerr-Sauer model to provide a method for calculating the mass transfer term.

[0088] Since the surface tension between air and water vapor is too small, considering only the surface tension between water and the other two phases, the surface tension coefficient is weighted and averaged based on the volume fraction α, and the expression for the surface tension coefficient is as follows:

[0089]

[0090] The volume fraction transport equations for each phase under the VoF method are shown below:

[0091]

[0092] Where α l α v α nc These are the volume fractions of water, water vapor, and air, respectively; U clv U clnc The artificial compression rates are for water-water vapor and water-air, respectively. This is to address the interface ambiguity caused by numerical diffusion inherent in the VoF method by performing artificial compression at the interface. At the water vapor-air interface, since the two phases have a certain degree of physical miscibility, no artificial compression is required.

[0093] The mass transfer term used to describe cavitation is derived from the governing equations and phase volume transport equations mentioned above. It represents the mass transfer rate between the liquid and vapor phases per unit time.

[0094]

[0095] Where ρl, ρv, and ρnc are the densities of the liquid phase, vapor phase, and non-condensable phase, respectively.

[0096] S13. Select the discretization method and construct the pressure-velocity decoupling algorithm according to the problem to be solved.

[0097] Based on the fundamental assumptions of the Schnerr-Sauer three-phase cavitation model, namely that a unit volume of liquid phase contains a fixed number of spherical gas nuclei, and the vapor phase is composed of gas nuclei growing to different radii, the equations are closed as follows:

[0098]

[0099] Where R is the bubble radius, n0 is the number of bubbles per unit volume of liquid phase, and DR / Dt is the rate of change of bubble radius, obtained from the Rayleigh-Plesset equation, neglecting the interaction between bubbles:

[0100]

[0101] Where pthreshold is the phase transition pressure threshold. It should be noted that the mass derivative of αv mentioned above... The reason for writing it in two forms is to improve computational stability during later calculations.

[0102] To avoid computational bias caused by the selection of the phase transition coefficient, the phase transition coefficient is adjusted based on the turbulence timescale. According to the linear turbulence assumption, the turbulence timescale is inversely proportional to the strain rate; therefore, the strain rate is used for direct correction, and normalization is performed using the mainstream timescale.

[0103]

[0104] Where cv is the original evaporation phase transition coefficient, and t∞ is the mainstream time scale.

[0105] After the model is built, the finite volume method (FVM) is selected to discretize the equations, and a co-located grid is selected to store the physical quantities.

[0106] To address the transient nature of the three-phase cavitation problem, the PIMPLE algorithm is used to decouple pressure and velocity. The calculation process is as follows: Figure 2 .

[0107] In step S2, a base class for the Schnerr-Sauer three-phase cavitation mixture model is constructed based on the library files and three-phase flow and interface processing classes provided by OpenFOAM. The three-phase cavitation mixture is a core class in OpenFOAM designed specifically for the physical property management of three-phase cavitation flow, which is responsible for processing the physical properties, phase change processes and mixture properties of the three-phase system of water, steam and air.

[0108] This includes step S21, which extends the base class of the two-phase cavitation mixture to the base class of the three-phase cavitation mixture. Specifically, the base class of the two-phase cavitation mixture (phaseChangeTwoPhaseMixtures) is extended to the three-phase to construct the threePhaseCavitationMixtures class of the three-phase cavitation model.

[0109] in, phaseChangeTwoPhaseMixtures is a base class of OpenFOAM used to define the basic framework for two-phase flow cavitation models. It supports a general interface for mass transfer (evaporation / condensation) and is designed for incompressible two-phase flows, but can be extended to support compressible flows.

[0110] Step S22: Make the base class of the three-phase cavitation mixture (threePhaseCavitationMixtures) inherit the interface properties module (threePhaseInterfaceProperties) and the basic mixture model (incompressibleThreePhaseMixture) of OpenFOAM.

[0111] After inheriting the threePhaseInterfaceProperties class for three-phase calculation and three-phase interface processing and the incompressibleThreePhaseMixture class, the constructor is modified synchronously and the original correct() function is divided into viscosity coefficient correction and interface tension correction.

[0112] Since this model requires the calculation of the three-phase volume fraction transport equation, an additional source term function vDotAlphav for the vapor phase volume fraction transport equation needs to be added. Its calculation method is as follows:

[0113]

[0114] It has no physical meaning in itself; it is merely an intermediate quantity in the calculation process.

[0115] To implement the runtime selection functionality, declare and implement the runtime selection table in the class, and register the base class on the runtime selection table.

[0116] It is also important to note that the original private member in the incompressibleThreePhaseMixture class should be changed to protected to ensure that it can be used by inherited classes.

[0117] The specific method in step S3 is as follows:

[0118] In the createFields.H header file used by the solver to initialize the fields, add the autoPtr... <phasechangethreephasemixture>The `mixture` object is used, and a dynamic selection mechanism is employed to point to a subclass object using the `New` function. Later at runtime, the `SchnerSauer` class can be selected using `TransportProperties`. In OpenFOAM, `mixture` is the core framework for handling the physical properties of multi-component mixtures. It calculates thermodynamic and transport properties through mass fraction weighted averaging or implicit polynomial coefficient superposition.

[0119] This application utilizes a mixture to construct the volume fraction fields and density constants of the three phases. Furthermore, it replaces the interface pointers in the original two-phase cavitation model with a mixture.

[0120] In the alphaEqn.H header file used by the solver to calculate the phase volume fraction transport equations, the parameters for obtaining the mixture, including the compressibility factor, vDotAlphal, and vDotAlphav, are rewritten. Next, the fluxes for each phase volume fraction transport equation are constructed and uniformly stored in the PtrList class (a list of pointers to dynamically allocated objects) object alphaPhiCorrs (used for volume fraction or phase fraction flux correction in multiphase flow simulations). This facilitates subsequent MULES correction. MULES is an algorithm specifically designed for numerical calculations of volume fractions (or scalar fields) in multiphase flows. Its core function is to ensure the boundedness of physical quantities during calculation (e.g., the volume fraction is always between 0 and 1), while suppressing numerical oscillations and improving the accuracy and stability of interface capture.

[0121] Because three-phase calculations are involved, the method for limiting fluxes using the MULES algorithm differs significantly from that used in two-phase cavitation calculations. First, fluxes for all three phases need to be constructed, including both the initial flux and the artificial compression flux, with the artificial compression flux added only between the water-vapor and water-air phases.

[0122] For the aqueous phase, the flux construction code is as follows:

[0123]

[0124] The other two phases have similar structures to the aqueous phase, so we will not impose too many restrictions on them here.

[0125] After constructing the three-phase fluxes, flux limits need to be imposed on each phase based on global and local extrema: First, calculate the maximum absolute value of the allowable inflow and outflow anti-diffusion flux for each unit. For the aqueous phase in any unit limiter, the absolute value of its maximum inflow flux is Q. - and the absolute value of the maximum outflow flux Q + As shown below:

[0126]

[0127] For the vapor phase in any unit limiter, the absolute value of its maximum inflow flux Q - and the absolute value of the maximum outflow flux Q + As shown below:

[0128]

[0129] For the air phase in any unit limiter, the absolute value of its maximum inflow flux Q - and the absolute value of the maximum outflow flux Q + As shown below:

[0130]

[0131] Where α l,max α l,min α v,max α v,min α nc,max α nc,min Let α represent the maximum and minimum allowable values ​​for the liquid phase, vapor phase, and non-condensable phase, respectively. l,o α v,o α nc,o Vc represents the volume fractions of the liquid, vapor, and non-condensable phases in the previous time step, where Vc is the control volume and F is the volume fraction of the control volume. BD,o This refers to the low-order format flux of the previous time step. These are the data parameters defined in the base class of OpenFOAM for three-phase cavitation mixtures, specifically as defined in the OpenFOAM base class threePhaseCavitationMixtures.

[0132] Next, based on the calculated maximum anti-diffusion flux, the limiter reduction factor λ for each cell is calculated. Within the OpenFOAM framework, an iterative method is used to calculate the outflow and inflow flux limiters for each cell, resulting in a larger λ compared to the traditional FCT method, thus achieving higher accuracy. - and λ + These are limiters for the inflow and outflow fluxes on each face, used to control the fluxes generated by the higher-order scheme. Their iterative formulas are as follows, where λ0 = 1:

[0133]

[0134] Where F AD,o For high-order flux formats.

[0135] Finally, the reduction factor λ of each unit is discretized onto each interface according to the following formula, thereby limiting the flux at each interface:

[0136]

[0137] Preferably, after all three phases have been solved, further adjustments to the correction flux are needed to ensure that each phase does not exceed its limits while maintaining the sum of the volume fractions of the three phases within each cell at 1. For each surface, the sum of the correction fluxes of each phase should be 0 to satisfy the physical definition that the sum of the volume fraction fluxes of each phase on this surface is 1, and this also ensures that the sum of the volume fractions within the cell is 1.

[0138] Therefore, it is necessary to calculate the sum of the corrected flux on each surface:

[0139] sum=ΣΦ i

[0140] Calculate its residual R:

[0141] R = sum - 0

[0142] Based on the sign and magnitude of the residuals, the inflow correction flux (negative flux) or outflow correction flux (positive flux) on this surface is reduced by a certain scaling factor λ.

[0143] If R is positive, the outflow correction flux is reduced; if R is negative, the inflow correction flux is reduced.

[0144] The reduction factor λ is calculated using the following formula:

[0145]

[0146] During this process, since only the modified flux is reduced, the original limiting effect is not destroyed. Instead, the modification becomes more conservative, resulting in the sum of the modified fluxes of each phase at each interface being 0 and the sum of the volume fractions inside the unit limiter being 1.

[0147] Finally, the base flux, i.e., the first-order upwind flux, is added to the corrected flux to form the complete flux. The following three equations are then explicitly solved using the MULES::explicitSolve function:

[0148]

[0149] These three equations are the equations that include the basic flux (low-order flux) and the anti-diffusion flux. Their calculation results meet certain accuracy requirements and ensure that they do not exceed the limits.

[0150] A new solver main program, 3PhaseChangeFoam.C, was then constructed. Through the construction of three-phase MULES, the solution of the phase volume fraction transport equation can avoid obtaining non-physical results, i.e., results greater than 0 or less than 1, while selecting a high-precision discretization scheme.

[0151] In another alternative embodiment, step S4 reconstructs the coefficients mDotAlpha and mDotP for the source terms of the phase equation and pressure equation based on the original SchnerrSauer class.

[0152] First, the rRb function in the Schnerr-Sauer source code is rewritten according to the following equation, which is derived based on the assumptions of the Schnerr-Sauer model:

[0153]

[0154] The pCoeff function is then rewritten according to the following equation:

[0155]

[0156] Here, psat refers to the phase change pressure threshold pthreshold mentioned earlier. This equation is from the aforementioned... Part of it is written as a function for the convenience of constructing equations later.

[0157] Furthermore, the mDotAlphal function used for the source terms of the phase volume transport equation and the mDotP function used for the source terms of the pressure equation are rewritten respectively.

[0158] The mDotAlphal function is rewritten according to the following equation, which is based on... The derivation yields that represent the condensation mass transfer rate and the evaporation mass transfer rate per unit phase volume fraction, respectively:

[0159]

[0160] The mDotP function is rewritten based on the following equation, which is based on The derivation yields that represent the condensation mass transfer rate and the evaporation mass transfer rate under a unit pressure difference, respectively:

[0161]

[0162]

[0163] The above as well as None of them have any actual physical meaning; they are constructed for the convenience of building equations and matrices and to improve computational stability.

[0164] This completes the construction of the three-phase cavitation model. The UML framework diagram of the specific cavitation model can be found here. Figure 3 .

[0165] In step S5, CMake is used for compilation. First, CMakeLists.txt is written, and the OpenFOAM version environment is set to v2312.

[0166] The threePhaseInterfaceProperties.C, incompressibleThreePhaseMixture.C, Schnerr-Sauer.C (the three-phase cavitation model built in steps S2 and S4), and threePhaseCavitationMixtures.C (the base class of the three-phase cavitation model) are compiled into threePhaseCavitationMixtures.so using CMake for subsequent linking. Similarly, threePhaseChangeFoam.C (built in step S3) is compiled into a new three-phase cavitation solver, threePhaseChangeFoam, by using CMake and linking it to the dynamic link library threePhaseCavitationMixtures.so.

[0167] To test the solver performance, we chose to use... Figure 4 The OpenFOAM native two-dimensional example shown uses a throttle mesh, and the mesh is refined in three layers on both sides of the throttle channel to ensure the accuracy of physical quantity calculation in the high gradient region. The calculation results are compared with those of interPhaseChangeFoam and interFoam under two-phase cavitation and no-cavitation conditions, respectively. Furthermore, the effectiveness of the phase transition coefficient correction is tested under three-phase cavitation conditions.

[0168] Under two-phase airflow conditions, αl = 1 in the initial computational domain, which is the water phase; the left inlet is a pressure inlet with a pressure of 100 bar, and the right outlet is a pressure outlet with a pressure of 1 bar, the same as the initial field; the k-ωSST model is selected as the turbulence model to ensure that the turbulence characteristics can be well described in both the boundary layer and free flow.

[0169] Figure 5 To stabilize the pressure distribution along the throttle valve axis within the range of 1mm before and 2mm after the valve orifice in the throttle channel, it can be seen that the calculation results of the two solvers are basically consistent. The pressure is maintained at a relatively high level of 100 bar, close to the inlet pressure, about 1mm before the valve orifice, with almost no pressure drop.

[0170] Under the two-phase cavitation-free condition, the phase transition coefficient of the three-phase cavitation solver is set to 0, the initial field αl = 1 in the region before the valve orifice, and αnc = 1 in other regions, with other boundary conditions remaining essentially unchanged. Since this example does not directly calculate cavitation, cavitation numbers are used for comparison and verification. The results are as follows: Figure 6 As shown: the sampling location is at the upper edge of the throttling channel, because this is a high-incidence area of ​​cavitation. The calculation results generally show good consistency, with the cavitation number distribution along the channel being basically the same, with only slight differences in the inlet area.

[0171] Under three-phase cavitation conditions, the boundary conditions are basically the same as those under two-phase non-cavitation conditions. Figure 7 The vapor phase volume fraction distribution after the three-phase cavitation flow reaches a steady state is calculated for both the original and modified models. In the figure, black represents the contour lines for αv = 0.9, and green represents the contour lines for αv = 0.5. It can be seen that the cavitation region of the original model is much longer longitudinally than that of the modified model. The latter captures the air backflow phenomenon at the outlet of the throttle valve pipe. This is because mDotP acts as a source term in the pressure correction equation, and the evaporation phase change coefficient increases under the influence of the large strain rate at the outlet, allowing the modified model to capture this feature. Simultaneously, radially, the cavitation cavity captured by the modified model has a greater thickness than that of the original model. This is because model A ignores the influence of the high strain rate near the channel wall on cavitation, thus underestimating the growth of the cavitation cavity. The pressure distribution curves at the midline of the first half of the channel in both models are shown below. Figure 8 As shown, it can be seen that all three models exhibited significant pressure fluctuations at 0.3 mm within the channel. However, the fluctuation amplitude of the original model was much larger than that of the corrected model. This overestimation of the pressure fluctuations further led to the compression of the cavitation cavity thickness at the front of the channel.

[0172] In summary, this patent proposes a method for constructing a three-phase cavitation solver based on the OpenFOAM framework, while also providing a convenient interface for subsequent user customization.

[0173] The above description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention shall be within the scope of protection of the pending claims of the present invention.< / phasechangethreephasemixture>

Claims

1. A method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework, characterized in that, Includes the following steps: Step S1: Construct a Schnerr-Sauer three-phase cavitation model based on the multiphase flow interface capture method and the two-phase Schnerr-Sauer cavitation model; Step S2: Construct a base class for the Schnerr-Sauer three-phase cavitation model; Step S3: Construct the three-phase phase volume fraction transport equation as the main solution program based on the two-phase Schnerr-Sauer cavitation solver. The process of constructing the volume fraction transport equations for the three phases includes: Add a mixture model to the header file used by the two-phase cavitation solver to initialize the field, and use the mixture model to construct the three-phase volume fraction field and density constant; Obtain the parameters of the mixture model from the header file of the two-phase cavitation solver used to calculate the phase volume fraction transport equation; The three-phase fluxes of the phase volume fraction transport equation are constructed. After the three-phase fluxes are constructed, flux limits are imposed on each phase flux based on global and local extrema. The maximum absolute value of the allowable inflow and outflow anti-diffusion fluxes of each unit is calculated. The reduction factor λ of each unit limiter is calculated based on the calculated maximum anti-diffusion flux. The outflow and inflow flux limiters of each unit limiter are calculated by iterative method to obtain the iterative formula. The reduction factor λ of each unit limiter is discretized to each interface, and then the flux of each interface is limited so that the sum of the corrected fluxes of each phase at each interface is 0 and the sum of the volume fractions inside the unit limiter is 1. The base flux is added to the modified flux to form the complete flux, which is then solved explicitly using the explicitSolve function of the MULES algorithm, thereby constructing the three-phase volume fraction transport equation. Step S4: Construct a Schnerr-Sauer three-phase cavitation derived class model based on the two-phase Schnerr-Sauer cavitation model. Step S5: Compile the source files of the Schnerr-Sauer three-phase cavitation derived class model and the corresponding base class into a dynamic library, compile the volume fraction transport equations of the three phases and link them with the dynamic library to build a new Schnerr-Sauer three-phase cavitation solver. The usability of the new Schnerr-Sauer three-phase cavitation solver and the effectiveness of model corrections were tested based on native two-dimensional example meshes.

2. The method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework according to claim 1, characterized in that, Step S1 includes: Step S11: Construct the governing equations for the three-phase cavitation of water, steam, and air; Step S12: Construct a three-phase cavitation model based on the two-phase Schnerr-Sauer model to provide a method for calculating the mass transfer term; Step S13: Select the discretization method and construct the pressure-velocity decoupling algorithm according to the problem to be solved.

3. The method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework according to claim 1, characterized in that, Step S2 includes: Step S21: Extend the base class of the two-phase cavitation mixture to the base class of the three-phase cavitation mixture; Step S22: Make the base class of the three-phase cavitation mixture inherit the interface feature module and the basic mixture model of OpenFOAM; The constructor is modified synchronously, and the original correct() function is divided into viscosity coefficient correction and interfacial tension correction; Add the source term function vDotAlphav for the vapor phase volume fraction transport equation and perform the calculation.

4. The method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework according to claim 1, characterized in that, In step S3, When limiting flux, first calculate the maximum absolute value of the allowed inflow and outflow of the back diffusion flux for each unit limiter.

5. The method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework according to claim 4, characterized in that, The iterative formula for the reduction factor λ of the unit limiter is as follows, where λ 0 =1: ； ； By discretizing the reduction coefficient λ of the unit limiter onto each interface, the flux at each interface is limited.

6. The method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework according to claim 5, characterized in that, Calculate the sum of the corrected fluxes on each facet: ； Calculate its residual R: ； Based on the sign and magnitude of the residuals, the inflow correction flux or outflow correction flux at the interface is reduced. If R is positive, the outflow correction flux is reduced; if R is negative, the inflow correction flux is reduced. This results in the sum of the modified fluxes of each phase at each interface being 0, and the sum of the volume fractions inside the cell confinement being 1.

7. The method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework according to claim 1, characterized in that, In step S4, the new Schnerr-Sauer three-phase cavitation derived class model is constructed using the following method: First, the rRb function in the Schnerr-Sauer source code was rewritten accordingly; Then the mDotAlphal function used for the source terms of the phase volume transport equation and the mDotP function used for the source terms of the pressure equation are rewritten respectively. This completes the construction of the cavitation model.

8. The method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework according to claim 2, characterized in that, In step S11, In the governing equations for three-phase cavitation of water, steam, and air, the continuity equation is: ； According to the assumptions of the VoF model, The average density is obtained by weighting the volume fractions of each phase within a control body. In step S12, since the surface tension between air and water vapor is too small, only the surface tension between water and the other two phases is considered.

9. The method for constructing a Schnerr-Sauer three-phase cavitation solver based on the OpenFOAM framework according to claim 1, characterized in that, In step S13, the equations are closed using the Schnerr-Sauer three-phase cavitation model; The phase transition coefficient is adjusted based on the turbulence timescale and corrected using strain rate according to the turbulence model assumptions. After the model is built, the finite volume method (FVM) is selected to discretize the equations, and a co-located grid is selected to store the physical quantities.

Images