A multiphase flow interaction simulation method based on MPM and phase field model

Abstract

The application discloses a multiphase flow interaction simulation method based on MPM and a phase field model, and comprises the following steps: step S1, storing mass parameters, phase field parameters and momentum parameters of Lagrange particles in an Euler grid; step S2, updating the phase field of the Euler grid according to the chemical potential to control the energy of the whole system; step S3, updating the momentum of the Euler grid according to gravity, elastic force and the surface tension of the updated equivalent phase field form; step S4, after the grid momentum is updated, processing the boundary on the Euler grid by using a slip boundary condition; step S5, after the phase field and the momentum are updated, converting the momentum parameters and the phase field parameters from the Euler grid variables to Lagrange variables; and step S6, updating the state of the Lagrange particles: updating the position and mass state of the Lagrange particles according to the information returned to the Lagrange particles.

Description

Technical Field

[0001] This invention relates to the field of computer fluid simulation, and in particular to a multiphase flow interactive simulation method based on MPM and phase field model. Background Technology

[0002] Physics-based simulation is a hot research area in computer graphics, widely used in film and animation, virtual reality, and other fields. Multiphase flow interaction is a very common phenomenon. Different fluids produce different effects during interaction. For example, oil droplets falling into water create a clear stratification effect, while ink droplets gradually dissolve and diffuse. More complex interactions exist, such as extraction and partial dissolution. Currently, multiphase flow simulation methods can be broadly divided into two categories: the Lagrange method and the Eulerian grid method.

[0003] The Lagrangian method focuses on capturing particle details and is mostly extended to particle solvers (such as SPH and PBF). The literature [Ren B, Li C, Yan X, et al. Multiple-fluid SPH simulation using a mixture model[J]. ACM Transactions on Graphics (TOG), 2014, 33(5): 1-11.] uses volume fractions to represent the distribution of different phases or components and introduces the concept of drift velocity to handle diffusion between components. However, this method does not consider energy and cannot simulate complex multiphase flow interactions such as extraction and partial dissolution. The literature [Yang T, Chang J, Ren B, et al. Fast multiple-fluid simulation using Helmholtz free energy[J]. ACM Transactions on Graphics (TOG), 2015, 34(6): 1-11.] first introduced Helmholtz free energy into multiphase flow interaction systems, extending the simulation of multifluid interactions. However, Lagrangian-based solvers introduce time step constraints and neglect the simulation of the gas-liquid interface. The literature [Ren B, Xu B, Li C. Unified particle system for multiple-fluid flow and porous material[J]. ACM Transactions on Graphics (TOG), 2021, 40(4): 1-14.] introduces a new concept of virtual phase and proposes a unified particle framework for simulating multiple fluid flows and porous materials.

[0004] Eulerian meshes focus on the processing of interfaces at different stages, achieving reasonable results by tracking the interfaces. The literature [Kang N, Park J, Noh J, et al. A hybrid approach to multiple fluid simulation using volume fractions[C] / / Computer Graphics Forum. Oxford, UK:Blackwell Publishing Ltd, 2010, 29(2): 685-694.] combines level set and volume fractions, capturing not only the discontinuous interfaces between immiscible fluids but also the smooth transitions between miscible fluids. The literature [Nielsen MB, Østerby O. A two-continua approach to Eulerian simulation of water spray[J]. ACM Transactions on Graphics (TOG), 2013, 32(4): 1-10.] simulates realistic water mist phenomena by equating each node in the Eulerian space to a bidirectionally coupled continuum of gas and liquid phases.

[0005] The Material Point Method (MPM), as a hybrid Euler-Lagrange solver, has demonstrated good scalability in the simulation of multiphase flow interactions. The paper [Yan X, Li CF, Chen XS, et al. MPM simulation of interacting fluids and solids[C] / / Computer Graphics Forum. 2018, 37(8): 183-193.] proposes a unified MPM framework for simulating gas and solid phases. However, this paper does not address complex multiphase flow interactions, such as extraction.

[0006] In summary, the simulation of multiphase flow interaction still has the following problems at this stage:

[0007] (1) The introduction of energy makes the unified solver based on Lagrange limited by time steps, resulting in very slow phase diffusion, and even requiring the addition of external forces.

[0008] (2) In the current unified solution of miscible and immiscible phases, the gas-liquid interface is mostly ignored, which reduces the realism of the simulation. Summary of the Invention

[0009] The technical problem to be solved by the present invention is to provide a multiphase flow interaction simulation method based on MPM and phase field model, which can solve the technical problems mentioned in the background art above.

