CN105956262A - Multi-component solid and fluid simulation method and system based on SPH (Smoothed Particle Hydrodynamics) method - Google Patents

Abstract

The invention discloses a multi-component solid and fluid simulation method and system based on a SPH (Smoothed Particle Hydrodynamics) method, can simulate solid and multi-component fluid, the coupling of the solid and the multi-component fluid, porous medium flow and dissolving. The method comprises the following steps: adding the divergence part of the deviatoric stress tensor of the solid into the momentum conservation formula of the multi-component fluid, retaining the part of intensity of pressure, and calculating the integral intensity of pressure of each particle; utilizing a standard multi-component fluid calculation method to update a volume fraction, and calculating the velocity gradient of all components of each particle according to the volume fraction of the solid and the fluid in the space and relative speed; according to the volume fraction, utilizing the constitutive equation of the solid to calculate the deviatoric stress tensor of the solid, and correcting the deviatoric stress tensor according to a yield criterion; utilizing the volume fraction of the solid to correct the deviatoric stress tensor; calculating viscous force among different particles; and according to the deviatoric stress tensor and the viscous force, utilizing the momentum conservation formula to carry out multi-component physical simulation.

Description

Multicomponent solid based on SPH method and fluid simulation method and system

Technical field

The present invention relates to computer graphics physical modeling and Rendering field, be specifically related to a kind of based on SPH method Multicomponent solid and fluid simulation method and system.

Background technology

Nearly ten years, the simulation based on physics relating to fluid and solid becomes more and more popular in figure educational circles.And fluid and Interaction between solid is also one of the difficult point being worth research in graphics field.

Work before relates to intercoupling of fluid and solid, the paper " Versatile delivered such as N.Akinci et al. Rigid-fluid coupling for incompressible SPH. " process the interaction of fluid and rigid body, afterwards Method is expanded to elastomer by " Coupling elastic solids with SPH fluids. " that they deliver. T.Lenaerts et al. publish thesis " Mixing fluids and granular materials. " processed sand and water it Between porous media flow problem.Their method can only corresponding specific scene and the unified framework of neither one processes these All of situation.Additionally, the effect needing volume fraction or mass fraction to state for dissolving etc. also is difficult to be simulated.

For multicomponent fluid, paper " the Multiple-fluid SPH simulation that B.Ren et al. delivers Using a mixture model " utilize mixed model, the simulation of the fluid that can mix and can not mix can be processed.T.Yang etc. People utilizes the method for energy can simulate similar extraction, the effect such as characteristic matching and mobility mixing.But they only consider Fluid, for solid simulation they all without reference to.

Generally speaking, method or scene before are the most specific, although or between processing and can mixing and can not mix Fluid, but do not expand to the situation of solid.

Summary of the invention

The technical issues that need to address of the present invention are, how to simulate the simulation of interaction between solid and fluid, main Solid to be related to has elasticoplastic body and particulate matter.Additionally, this technology maintains a unified framework so that solid, many The coupling of component flow, solid and multicomponent fluid, porous media flow, the simulation of dissolving, and between solid and fluid components The simulation mixing with can not mixing can be carried out under a unified framework, there is robustness, and be easily achieved.

On the one hand, the embodiment of the present invention proposes a kind of multicomponent solid based on SPH method and fluid simulation method, bag Include:

S1, the conservation of momentum formula of multicomponent fluid will add the divergence part of the deviatoric stress tensor of solid, retain pressure Strong part, and calculate the overall pressure of each particle；

S2, utilize standard multicomponent fluid computational methods update volume fraction, according to solid in space and the body of fluid Fraction and relative velocity calculate the velocity gradient of each particle all components；

S3, the constitutive equation of solid is utilized to calculate the deviatoric stress tensor of solid, and according to surrender according to described velocity gradient Deviatoric stress tensor described in criterion modification；

S4, utilize deviatoric stress tensor described in the volume fraction correction of solid；

S5, the viscous force size calculated between different particle；

S6, according to described deviatoric stress tensor and viscous force size, utilize described conservation of momentum formula to carry out multicomponent physics Simulation.

On the other hand, the embodiment of the present invention proposes a kind of multicomponent solid based on SPH method and fluid simulation system, bag Include:

First computing unit, for dissipating of the deviatoric stress tensor by adding solid in the conservation of momentum formula of multicomponent fluid Degree part, retains pressure part, and calculates the overall pressure of each particle；

Second computing unit, for utilizing the multicomponent fluid computational methods of standard to update volume fraction, according in space Solid and the volume fraction of fluid and relative velocity calculate the velocity gradient of each particle all components；

3rd computing unit, opens for the deviatoric stress utilizing the constitutive equation of solid to calculate solid according to described velocity gradient Amount, and according to deviatoric stress tensor described in yield criterion correction；

First amending unit, for utilizing deviatoric stress tensor described in the volume fraction correction of solid；

Second amending unit, for calculating the viscous force size between different particle；

Analogue unit, for according to described deviatoric stress tensor and viscous force size, utilizes described conservation of momentum formula to carry out Multicomponent physical modeling.

The multicomponent solid based on SPH method of embodiment of the present invention offer and fluid simulation method and system, based on thing Reason model in the conservation of mass, the conservation of momentum and this structure of solid, can with the solid in mould real world, multicomponent fluid and Interaction between multicomponent fluid and solid, having only to arrange the volume fraction of solid in particle when simulating solid is 1 , and the simulation of the interaction between multicomponent fluid and solid includes that solid and the coupling of multicomponent fluid, porous are situated between The simulation of mass flow, and the simulation mixing Yu can not mixing between solid and fluid components.For the simulation dissolved, it is only necessary to repair Changing the diffusion part based on volume fraction in diffusion equation, constraint makes volume fraction be less than saturated concentration.It addition, this Invention is because pertain only to extra deviatoric stress tensor, so extra expense is it is also ensured that will not be excessive.

Accompanying drawing explanation

Fig. 1 is present invention multicomponent based on SPH method solid and the schematic flow sheet of fluid simulation method one embodiment；

Fig. 2 is present invention multicomponent based on SPH method solid and the structural representation of fluid simulation system one embodiment.

Detailed description of the invention

For making the purpose of the embodiment of the present invention, technical scheme and advantage clearer, below in conjunction with the embodiment of the present invention In accompanying drawing, the technical scheme in the embodiment of the present invention is explicitly described, it is clear that described embodiment is the present invention A part of embodiment rather than whole embodiments.Based on the embodiment in the present invention, those of ordinary skill in the art are not having Make the every other embodiment obtained under creative work premise, broadly fall into the scope of protection of the invention.

As it is shown in figure 1, the present embodiment discloses a kind of multicomponent solid based on SPH method and fluid simulation method, including:

S1, the conservation of momentum formula of multicomponent fluid will add the divergence part of the deviatoric stress tensor of solid, retain pressure Strong part, and calculate the overall pressure of each particle；

S2, utilize standard multicomponent fluid computational methods update volume fraction, according to solid in space and the body of fluid Fraction and relative velocity calculate the velocity gradient of each particle all components；

S3, the constitutive equation of solid is utilized to calculate the deviatoric stress tensor of solid, and according to surrender according to described velocity gradient Deviatoric stress tensor described in criterion modification；

S4, utilize deviatoric stress tensor described in the volume fraction correction of solid；

S5, the viscous force size calculated between different particle；

S6, according to described deviatoric stress tensor and viscous force size, utilize described conservation of momentum formula to carry out multicomponent physics Simulation.

The present invention modifies first against the momentum conservation equation of the Navier-Stokes equation of fluid-mixing:

$\frac{D\left({\mathrm{\ρ}}_m{u}_m\right)}{D\left(t\right)}=-\▿p+\▿\·({\mathrm{\τ}}_m+{\mathrm{\τ}}_{Dm})+{\mathrm{\ρ}}_{m}g---\left(1\right)$

Formula (1) is the momentum conservation equation of multicomponent fluid, and wherein D/D (t) represents material derivative, ρ_mMixing for component Close density, u_mFor mixing velocity, t is the time, and p is pressure, τ_mFor viscous force tensor, τ_DmFor relative velocity tensor, g is that gravity adds Speed.In order to the simulation of solid can be expanded to, simultaneously in order to ensure the seriality of multi-component pressure, present invention preserves pressure Strong component, adds the hybrid solid deviatoric stress tensor τ of a similar viscid tensor on the right of equation_SmSo that calculating each group The solid property when of interaction within consideration divided.The most whole momentum conservation equation reforms into:

$\frac{D\left({\mathrm{\ρ}}_m{u}_m\right)}{D\left(t\right)}=-\▿p+\▿\·({\mathrm{\τ}}_{Sm}+{\mathrm{\τ}}_m+{\mathrm{\τ}}_{Dm})+{\mathrm{\ρ}}_{m}g---\left(2\right)$

For this hybrid solid deviatoric stress tensor τ_SmCalculating, the present invention is directed to elasticoplastic body and for simulating particulate matter The constitutive equation of the model Asia plastic body model of matter extends and calculates corresponding deviatoric stress tensor.

Elasticoplastic body obtains corresponding deviatoric stress tensor according to generalized Hooke law:

${\stackrel{\·}{\mathrm{\σ}}}^{\′}=\stackrel{\·}{\mathrm{\ω}}\·{\mathrm{\σ}}^{\′}-{\mathrm{\σ}}^{\′}\·\stackrel{\·}{\mathrm{\ω}}+2G({\stackrel{\·}{\mathrm{\ϵ}}}^{\′}-{\stackrel{\·}{\mathrm{\ϵ}}}^{p^{\′}})---\left(3\right)$

Wherein σ is stress tensor, and G is modulus of shearing, and ω is rotation tensor, and ε is strain tensor, ε^pOpen for plastic strain Amount, symbolRepresent that a represents the deviatoric stress part of tensor b to time derivation, symbol " b ' ".ε and ω the two tensor all may be used Directly calculate with Negotiation speed gradient.For the part of plastic stress, this method directly uses von Mises criterion to partially Stress tensor is modified:

σ ' :=σ '/Y (4)

Wherein Y is surrender radius.So for the deviatoric stress part of elasticoplastic body, first calculate its elastic part, then root It is modified from the deviatoric stress tensor being accurately obtained elasticoplastic body according to von Mises criterion.

Sub-plastic body obtains corresponding deviatoric stress tensor according to sub-Plastic Constitutive:

Wherein, c₁、c₂And c₃Be the parameter that is manually set to adjust the character of particulate matter, Tr () represents the mark of tensor, Represent Jaumann derivative.The present invention utilizes Drucker-Prager criterion to be modified deviatoric stress tensor:

${\mathrm{\σ}}^{\′}:=({\mathrm{\α}}_{\mathrm{\φ}}p+k_c){\mathrm{\σ}}^{\′}/\sqrt{J_2^{\left({\mathrm{\σ}}^{\′}\right)}}---\left(6\right)$

Wherein, α_φAnd k_cFor the parameter being manually set, it is used for adjusting the character of surrender.For σ ' second is main constant Amount.

Owing to this structure of elasticoplastic body and sub-plastic body is required for using velocity gradient, for the SPH particle that can mix, this Method utilizes average method to calculate the velocity gradient of each component:

$\▿{\left(u_k\right)}_i=(\underset{j}{\Σ}\frac{m_j}{{\mathrm{\ρ}}_j}\frac{2{\left({\mathrm{\α}}_k\right)}_i{\left({\mathrm{\α}}_k\right)}_j}{{\left({\mathrm{\α}}_k\right)}_i+{\left({\mathrm{\α}}_k\right)}_j}\×---\left(7\right)$

$({\left(u_{\mathrm{mk}}\right)}_j-{\left(u_{\mathrm{mk}}\right)}_i+{\left(u_m\right)}_j-{\left(u_m\right)}_i)\▿W_{\mathrm{ij}}$

Wherein, i and j represents particle subscript, u_kFor the absolute velocity of component k, α_kFor the volume fraction of component k, m_jFor matter Amount, ρ_jRepresent the density of particle j, W_ijFor the kernel function of SPH, u_mkFor component relative to the relative velocity of particle centre, u_mFor grain The speed of subcenter, u_k=u_mk+u_m。

For stability, this method uses Moving Least to carry out discrete:

$\begin{array}{c}\▿{\left(u_k\right)}_i=\left(\underset{j}{\Σ}\frac{m_j}{{\mathrm{\ρ}}_j}\frac{2{\left({\mathrm{\α}}_k\right)}_i{\left({\mathrm{\α}}_k\right)}_j}{{\left({\mathrm{\α}}_k\right)}_i+{\left({\mathrm{\α}}_k\right)}_j}\right({\left(u_k\right)}_j-{\left(u_k\right)}_i)x_{ij}{W}_{ij})\\ {\left(\underset{j}{\Σ}\frac{m_j}{{\mathrm{\ρ}}_j}\frac{2{\left({\mathrm{\α}}_k\right)}_i{\left({\mathrm{\α}}_k\right)}_j}{{\left({\mathrm{\α}}_k\right)}_i+{\left({\mathrm{\α}}_k\right)}_j}{x}_{ij}{x}_{ij}{W}_{ij}\right)}^{-1}\end{array}---\left(8\right)$

x_ijFor particle i and the distance vector of particle j, the deviatoric stress tensor of the most calculated different component is only suitable for solely The calculating of vertical component, in the case of can mixing, the when of calculating the divergence of mixing deviatoric stress tensor, this method needs also exist for basis Volume fraction averages:

$F_i^{{\mathrm{\τ}}_{Sm}}=m_i\underset{k}{\Σ}\underset{j}{\Σ}\frac{m_j}{{\mathrm{\ρ}}_j}\frac{2{\left({\mathrm{\α}}_k\right)}_i{\left({\mathrm{\α}}_k\right)}_j}{{\left({\mathrm{\α}}_k\right)}_i+{\left({\mathrm{\α}}_k\right)}_j}(\frac{{\left({\mathrm{\σ}}_k\right)}_i^{\′}}{{\mathrm{\ρ}}_i^2}+\frac{{\left({\mathrm{\σ}}_k\right)}_j^{\′}}{{\mathrm{\ρ}}_j^2})\·\▿W_{ij}---\left(9\right)$

Wherein k travels through all possible component,Represent mixing power corresponding to inclined gravitatonal tensor, (σ_k) ' represent the The deviatoric stress tensor that k component calculates.The divergence of the deviatoric stress tensor being thus calculated mixing (is exactly substantially this The power form of individual deviatoric stress tensor:)。

Overall viscous force is calculated, this method for mixing component and component can not be mixed individually analyze, if Two kinds of materials can not mix, then needs to strengthen viscous force and prevents from penetrating, if can mix, utilize average method average Viscous force:

Π_ijFor i particle and the mutual viscous force of j particle, γ_bFor the coefficient of viscosity of situation, γ can not be mixed_kFor situation can be mixed The coefficient of viscosity of component k, π_ijFor viscid variable, the distance between it and two particles is inversely proportional to, and between two particles Relative velocity is directly proportional.

For the simulation dissolved, this method have modified diffusion equation so that diffusion equation will not make the concentration of solute surpass Supersaturation concentration, and for the average weighted mode of velocity gradient, this method is also configured as when the volume fraction of solid is little In saturated concentration when, its weight coefficient becomes zero:

${\mathrm{\μ}}_k(\frac{\▿{\mathrm{\β}}_k}{{\mathrm{\β}}_k}-\underset{k^{\′}}{\Σ}\frac{{\mathrm{\ρ}}_{k^{\′}}}{{\mathrm{\ρ}}_m}\▿{\mathrm{\β}}_{k^{\′}})---\left(10\right)$

Expression formula (10) is the amendment to diffusion velocity, wherein, μ_kFor the diffusion coefficient of component k, β_kThe volume integral of component k Number and the reckling of saturated concentration, i.e. β_k=(α_s,α_k), wherein α_sFor saturated concentration.

Average for the volume fraction for velocity gradient, this method is by formula (9)Amendment For:

$\frac{2({\left({\mathrm{\α}}_k\right)}_i-{\mathrm{\α}}_s)({\left({\mathrm{\α}}_k\right)}_j-{\mathrm{\α}}_s)}{({\left({\mathrm{\α}}_k\right)}_i-{\mathrm{\α}}_s)+({\left({\mathrm{\α}}_k\right)}_j-{\mathrm{\α}}_s)}---\left(11\right)$

If (α in this expression formula_k)_i-α_sOr (α_k)_j-α_sLess than zero, then become zero with regard to another expression formula (11).

The present invention obtains according to physical conservation equation, can a unified framework Imitating elasticoplastic body, Grain material, multicomponent fluid, solid and the coupling of multicomponent fluid, the phenomenon such as porous media fluid and dissolving.This method can To extend multiple components and ensure that robustness and need not the biggest expense.

The multicomponent solid based on SPH method of embodiment of the present invention offer and fluid simulation method, based on physical model In the conservation of mass, the conservation of momentum and this structure of solid, can be with the solid in mould real world, multicomponent fluid and multicomponent Interaction between fluid and solid, having only to arrange the volume fraction of solid in particle when simulating solid is 1, and The simulation of the interaction between multicomponent fluid and solid includes the coupling of multicomponent fluid and the simulation of porous media flow.Right In the simulation dissolved, it is only necessary to the diffusion part based on volume fraction in amendment diffusion equation, constraint makes volume fraction not surpass Supersaturation concentration.

Alternatively, in another embodiment of present invention multicomponent based on SPH method solid and fluid simulation method, institute State S1, including:

The overall pressure of each particle is calculated by the state equation of standard SPH.

Alternatively, in another embodiment of present invention multicomponent based on SPH method solid and fluid simulation method, institute State S2, including:

Add the overall rate of particle relative to the relative velocity of overall particle speed according to each component in SPH particle The absolute velocity of each component obtained, and the volume fraction of solid utilize Moving Least calculate component speed ladder Degree.

Alternatively, in another embodiment of present invention multicomponent based on SPH method solid and fluid simulation method, institute State S3, including:

According to the velocity gradient of each component, utilize the constitutive equation of elasticoplastic body to be calculated stress tensor, ask afterwards Solve the deviator part of described stress tensor, obtain the deviatoric stress tensor of each component of elasticoplastic body, finally use von Mises accurate Then the deviatoric stress tensor of component each to described elasticoplastic body is modified；Or

According to the velocity gradient of each component, utilize the constitutive equation of sub-plastic body to be calculated stress tensor, ask afterwards Solve the deviator part of described stress tensor, obtain the deviatoric stress tensor of each component of particulate matter, finally use Drucker- The deviatoric stress tensor of Prager criterion component each to described particulate matter is modified.

Alternatively, in another embodiment of present invention multicomponent based on SPH method solid and fluid simulation method, institute State S4, including:

The deviatoric stress tensor that the solid that obtains for step S3 is exclusive, obtains SPH grain by the weighted average of volume fraction The deviatoric stress tensor of son.

Alternatively, in another embodiment of present invention multicomponent based on SPH method solid and fluid simulation method, institute State S5, including:

For particle can be mixed, by the frictional force between volume fraction weighted calculation particle and viscous force；Or

For particle can not be mixed, the coefficient of viscosity set directly is used to calculate the frictional force between particle and viscous force；

Wherein, described S6, including:

Multicomponent physical modeling is carried out according to described deviatoric stress tensor, frictional force and viscous force size.

As in figure 2 it is shown, the present embodiment discloses a kind of multicomponent solid based on SPH method and fluid simulation system, including:

First computing unit 1, for the deviatoric stress tensor by adding solid in the conservation of momentum formula of multicomponent fluid Divergence part, retains pressure part, and calculates the overall pressure of each particle；

Second computing unit 2, for utilizing the multicomponent fluid computational methods of standard to update volume fraction, according in space Solid and the volume fraction of fluid and relative velocity calculate the velocity gradient of each particle all components；

3rd computing unit 3, for utilizing the constitutive equation of solid to calculate the deviatoric stress of solid according to described velocity gradient Tensor, and according to deviatoric stress tensor described in yield criterion correction；

First amending unit 4, for utilizing deviatoric stress tensor described in the volume fraction correction of solid；

Second amending unit 5, for calculating the viscous force size between different particle；

Analogue unit 6, for according to described deviatoric stress tensor and viscous force size, utilizes described conservation of momentum formula to carry out Multicomponent physical modeling.

Although being described in conjunction with the accompanying embodiments of the present invention, but those skilled in the art can be without departing from this Making various modifications and variations in the case of bright spirit and scope, such amendment and modification each fall within by claims Within limited range.

Claims (7)

1. a multicomponent solid based on SPH method and fluid simulation method, it is characterised in that including:

S1, the conservation of momentum formula of multicomponent fluid will add the divergence part of the deviatoric stress tensor of solid, retain pressure portion Point, and calculate the overall pressure of each particle；

S2, utilize standard multicomponent fluid computational methods update volume fraction, according to solid in space and the volume integral of fluid Number and relative velocity calculate the velocity gradient of each particle all components；

S3, the constitutive equation of solid is utilized to calculate the deviatoric stress tensor of solid, and according to yield criterion according to described velocity gradient Revise described deviatoric stress tensor；

S4, utilize deviatoric stress tensor described in the volume fraction correction of solid；

S5, the viscous force size calculated between different particle；

S6, according to described deviatoric stress tensor and viscous force size, utilize described conservation of momentum formula to carry out multicomponent physical modeling.

Multicomponent solid based on SPH method the most according to claim 1 and fluid simulation method, it is characterised in that institute State S1, including:

The overall pressure of each particle is calculated by the state equation of standard SPH.

Multicomponent solid based on SPH method the most according to claim 1 and fluid simulation method, it is characterised in that institute State S2, including:

Add that the overall rate of particle obtains according to each component in SPH particle relative to the relative velocity of overall particle speed The absolute velocity of each component, and the volume fraction of solid utilizes Moving Least to calculate the velocity gradient of component.

Multicomponent solid based on SPH method the most according to claim 1 and fluid simulation method, it is characterised in that institute State S3, including:

According to the velocity gradient of each component, utilize the constitutive equation of elasticoplastic body to be calculated stress tensor, solve institute afterwards State the deviator part of stress tensor, obtain the deviatoric stress tensor of each component of elasticoplastic body, finally use von Mises criterion pair The deviatoric stress tensor of each component of described elasticoplastic body is modified；Or

According to the velocity gradient of each component, utilize the constitutive equation of sub-plastic body to be calculated stress tensor, solve institute afterwards State the deviator part of stress tensor, obtain each component of particulate matter deviatoric stress tensor, finally use Drucker- The deviatoric stress tensor of Prager criterion component each to described particulate matter is modified.

Multicomponent solid based on SPH method the most according to claim 1 and fluid simulation method, it is characterised in that institute State S4, including:

The deviatoric stress tensor that the solid that obtains for step S3 is exclusive, obtains SPH particle by the weighted average of volume fraction Deviatoric stress tensor.

Multicomponent solid based on SPH method the most according to claim 1 and fluid simulation method, it is characterised in that institute State S5, including:

For particle can be mixed, by the frictional force between volume fraction weighted calculation particle and viscous force；Or

For particle can not be mixed, the coefficient of viscosity set directly is used to calculate the frictional force between particle and viscous force；

Wherein, described S6, including:

According to described deviatoric stress tensor, frictional force and viscous force size, described conservation of momentum formula is utilized to carry out multicomponent physics Simulation.

7. a multicomponent solid based on SPH method and fluid simulation system, it is characterised in that including:

First computing unit, the divergence portion of the deviatoric stress tensor for solid will be added in the conservation of momentum formula of multicomponent fluid Point, retain pressure part, and calculate the overall pressure of each particle；

Second computing unit, for utilizing the multicomponent fluid computational methods of standard to update volume fraction, according to solid in space With the velocity gradient that the volume fraction of fluid and relative velocity calculate each particle all components；

3rd computing unit, for utilizing the constitutive equation of solid to calculate the deviatoric stress tensor of solid according to described velocity gradient, And according to deviatoric stress tensor described in yield criterion correction；

First amending unit, for utilizing deviatoric stress tensor described in the volume fraction correction of solid；

Second amending unit, for calculating the viscous force size between different particle；

Analogue unit, for according to described deviatoric stress tensor and viscous force size, utilizes described conservation of momentum formula to carry out many groups Divide physical modeling.