[0010] To solve the above-mentioned technical problems, the technical solution adopted by the present invention is: a multiphase flow interaction simulation method based on MPM and phase field model, comprising the following steps:

[0011] Step S1: Transfer the mass parameters, phase field parameters, and momentum parameters stored in the Lagrange particles to the Eulerian grid, respectively.

[0012] Step S2: Update the Eulerian grid phase field by controlling the energy of the entire system based on the chemical potential;

[0013] Step S3: Update the momentum of the Euler mesh based on gravity, elastic force, and surface tension in the updated equivalent phase field form;

[0014] Step S4: After updating the mesh momentum, apply slip boundary conditions to the boundary on the Eulerian mesh;

[0015] Step S5: After updating the phase field and momentum, convert the momentum parameters and phase field parameters from Eulerian grid variables to Lagrangian variables;

[0016] Step S6: Update the state of the Lagrange particle: Update the position and mass state of the Lagrange particle based on the information transmitted back.

[0017] A further improvement to the technical solution of the present invention is that the specific process of step S1 is as follows:

[0018] Step S11: Transfer the mass parameters stored in the Lagrangian particles to the Eulerian grid: Convert the Lagrangian mass to the Eulerian grid mass using APIC as the transfer method. This process satisfies the following formula:

[0019] ,

[0020] in, Represents Euler nodes quality Represents the Lagrange particle quality These are the weights in the continuous quadratic B-spline interpolation function of the Lagrange particle and Eulerian grid mapping;

[0021] Step S12: Transfer the momentum parameters stored in the Lagrange particles to the Eulerian grid: The transfer process satisfies momentum conservation, using APIC as the transfer method, and the process satisfies the following formula:

[0022] ,

[0023] in, For Lagrange particles affine velocity, Represents the Lagrange particle volume, Euler nodes Location, For Lagrange particles Location;

[0024] Step S13: Transfer the phase field parameters stored in the Lagrange particles to the Eulerian grid: Use mass fraction as the phase field parameter to reflect the compatibility of the multiphase fluid. Interactive systems only need to solve The phase field equations, and the transfer process only requires the transfer of... For each phase field parameter, selecting the MPM transfer method, the phase field parameter transfer satisfies the following formula:

[0025] ,

[0026] in, This represents the phase field parameters of the Eulerian grid. This represents the phase field parameters of the Lagrange particle.

[0027] A further improvement to the technical solution of the present invention is that the process of updating the Eulerian grid phase field according to the chemical potential in step S2 is as follows:

[0028] Step S21: Solve for chemical potential: Chemical potential characterizes the transfer trend between phase components and controls the development of the entire system from high chemical potential to low chemical potential. After interpolating the phase field parameters to the Eulerian grid, the central grid and central difference are used to discretize them.

[0029] The governing equations for the phase field are:

[0030] ,

[0031] in, Represents the spatial region of the system. It is a constant used to control the diffusion rate. This represents the interfacial energy, used to reduce the gradient of the mass fraction field, making surface changes smoother. It represents physical capabilities and controls how fluids interact. This represents the phase field value at position x at time t;

[0032] The governing equation for the chemical potential is:

[0033] ,

[0034] Represents chemical potential, Indicates physical fitness. It is a constant. Represents phase field parameters, Represents the Laplace operator;

[0035] Step S22: Selecting the volumetric energy function: The local minimum value of the volumetric energy function is used to control the low-energy state of the fluid. The process of the fluid evolving from a high-energy state to a low-energy state is simulated. Different volumetric energy functions are designed to achieve different interactive effects. Two-phase interaction is displayed using a volumetric energy function graph, and the minimum value in the effective region is marked. Three-phase interaction uses a Gibbs triangle as the energy display. The design of the volumetric energy function follows the following principles: the minimum energy value is located within the effective region to ensure that the energy minimization state can be achieved in the end; the partial derivative of the chemical potential with respect to a certain phase should only be related to the mass ratio of that phase.

[0036] Step S23: Update the phase field: After designing the corresponding bulk energy function for different multiphase flow interactions, the phase field is updated using the chemical potential. The degenerate mobility is set as a constant, and the phase field is discretized using a central grid and central difference. The phase field update satisfies the following formula:

[0037] ,

[0038]

[0039] in, Represents the mass derivative. Indicates degenerate mobility. It represents chemical potential.

[0040] A further improvement of the technical solution of the present invention lies in: in step S3, an updated phase field is used to capture the interface between different components, surface tension is applied, and momentum is further updated. The surface tension control equation satisfies:

[0041] ,

[0042] in, It uses the equivalent surface tension in the form of phase field. express , The surface tension coefficient between phases express Phase field parameters of the phase, express Phase field parameters of the phase, Coefficients related to the diffusion interface.

[0043] A further improvement to the technical solution of the present invention is that: in step S3, the momentum of the Euler grid is updated according to gravity and elastic force according to the following formula:

[0044] ,

[0045] in, Indicates elastic force. This represents the grid spacing of the Euler mesh. Let be the initial volume of the Lagrange particle. Indicates Piolakirchhoff stress, express Lagrange particles at time Deformation gradient, Represents the Lagrange particle The transpose of the deformation gradient matrix, where T denotes the transpose. Represents gravity. Represents gravitational acceleration. Represents the Euler mesh The sum of forces acting on the body Indicates the current moment. Indicates the next moment, Indicates a time step, Represents the Euler mesh The speed in the next moment, Represents the Euler mesh The speed at the current moment.

[0046] A further improvement to the technical solution of this invention lies in the following: after updating the mesh momentum in step S4, a slip boundary condition is applied to the boundary on the Eulerian mesh. The slip boundary condition is as follows:

[0047] ,

[0048] in, This indicates the mesh velocity where the collision has not yet ended. Collision processing is applied to Lagrange particle pairs multiple times before updating their positions. Represents the surface normal.

[0049] A further improvement to the technical solution of the present invention is that the specific process of step S6 is as follows:

[0050] Step S51: After updating the momentum parameters on the Eulerian grid, transfer the momentum parameters back from the Eulerian grid to the Lagrange particle: Using APIC as the transfer method, the momentum parameter transfer satisfies:

[0051]

[0052] in, Represents the Lagrange particle speed, Represents the Euler mesh speed;

[0053] Step S52: After updating the phase field parameters on the Eulerian grid, transfer the phase field parameters back from the Eulerian grid to the Lagrangian particles. The phase field information transferred back to the Lagrangian particles is the mass fraction of each phase component. This process follows the formula below:

[0054] .

[0055] A further improvement to the technical solution of the present invention is that the specific process of step S7 is as follows:

[0056] Step S61: After transmitting the momentum information back to the Lagrange particle, update the position of the Lagrange particle based on the momentum information. The updated position satisfies:

[0057] ,

[0058] in, Represents the Lagrange particle at the current moment. Location, Indicates the Lagrange particle at the next moment Location, Indicates a time step, This indicates the velocity of the Lagrange particle at the next instant;

[0059] Step S62: Update the Lagrange particle mass: Update the mass of the Lagrange variables according to the phase field parameters. The mass is not updated during the transfer process. After one time step, the mass is updated, and the updated mass satisfies the following formula:

[0060] ,

[0061] in, Indicates the Lagrange particle at the next moment quality express Phase mass fraction, express Phase quality.

[0062] The technological advancements achieved by this invention due to the adoption of the above technical solutions are as follows:

[0063] 1. This invention proposes a solver that uniformly handles diffusion effects and chemical reactions between immiscible and immiscible phases. Multiphase flow interactions produce a wealth of visual effects; processing each effect separately would be a huge workload and would severely limit the system's scalability and compatibility. The system proposed in this invention highly unifies the characteristics of multiphase flow interactions. This system is applicable to the simulation of almost all multiphase flow interactions and can be easily extended without much modification.

[0064] 2. This invention provides an energy-based solver that allows users to define their own volumetric energy functions to control the evolution of the entire system. In nature, systems evolve towards energy minimization. Based on this principle, this invention proposes an energy-based solver to simulate phenomena such as partial dissolution and extraction, overcoming the limitations of traditional multiphase flow solvers. Different volumetric energy functions produce completely different visual effects, allowing users to define their own volumetric energy functions and expand the simulation of multiphase flow interactions.

[0065] 3. This invention combines MPM and PFM to model the entire MPM mesh, enabling rapid capture of interfaces between different phases, including the gas-liquid interface. The composite model of MPM and Phase Field PFM effectively combines the advantages of the Lagrangian and Eulerian methods. MPM uses Eulerian meshes and Lagrangian (material) points to represent material properties. During motion, the Eulerian mesh remains unchanged while the Lagrangian points move, preserving motion details. PFM is introduced to model the background mesh of the MPM, allowing for rapid capture of interfaces between different phases, including fluid-gas interfaces, without the need for interface tracking. This method has physical significance not only at the interface but also in the bulk phase. Attached Figure Description

[0066] Figure 1 This is a flowchart of the present invention;

[0067] Figure 2 This is a schematic diagram of the discretization of the Hamiltonian operator of this invention;

[0068] Figure 3(a) is a volume energy function image of the two-phase immiscible system of the present invention;

[0069] Figure 3(b) is a volume energy function image of partial dissolution of two phases in the two-phase flow system of the present invention;

[0070] Figure 4(a) is a volume energy function graph of the three-phase flow system of the present invention, illustrated by the Gibbs triangle;

[0071] Figure 4(b) is a volume energy function image of the three-phase flow system of the present invention, using the Gibbs triangle energy contour lines as an example;

[0072] Figure 4(c) is a comparison diagram of the energy magnitude and color of the Gibbs triangle of the present invention;

[0073] Figure 5 This is a diagram illustrating the effect of gas-liquid interface capture in this invention;

[0074] Figure 6 This is a diagram illustrating the effect of two-phase immiscibility in the two-phase flow system of this invention;

[0075] Figure 7 This is a diagram illustrating the effect of partial dissolution of two phases in the two-phase flow system of this invention;

[0076] Figure 8 This is a diagram showing the extraction effect in the three-phase flow system of this invention. Detailed Implementation

[0077] The present invention will be further described in detail below with reference to embodiments:

[0078] like Figure 1 As shown, a multiphase flow interaction simulation method based on MPM and phase field model includes the following steps:

[0079] Step S1: Transfer the mass parameters, phase field parameters, and momentum parameters stored in the Lagrange particles to the Eulerian grid, respectively.

[0080] The specific process of step S1 is as follows:

[0081] Step S11: Transfer the mass parameters stored in the Lagrangian particles to the Eulerian grid: Convert the Lagrangian mass to the Eulerian grid mass using APIC as the transfer method. This process satisfies the following formula:

[0082] ,

[0083] in, Represents Euler nodes quality Represents the Lagrange particle quality These are the weights in the continuous quadratic B-spline interpolation function of the Lagrange particle and Eulerian grid mapping.

[0084] In the hybrid Eulerian-Lagrange method, Lagrange particles are typically used to store all state information, while Eulerian grids are used to assist in calculations. This process requires multiple information transfers between Lagrange particles and Eulerian grids. Therefore, choosing an information transfer method with good numerical properties is crucial.

[0085] APIC is a method for transferring information between Lagrange particles and Eulerian grids. This method eliminates particle dissipation, prevents information loss, and ensures the precise conservation of transferred angular momentum. In addition to position and velocity, the APIC method adds an Affine matrix to each Lagrange particle to record its divergence and angular velocity.

[0086] Step S12: Transfer the momentum parameters stored in the Lagrange particles to the Eulerian grid: The transfer process satisfies momentum conservation, using APIC as the transfer method, and the process satisfies the following formula:

[0087] ,

[0088] in, For Lagrange particles affine velocity, Represents the Lagrange particle volume, Euler nodes Location, For Lagrange particles The location.

[0089] Step S13: Transfer the phase field parameters stored in the Lagrange particles to the Eulerian grid: Solving the phase field equations requires multiple calculations of the Hamiltonian operator. To more quickly discretize it, the phase field update will be performed on the Eulerian grid. To better illustrate the transfer of phase components, mass fraction is used as the phase field parameter to reflect the compatibility of the multiphase fluid. Interactive systems only need to solve The phase field equations, and the transfer process only requires the transfer of... For each phase field parameter, selecting the MPM transfer method, the phase field parameter transfer satisfies the following formula:

[0090] ,

[0091] in, This represents the phase field parameters of the Eulerian grid. This represents the phase field parameters of the Lagrange particle.

[0092] The Material Point Method (MPM) is a Galerkin method where each MPM particle represents a continuous block of material. It uses a continuous description of the governing equations. The MPM uses Eulerian meshes for collision detection and momentum updates, and Lagrangian particles for state tracking (e.g., advection and deformation). To eliminate Lagrangian particle dissipation, prevent information loss, and ensure accurate conservation of transferred angular momentum, Affine Particle-In-Cell (APIC) is used for information transfer between the Eulerian mesh and the Lagrangian particles. To prevent element cross-instability, a continuous quadratic B-spline interpolation function is used.

[0093] Step S2: Update the Eulerian grid phase field by controlling the energy of the entire system based on the chemical potential;

[0094] The phase-field method (PFM) is a mathematical approach for handling interfaces. Solving the phase-field equations involves a fourth-order spatial derivative, and the sensitivity of higher-order derivatives introduces time-step limitations in particle systems. This negates the computational advantage of the fast solution of the phase-field method. Furthermore, because particle systems use kernel functions to interpolate surrounding particles, it is difficult to model air, requiring further separate processing for operations such as capturing gas-liquid interfaces. Therefore, updating the phase field on an Eulerian grid is chosen.

[0095] The specific process of updating the phase field of the Eulerian mesh in step S2 is as follows:

[0096] Step S21: Solve for the chemical potential: The chemical potential characterizes the transfer trend between phase components and controls the entire system's evolution from high to low chemical potential. After interpolating the phase field parameters to an Eulerian grid, discretize them using a central grid and central difference. A schematic diagram of the discretization is attached. Figure 2 As shown.

[0097] One of the core ideas of the phase-field method is to use the gradient of the order parameter to track the interface. The phase states of various fluids define a potential energy called free energy. The evolution of the phase states minimizes the energy of the entire system.

[0098] The governing equations for the phase field are:

[0099] ,

[0100] in, Represents the spatial region of the system. It is a constant used to control the diffusion rate. This represents the interfacial energy, used to reduce the gradient of the mass fraction field, making surface changes smoother. It represents physical capabilities and controls how fluids interact. This represents the phase field value at position x at time t;

[0101] The governing equation for the chemical potential is:

[0102] ,

[0103] Represents chemical potential, Indicates physical fitness. It is a constant. Represents phase field parameters, Represents the Laplace operator;

[0104] Step S22: Selecting the volumetric energy function: The local minimum value of the volumetric energy function is used to control the low-energy state of the fluid. The process of the fluid evolving from a high-energy state to a low-energy state is simulated. Different volumetric energy functions are designed to achieve different interactive effects. Two-phase interaction is displayed using a volumetric energy function graph, and the minimum value in the effective region is marked. Three-phase interaction uses a Gibbs triangle as the energy display. The design of the volumetric energy function follows the following principles: the minimum energy value is located within the effective region to ensure that the energy minimization state can be achieved in the end; the partial derivative of the chemical potential with respect to a certain phase should only be related to the mass ratio of that phase.

[0105] Different volumetric energy functions can achieve different interactive effects. Figure 3(a) shows the volumetric energy function graphs for two-phase immiscible systems. The minimum values ​​of this volumetric energy function are respectively located in... The system evolves towards energy minimization, eventually reaching a local minimum in the effective region—the pure state where the fluid has not dissolved. This phenomenon is captured using the local minimum of the volume energy function. Figure 6 This is a diagram illustrating the effect of two-phase immiscibility in a two-phase flow system. Blue and green liquids are simultaneously injected into the container, and the entire evolution process is driven by a volume energy function. As shown in Figure 3(a), the minimum value of the effective region is obtained in the pure state where the two phases have not dissolved. Therefore, the effect of two-phase immiscibility is simulated under the drive of this volume energy function. Although the two phases are immiscible, there is still interaction in the dynamics. The surface tension in the phase field equivalent form is used to capture the interface between the two phases. The experimental results show that the blue and green liquids exhibit good stratification throughout the injection process. After stabilization, a clear interface can be observed. Because APIC is used as the information transmission method and a quadratic B-spline function as the interpolation function, it can be observed that the fluid is very "active" throughout the interaction process, without any phenomenon of fluid decreasing, and the energy loss is very small.

[0106] Figure 3(b) shows a two-phase flow system. The volume energy function of partial dissolution. Two local minima, located at... The fluid is obtained from [location]. Assume the fluid is initially a pure phase. Then, for pure [phase]... Phase, reaching the minimum value For pure Phase, reaching the minimum value The free energy function drives the system to minimize energy, and partial dissolution is captured by utilizing local minima. (Appendix) Figure 7 This is a diagram showing the effect of partial dissolution in a two-phase flow system.

[0107] Equal volumes of red and green liquids partially dissolve. Using the red liquid as a pure solvent... Phase, green liquid as pure Phase. The interaction process is driven by the volumetric energy function given in Figure 3(b). For the red liquid, the regional minimum is in The value was obtained from [location], therefore, after reaching a steady state, it shows a red proportion. Green percentage A reddish saturated substance. For green liquids, the regional minimum is at... The value was obtained from [location], therefore, after reaching a steady state, it shows a red proportion. Green percentage The greenish saturated phase eventually produces a mixture of two fully saturated phases. Applying surface tension using phase-field equivalence reveals a clear interface between the two phases. (Appendix) Figure 7 With appendix Figure 6 This invention realizes two different states of biphase interaction, the difference between which is entirely driven by the volume function. Therefore, the solver proposed in this invention allows users to design different volume functions to simulate the rich interaction effects of multiphase flow. For miscible and immiscible phases, separate processing is not required; only appropriate volume functions need to be designed. Therefore, the solver proposed in this invention is a solver that handles miscible and immiscible phases uniformly.

[0108] Figure 4(a) shows the volume energy function represented by the Gibbs triangle in the extraction simulation of the three-phase flow system (A, B, C). The Gibbs triangle was chosen for energy representation. The vertices of the triangle represent the pure fluid state (…). The mass fractions of parallel lines (of fluids A, B, or C) are the same. See Figure 4(a) for " "This indicates" Figure 4(b) shows the volumetric energy function represented by Gibbs triangle energy contour lines in the extraction simulation of the three-phase flow system (A, B, C). The volumetric energy is plotted in color, with white dashed lines representing contour lines. The energy magnitude corresponding to each color is shown in Figure 4(c). From the energy diagram, it can be observed that phase B has a lower minimum energy to dissolve in phase A than it does in phase C. Therefore, when phase A is added to the mixture of phases B and C, phase B, which was originally dissolved in phase C, will redissolve in phase A, achieving a chemical extraction effect. Figure 8This diagram illustrates the extraction effect in a three-phase flow system. Initially, the upper layer is a mixture of red and green liquids, while the upper layer, blue liquid, is the extractant. The separatory funnel is shaken thoroughly. The entire interaction process is driven by the volumetric energy function shown in Figure 4(b). The separatory funnel is shaken thoroughly to ensure complete mixing of the blue liquid with the mixture. Due to the greater solubility of the red liquid in the blue liquid, the blue liquid extracts the red liquid from the mixture. After settling, the upper layer is the remaining liquid from the initial red liquid extraction, and the lower layer is the mixture after the red liquid has redissolved in the blue liquid. This achieves the effect of chemical extraction.

[0109] Step S23: Update the phase field: After designing the corresponding bulk energy function for different multiphase flow interactions, the phase field is updated using the chemical potential. The degenerate mobility is set as a constant, and the phase field is discretized using a central grid and central difference. The phase field update satisfies the following formula:

[0110] ,

[0111]

[0112] in, Represents the mass derivative. Indicates degenerate mobility. It represents chemical potential.

[0113] Step S3: Discretization of MPM momentum conservation is based on the material's potential energy. Surface tension is considered during the evolution process, and the surface tension formula based on the level set is replaced with an equivalent phase field form. It is assumed that the material's internal energy includes potential energy related to surface tension, elastic force, and gravity.

[0114] The momentum of the Euler mesh is updated based on gravity and elastic force according to the following formula:

[0115] ,

[0116] in, Indicates elastic force. This represents the grid spacing of the Euler mesh. Let be the initial volume of the Lagrange particle. Indicates Piolakirchhoff stress, express Lagrange particles at time Deformation gradient, Represents the Lagrange particle The transpose of the deformation gradient matrix, where T denotes the transpose. Represents gravity. Represents gravitational acceleration. Represents the Euler mesh The sum of forces acting on the body Indicates the current moment. Indicates the next moment, Indicates a time step, Represents the Euler mesh The speed in the next moment, Represents the Euler mesh The speed at the current moment.

[0117] Surface tension in the equivalent phase field form refers to the tension acting along any boundary line on the surface of a liquid due to the imbalance of molecular attraction. In the phase field method, the classic sharp interface between two immiscible fluids is replaced by a transition region. In this transition region, the fluid properties continuously change, and the equivalent phase field replaces surface tension. Figure 5 This image shows the effect of using an equivalent phase field to replace surface tension in capturing the gas-liquid interface. Traditional multiphase flow solvers based on Lagrange particles struggle to capture the gas-liquid interface. The solver proposed in this invention models the entire MPM mesh, treating the gas as one phase in a two-phase flow. An equivalent phase field is used to replace surface tension in capturing the gas-liquid interface. The image shows that the movement of surface layer particles gradually reduces the surface area, consistent with the effect of surface tension. The Eulerian mesh momentum is updated based on gravity, elastic force, and the surface tension in the updated equivalent phase field form.

[0118] Using an updated phase field, the interface between different components is captured, surface tension is applied, and the momentum is further updated. The surface tension governing equation satisfies:

[0119] ,

[0120] in, It uses the equivalent surface tension in the form of phase field. express , The surface tension coefficient between phases express Phase field parameters of the phase, express Phase field parameters of the phase, Coefficients related to the diffusion interface.

[0121] Step S4: After updating the mesh momentum, apply slip boundary conditions to the boundaries on the Eulerian mesh;

[0122] Because particle-based collisions can cause inconsistencies in deformation gradients, artifacts can appear in the simulation. To avoid these artifacts, after updating the mesh momentum, slip boundary conditions are applied to the boundaries of the Eulerian mesh.

[0123]

[0124] in, This indicates the mesh velocity where the collision has not yet ended. Collision processing is applied to Lagrange particle pairs multiple times before updating their positions. Represents the surface normal.

[0125] Step S5: After updating the phase field and momentum, the momentum and phase field parameters are converted from Eulerian grid variables to Lagrangian variables. The basic process is consistent with the information transfer from Lagrangian particles to Eulerian grids and does not affect the mass. The transfer is performed because the system assumes that the mass of the Lagrange particle remains constant within a time step.

[0126] The specific process of step S5 is as follows:

[0127] Step S51: After updating the momentum parameters on the Eulerian grid, transfer the momentum parameters back from the Eulerian grid to the Lagrange particle: To maintain high consistency with the information transfer from the Lagrange particle to the Eulerian grid, APIC is used as the transfer method, and the momentum parameter transfer satisfies:

[0128]

[0129] in, Represents the Lagrange particle speed, Represents the Euler mesh speed;

[0130] Step S52: After updating the phase field parameters on the Eulerian grid, transfer the phase field parameters back from the Eulerian grid to the Lagrangian particles. The phase field information transferred back to the Lagrangian particles is the mass fraction of each phase component. This process follows the formula below:

[0131] .

[0132] Step S6: Update the state of the Lagrange particle: Update the position and mass state of the Lagrange particle based on the information transmitted back.

[0133] The specific process of step S6 is as follows:

[0134] Step S61: After transmitting momentum information back to the Lagrange particle, update the Lagrange particle's position based on the momentum information. Since the boundaries have already been constrained on the Eulerian mesh, it is not necessary to constrain the Lagrange particle's position in this step. The updated position satisfies:

[0135] ,

[0136] in, Represents the Lagrange particle at the current moment. Location, Indicates the Lagrange particle at the next moment Location, Indicates a time step, This indicates the velocity of the Lagrange particle at the next instant;

[0137] Step S62: Update Lagrange particle mass: Update the mass of the Lagrange variable according to the phase field parameters. Since the properties (such as density) of each fluid substance are different, different parameters can be set for different phases. In this invention, each phase is allowed to have a different mass, and the mass of the Lagrange variable is updated according to the phase field parameters. To achieve high compatibility with the MPM solver, the mass is not updated during the transfer process. After one time step, the mass is updated, and the updated mass satisfies the following formula:

[0138]

[0139] in, Indicates the Lagrange particle at the next moment quality express Phase mass fraction, express Phase quality.

Claims

1. A multiphase flow interaction simulation method based on MPM and phase field model, characterized in that: Includes the following steps: Step S1: Transfer the mass parameters, phase field parameters, and momentum parameters stored in the Lagrange particles to the Eulerian grid, respectively. The specific process of step S1 is as follows: Step S11: Transfer the mass parameters stored in the Lagrangian particles to the Eulerian grid: Convert the Lagrangian mass to the Eulerian grid mass using APIC as the transfer method. This process satisfies the following formula: ， in, Represents Euler nodes quality Represents the Lagrange particle quality These are the weights in the continuous quadratic B-spline interpolation function of the Lagrange particle and Eulerian grid mapping; Step S12: Transfer the momentum parameters stored in the Lagrange particles to the Eulerian grid: The transfer process satisfies momentum conservation, using APIC as the transfer method, and the process satisfies the following formula: ， in, For Lagrange particles affine velocity, Represents the Lagrange particle volume, Euler nodes Location, For Lagrange particles Location; Step S13: Transfer the phase field parameters stored in the Lagrange particles to the Eulerian grid: Use mass fraction as the phase field parameter to reflect the compatibility of the multiphase fluid. Interactive systems only need to solve The phase field equations, and the transfer process only requires the transfer of... For each phase field parameter, selecting the MPM transfer method, the phase field parameter transfer satisfies the following formula: ， in, This represents the phase field parameters of the Eulerian grid. Represents the phase field parameters of a Lagrange particle; Step S2: Update the Eulerian grid phase field by controlling the energy of the entire system based on the chemical potential; Step S3: Update the momentum of the Euler mesh based on gravity, elastic force, and surface tension in the updated equivalent phase field form; Step S4: After updating the mesh momentum, apply slip boundary conditions to the boundary on the Eulerian mesh; After updating the Eulerian mesh momentum in step S4, the boundary conditions are applied to the Eulerian mesh using slip boundary conditions: ， in, This indicates the mesh velocity where the collision has not yet ended. Collision processing is applied to Lagrange particle pairs multiple times before updating their positions. Represents the surface normal; Step S5: After updating the phase field and momentum, convert the momentum parameters and phase field parameters from Eulerian grid variables to Lagrangian variables; The specific process of step S5 is as follows: Step S51: After updating the momentum parameters on the Eulerian grid, transfer the momentum parameters back from the Eulerian grid to the Lagrange particle: Using APIC as the transfer method, the momentum parameter transfer satisfies: in, Represents the Lagrange particle speed, Represents the Euler mesh speed; Step S52: After updating the phase field parameters on the Eulerian grid, transfer the phase field parameters back from the Eulerian grid to the Lagrangian particles. The phase field information transferred back to the Lagrangian particles is the mass fraction of each phase component. This process follows the formula below: Step S6: Update the state of the Lagrange particle: Update the position and mass state of the Lagrange particle based on the information transmitted back.

2. The multiphase flow interaction simulation method based on MPM and phase field model according to claim 1, characterized in that: The process of updating the Eulerian grid phase field based on the chemical potential in step S2 is as follows: Step S21: Solve for chemical potential: Chemical potential characterizes the transfer trend between phase components and controls the development of the entire system from high chemical potential to low chemical potential. After interpolating the phase field parameters to the Eulerian grid, the central grid and central difference are used to discretize them. The governing equations for the phase field are: ， in, Represents the spatial region of the system. It is a constant used to control the diffusion rate. This represents the interfacial energy, used to reduce the gradient of the mass fraction field, making surface changes smoother. It represents physical capabilities and controls how fluids interact. This represents the phase field value at position x at time t; The governing equation for the chemical potential is: ， Represents chemical potential, Indicates physical fitness. It is a constant. Represents phase field parameters, Represents the Laplace operator; Step S22: Selecting the volumetric energy function: The local minimum value of the volumetric energy function is used to control the low-energy state of the fluid. The process of the fluid evolving from a high-energy state to a low-energy state is simulated. Different volumetric energy functions are designed to achieve different interactive effects. Two-phase interaction is displayed using a volumetric energy function graph, and the minimum value in the effective region is marked. Three-phase interaction uses a Gibbs triangle as the energy display. The design of the volumetric energy function follows the following principles: the minimum energy value is located within the effective region to ensure that the energy minimization state is ultimately achieved; the partial derivative of the chemical potential with respect to a certain phase should only be related to the mass ratio of that phase. Step S23: Update the phase field: After designing the corresponding bulk energy function for different multiphase flow interactions, the phase field is updated using the chemical potential. The degenerate mobility is set as a constant, and the phase field is discretized using a central grid and central difference. The phase field update satisfies the following formula: ， in, Represents the mass derivative. Indicates degenerate mobility. It represents chemical potential.

3. The multiphase flow interaction simulation method based on MPM and phase field model according to claim 2, characterized in that: In step S3, the momentum of the Euler mesh is updated based on gravity and elastic force according to the following formula: ， in, Indicates elastic force. This represents the grid spacing of the Euler mesh. Let be the initial volume of the Lagrange particle. Indicates Piolakirchhoff stress, express Lagrange particles at time Deformation gradient, Represents the Lagrange particle The transpose of the deformation gradient matrix, where T denotes the transpose. Represents gravity. Represents gravitational acceleration. Represents the Euler mesh The sum of forces acting on the body Indicates the current moment. Indicates the next moment, Indicates a time step, Represents the Euler mesh The speed in the next moment, Represents the Euler mesh The speed at the current moment.

4. The multiphase flow interaction simulation method based on MPM and phase field model according to claim 3, characterized in that: In step S3, the updated phase field is used to capture the interface between different components, surface tension is applied, and momentum is further updated. The surface tension governing equation satisfies: ， in, It uses the equivalent surface tension in the form of phase field. express , The surface tension coefficient between phases express Phase field parameters of the phase, express Phase field parameters of the phase, Coefficients related to the diffusion interface.

5. The multiphase flow interaction simulation method based on MPM and phase field model according to claim 1, characterized in that: The specific process of step S6 is as follows: Step S61: After transmitting the momentum information back to the Lagrange particle, update the position of the Lagrange particle based on the momentum information. The updated position satisfies: ， in, Represents the Lagrange particle at the current moment. Location, Indicates the Lagrange particle at the next moment Location, Indicates a time step, This indicates the velocity of the Lagrange particle at the next instant; Step S62: Update the Lagrange particle mass: Update the mass of the Lagrange variables according to the phase field parameters. The mass is not updated during the transfer process. After one time step, the mass is updated, and the updated mass satisfies the following formula: in, Indicates the Lagrange particle at the next moment quality express Phase mass fraction, express Phase quality.

Images