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

Multi-component solid and fluid simulation method and system based on SPH (Smoothed Particle Hydrodynamics) method Download PDF

Info

Publication number
CN105956262A
CN105956262A CN201610280000.0A CN201610280000A CN105956262A CN 105956262 A CN105956262 A CN 105956262A CN 201610280000 A CN201610280000 A CN 201610280000A CN 105956262 A CN105956262 A CN 105956262A
Authority
CN
China
Prior art keywords
solid
stress tensor
particle
component
fluid
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201610280000.0A
Other languages
Chinese (zh)
Other versions
CN105956262B (en
Inventor
胡事民
严枭
江云涛
李晨锋
拉尔夫·罗伯特·马丁
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Tsinghua University
Original Assignee
Tsinghua University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Tsinghua University filed Critical Tsinghua University
Priority to CN201610280000.0A priority Critical patent/CN105956262B/en
Publication of CN105956262A publication Critical patent/CN105956262A/en
Application granted granted Critical
Publication of CN105956262B publication Critical patent/CN105956262B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/30—Circuit design
    • G06F30/36—Circuit design at the analogue level
    • G06F30/367—Design verification, e.g. using simulation, simulation program with integrated circuit emphasis [SPICE], direct methods or relaxation methods

Landscapes

  • Engineering & Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Microelectronics & Electronic Packaging (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开一种基于SPH方法的多组分固体和流体模拟方法及系统，能模拟固体，多组分流体，固体和多组分流体的耦合，多孔介质流以及溶解。方法包括：将多组分流体的动量守恒公式中加入固体的偏应力张量的散度部分，保留压强部分，计算每个粒子的整体压强；利用标准的多组分流体计算方法更新体积分数，根据空间中固体和流体的体积分数以及相对速度计算出每个粒子所有组分的速度梯度；根据所述速度梯度利用固体的本构方程计算固体的偏应力张量，并根据屈服准则修正所述偏应力张量；利用固体的体积分数修正所述偏应力张量；计算不同粒子之间的粘滞力大小；根据所述偏应力张量和粘滞力大小，利用所述动量守恒公式进行多组分物理模拟。

The invention discloses a multi-component solid and fluid simulation method and system based on the SPH method, capable of simulating solid, multi-component fluid, coupling of solid and multi-component fluid, porous medium flow and dissolution. The method includes: adding the divergence part of the deviatoric stress tensor of the solid to the momentum conservation formula of the multi-component fluid, retaining the pressure part, and calculating the overall pressure of each particle; using the standard multi-component fluid calculation method to update the volume fraction, The velocity gradient of all components of each particle is calculated according to the volume fraction and relative velocity of solid and fluid in the space; the deviatoric stress tensor of the solid is calculated using the constitutive equation of the solid according to the velocity gradient, and the yield criterion is corrected Deviating stress tensor; using the volume fraction of solids to correct the deviating stress tensor; calculating the size of the viscous force between different particles; according to the deviating stress tensor and the size of the viscous force, using the momentum conservation formula to perform multiple Component physics simulation.

Description

基于SPH方法的多组分固体和流体模拟方法及系统Multi-component solid and fluid simulation method and system based on SPH method

技术领域technical field

本发明涉及计算机图形学物理模拟与渲染技术领域，具体涉及一种基于SPH方法的多组分固体和流体模拟方法及系统。The invention relates to the technical field of computer graphics physical simulation and rendering, in particular to a multi-component solid and fluid simulation method and system based on the SPH method.

背景技术Background technique

近十年来，涉及流体和固体的基于物理的模拟在图形学界越来越流行。而流体和固体之间的相互作用也是图形学领域的一个值得研究的难点之一。Physics-based simulations involving fluids and solids have grown in popularity in the graphics community over the past decade. And the interaction between fluid and solid is also one of the difficulties worth studying in the field of graphics.

之前的工作涉及流体和固体的相互耦合，如N.Akinci等人发表的论文“Versatilerigid-fluid coupling for incompressible SPH.”处理了流体和刚体的相互作用，之后他们发表的“Coupling elastic solids with SPH fluids.”将方法扩展到弹性体。T.Lenaerts等人发表论文“Mixing fluids and granular materials.”处理了沙子和水之间的多孔介质流问题。他们的方法只能对应特定的场景而没有一个统一的框架来处理这些所有的情况。此外，对于溶解等需要体积分数或质量分数来表述的效果也很难进行模拟。The previous work involved the mutual coupling of fluids and solids, such as the paper "Versatilerigid-fluid coupling for incompressible SPH." published by N.Akinci et al. dealt with the interaction of fluids and rigid bodies, and then they published "Coupling elastic solids with SPH fluids .” extends the method to elastic bodies. T. Lenaerts et al published a paper "Mixing fluids and granular materials." dealing with the flow of porous media between sand and water. Their methods can only correspond to specific scenarios without a unified framework to handle all these situations. In addition, it is difficult to simulate effects such as dissolution that need to be expressed by volume fraction or mass fraction.

对于多组分流体，B.Ren等人发表的论文“Multiple-fluid SPH simulationusing a mixture model”利用混合模型，可以处理可混与不可混的流体的模拟。T.Yang等人利用能量的方法可以模拟类似萃取，特征匹配以及流动性混合等效果。但是他们只考虑了流体，对于固体的模拟他们都没有涉及。For multi-component fluids, the paper "Multiple-fluid SPH simulation using a mixture model" published by B.Ren et al. uses a mixed model to handle the simulation of mixed and immiscible fluids. T. Yang et al. use energy methods to simulate effects such as extraction, feature matching, and fluidity mixing. But they only considered fluids, and they didn't involve in the simulation of solids.

总的来说，之前的方法要么场景过于特定，要么虽然能处理可混与不可混之间的流体，但是没有扩展到固体的情况。In general, the previous methods are either too scene-specific, or although they can handle fluids between mixable and immiscible, they have not been extended to solids.

发明内容Contents of the invention

本发明需要解决的技术问题是，如何模拟固体和流体之间的相互作用的模拟，主要涉及的固体有弹塑性体和颗粒物质。此外，本技术保持了一个统一的框架，使得固体、多组分流体、固体和多组分流体的耦合、多孔介质流、溶解的模拟，以及固体和流体组分之间的可混与不可混的模拟在一个统一的框架下都可以进行，具有鲁棒性，且易于实现。The technical problem to be solved in the present invention is how to simulate the interaction between solid and fluid, mainly related to solid, elastoplastic and granular matter. Furthermore, the technique maintains a unified framework that enables the simulation of solids, multicomponent fluids, coupling of solids and multicomponent fluids, porous media flow, dissolution, and miscibility and immiscibility between solid and fluid components All simulations can be carried out under a unified framework, which is robust and easy to implement.

一方面，本发明实施例提出一种基于SPH方法的多组分固体和流体模拟方法，包括：On the one hand, the embodiment of the present invention proposes a multi-component solid and fluid simulation method based on the SPH method, including:

S1、将多组分流体的动量守恒公式中加入固体的偏应力张量的散度部分，保留压强部分，并计算每个粒子的整体压强；S1. Add the divergence part of the deviatoric stress tensor of the solid to the momentum conservation formula of the multi-component fluid, retain the pressure part, and calculate the overall pressure of each particle;

S2、利用标准的多组分流体计算方法更新体积分数，根据空间中固体和流体的体积分数以及相对速度计算出每个粒子所有组分的速度梯度；S2. Use the standard multi-component fluid calculation method to update the volume fraction, and calculate the velocity gradient of all components of each particle according to the volume fraction and relative velocity of solid and fluid in the space;

S3、根据所述速度梯度利用固体的本构方程计算固体的偏应力张量，并根据屈服准则修正所述偏应力张量；S3. Calculate the deviatoric stress tensor of the solid by using the constitutive equation of the solid according to the velocity gradient, and correct the deviatoric stress tensor according to the yield criterion;

S4、利用固体的体积分数修正所述偏应力张量；S4. Using the volume fraction of solids to correct the deviatoric stress tensor;

S5、计算不同粒子之间的粘滞力大小；S5, calculating the size of the viscous force between different particles;

S6、根据所述偏应力张量和粘滞力大小，利用所述动量守恒公式进行多组分物理模拟。S6. According to the magnitude of the deviatoric stress tensor and the viscous force, use the momentum conservation formula to perform multi-component physical simulation.

另一方面，本发明实施例提出一种基于SPH方法的多组分固体和流体模拟系统，包括：On the other hand, an embodiment of the present invention proposes a multi-component solid and fluid simulation system based on the SPH method, including:

第一计算单元，用于将多组分流体的动量守恒公式中加入固体的偏应力张量的散度部分，保留压强部分，并计算每个粒子的整体压强；The first calculation unit is used to add the divergence part of the deviatoric stress tensor of the solid to the momentum conservation formula of the multi-component fluid, retain the pressure part, and calculate the overall pressure of each particle;

第二计算单元，用于利用标准的多组分流体计算方法更新体积分数，根据空间中固体和流体的体积分数以及相对速度计算出每个粒子所有组分的速度梯度；The second calculation unit is used to update the volume fraction using the standard multi-component fluid calculation method, and calculate the velocity gradient of all components of each particle according to the volume fraction and relative velocity of solid and fluid in the space;

第三计算单元，用于根据所述速度梯度利用固体的本构方程计算固体的偏应力张量，并根据屈服准则修正所述偏应力张量；The third calculation unit is used to calculate the deviatoric stress tensor of the solid by using the constitutive equation of the solid according to the velocity gradient, and correct the deviatoric stress tensor according to the yield criterion;

第一修正单元，用于利用固体的体积分数修正所述偏应力张量；a first correction unit, for correcting the deviatoric stress tensor by using the volume fraction of solid;

第二修正单元，用于计算不同粒子之间的粘滞力大小；The second correction unit is used to calculate the size of the viscous force between different particles;

模拟单元，用于根据所述偏应力张量和粘滞力大小，利用所述动量守恒公式进行多组分物理模拟。The simulation unit is used to perform multi-component physical simulation by using the momentum conservation formula according to the magnitude of the deviatoric stress tensor and the viscous force.

本发明实施例提供的基于SPH方法的多组分固体和流体模拟方法及系统，基于物理模型中的质量守恒、动量守恒和固体的本构，可以模真实世界中的固体，多组分流体以及多组分流体和固体之间的相互作用，在模拟固体时只需要设置粒子中固体的体积分数为1即可，而多组分流体和固体之间的相互作用的模拟包括固体和多组分流体的耦合、多孔介质流的模拟，以及固体和流体组分之间的可混与不可混的模拟。对于溶解的模拟，只需要修改扩散方程中的基于体积分数的扩散部分，约束使体积分数不超过饱和浓度即可。另外，本发明因为只涉及额外的偏应力张量，所以额外的开销也可以保证不会过大。The multi-component solid and fluid simulation method and system based on the SPH method provided by the embodiments of the present invention can simulate solids in the real world, multi-component fluids and For the interaction between multi-component fluids and solids, when simulating solids, you only need to set the volume fraction of solids in particles to 1, while the simulation of interactions between multi-component fluids and solids includes solids and multi-components Coupling of fluids, simulation of flow in porous media, and simulation of miscibility and immiscibility between solid and fluid components. For the simulation of dissolution, it is only necessary to modify the volume fraction-based diffusion part of the diffusion equation, and constrain the volume fraction not to exceed the saturation concentration. In addition, since the present invention only involves the additional deviatoric stress tensor, the additional overhead can also be guaranteed not to be too large.

附图说明Description of drawings

图1为本发明基于SPH方法的多组分固体和流体模拟方法一实施例的流程示意图；Fig. 1 is a schematic flow sheet of an embodiment of the multi-component solid and fluid simulation method based on the SPH method of the present invention;

图2为本发明基于SPH方法的多组分固体和流体模拟系统一实施例的结构示意图。Fig. 2 is a structural schematic diagram of an embodiment of the multi-component solid and fluid simulation system based on the SPH method of the present invention.

具体实施方式detailed description

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。In order to make the purpose, technical solutions and advantages of the embodiments of the present invention clearer, the technical solutions in the embodiments of the present invention will be clearly described below in conjunction with the accompanying drawings in the embodiments of the present invention. Obviously, the described embodiments are the Some, but not all, embodiments are invented. Based on the embodiments of the present invention, all other embodiments obtained by persons of ordinary skill in the art without making creative efforts belong to the protection scope of the present invention.

如图1所示，本实施例公开一种基于SPH方法的多组分固体和流体模拟方法，包括：As shown in Figure 1, this embodiment discloses a multi-component solid and fluid simulation method based on the SPH method, including:

S1、将多组分流体的动量守恒公式中加入固体的偏应力张量的散度部分，保留压强部分，并计算每个粒子的整体压强；S1. Add the divergence part of the deviatoric stress tensor of the solid to the momentum conservation formula of the multi-component fluid, retain the pressure part, and calculate the overall pressure of each particle;

S2、利用标准的多组分流体计算方法更新体积分数，根据空间中固体和流体的体积分数以及相对速度计算出每个粒子所有组分的速度梯度；S2. Use the standard multi-component fluid calculation method to update the volume fraction, and calculate the velocity gradient of all components of each particle according to the volume fraction and relative velocity of solid and fluid in the space;

S3、根据所述速度梯度利用固体的本构方程计算固体的偏应力张量，并根据屈服准则修正所述偏应力张量；S3. Calculate the deviatoric stress tensor of the solid by using the constitutive equation of the solid according to the velocity gradient, and correct the deviatoric stress tensor according to the yield criterion;

S4、利用固体的体积分数修正所述偏应力张量；S4. Using the volume fraction of solids to correct the deviatoric stress tensor;

S5、计算不同粒子之间的粘滞力大小；S5, calculating the size of the viscous force between different particles;

S6、根据所述偏应力张量和粘滞力大小，利用所述动量守恒公式进行多组分物理模拟。S6. According to the magnitude of the deviatoric stress tensor and the viscous force, use the momentum conservation formula to perform multi-component physical simulation.

本发明首先针对混合流体的Navier-Stokes方程的动量守恒方程进行修改：The present invention first modifies the momentum conservation equation of the Navier-Stokes equation of the mixed fluid:

DD. (( ρρ mm uu mm )) DD. (( tt )) == -- ▿▿ pp ++ ▿▿ ·&Center Dot; (( ττ mm ++ ττ DD. mm )) ++ ρρ mm gg -- -- -- (( 11 ))

公式(1)是多组分流体的动量守恒方程，其中D/D(t)表示随体导数，ρm为组分的混合密度，um为混合速度，t为时间，p为压强，τm为黏滞力张量，τDm为相对速度张量，g为重力加速度。为了可以扩展到固体的模拟，同时为了保证多组分的压强的连续性，本发明保留了压强分量，在等式右边加入一个类似黏滞张量的混合固体偏应力张量τSm，使得在计算各个组分的相互作用的时候考虑内部的固体特性。这样整个动量守恒方程就变成了：Equation (1) is the momentum conservation equation of a multi-component fluid, where D/D(t) represents the derivative with the body, ρ m is the mixing density of the components, u m is the mixing velocity, t is the time, p is the pressure, τ m is the viscous force tensor, τ Dm is the relative velocity tensor, and g is the gravitational acceleration. In order to be able to extend to the simulation of solids, and to ensure the continuity of multi-component pressures, the present invention retains the pressure component, and adds a mixed solid deviatoric stress tensor τ Sm similar to the viscous tensor to the right side of the equation, so that in Internal solid properties are taken into account when calculating the interactions of the individual components. The entire momentum conservation equation then becomes:

DD. (( ρρ mm uu mm )) DD. (( tt )) == -- ▿▿ pp ++ ▿▿ ·&Center Dot; (( ττ SS mm ++ ττ mm ++ ττ DD. mm )) ++ ρρ mm gg -- -- -- (( 22 ))

对于这个混合固体偏应力张量τSm的计算，本发明针对弹塑性体和用来模拟颗粒物质的模型——亚塑性体模型的本构方程扩展并计算相应的偏应力张量。For the calculation of the deviatoric stress tensor τ Sm of the mixed solid, the present invention expands and calculates the corresponding deviatoric stress tensor for the constitutive equation of the elastoplastic body and the hypoplastic body model, which is a model for simulating granular matter.

弹塑性体根据广义胡克定律得到相应的偏应力张量：According to the generalized Hooke's law, the elastic-plastic body obtains the corresponding deviatoric stress tensor:

σσ ·&Center Dot; ′′ == ωω ·· ·&Center Dot; σσ ′′ -- σσ ′′ ·· ωω ·&Center Dot; ++ 22 GG (( ϵϵ ·· ′′ -- ϵϵ ·· pp ′′ )) -- -- -- (( 33 ))

其中σ是应力张量，G是剪切模量，ω为旋转张量，ε为应变张量，εp为塑性应变张量，符号表示a对时间求导，符号“b′”表示张量b的偏应力部分。ε和ω这两个张量都可以通过速度梯度直接计算出来。对于塑性应力的部分，本方法直接采用von Mises准则对偏应力张量进行修正：where σ is the stress tensor, G is the shear modulus, ω is the rotation tensor, ε is the strain tensor, ε p is the plastic strain tensor, symbol Indicates the derivation of a with respect to time, and the symbol "b'" indicates the deviatoric stress part of the tensor b. Both tensors ε and ω can be directly computed from the velocity gradient. For the part of plastic stress, this method directly uses the von Mises criterion to correct the deviatoric stress tensor:

σ′:＝σ′/Y (4)σ′:=σ′/Y (4)

其中Y是屈服半径。这样对于弹塑性体的偏应力部分，首先计算其弹性部分，再根据von Mises准则进行修正从准确的得到弹塑性体的偏应力张量。where Y is the yield radius. In this way, for the deviatoric stress part of the elastoplastic body, the elastic part is calculated first, and then corrected according to the von Mises criterion to obtain the deviatoric stress tensor of the elastoplastic body accurately.

亚塑性体根据亚塑性本构得到相应的偏应力张量：The hypoplastic body obtains the corresponding deviatoric stress tensor according to the hypoplastic constitutive:

其中，c1、c2和c3是人为设定的参数以调整颗粒物质的性质，Tr()表示张量的迹， 表示Jaumann导数。本发明利用Drucker-Prager准则对偏应力张量进行修正：Among them, c 1 , c 2 and c 3 are artificially set parameters to adjust the properties of granular matter, Tr() represents the trace of the tensor, Denotes the Jaumann derivative. The present invention utilizes the Drucker-Prager criterion to correct the deviatoric stress tensor:

σσ ′′ :: == (( αα φφ pp ++ kk cc )) σσ ′′ // JJ 22 (( σσ ′′ )) -- -- -- (( 66 ))

其中，αφ和kc为人为设定的参数，用来调整屈服的性质。为σ′的第二主不变量。Among them, α φ and k c are artificially set parameters, which are used to adjust the yield properties. is the second main invariant of σ′.

由于弹塑性体和亚塑性体的本构都需要利用到速度梯度，对于可混的SPH粒子，本方法利用平均的方法来计算每个组分的速度梯度：Since both elastoplastic and hypoplastic constitutives need to use the velocity gradient, for the mixed SPH particles, this method uses the average method to calculate the velocity gradient of each component:

▿▿ (( uu kk )) ii == (( ΣΣ jj mm jj ρρ jj 22 (( αα kk )) ii (( αα kk )) jj (( αα kk )) ii ++ (( αα kk )) jj ×× -- -- -- (( 77 ))

(( (( uu mkmk )) jj -- (( uu mkmk )) ii ++ (( uu mm )) jj -- (( uu mm )) ii )) ▿▿ WW ijij

其中，i和j表示粒子下标，uk为组分k的绝对速度，αk为组分k的体积分数，mj为质量，ρj表示粒子j的密度，Wij为SPH的核函数，umk为组分相对于粒子中心的相对速度，um为粒子中心的速度，uk＝umk+um。where i and j represent particle subscripts, u k is the absolute velocity of component k, α k is the volume fraction of component k, m j is the mass, ρ j represents the density of particle j, W ij is the kernel function of SPH , u mk is the relative velocity of the component relative to the particle center, u m is the velocity of the particle center, u k = umk +u m .

为了稳定性，本方法使用移动最小二乘法进行离散：For stability, this method uses the moving least squares method for discretization:

▿▿ (( uu kk )) ii == (( ΣΣ jj mm jj ρρ jj 22 (( αα kk )) ii (( αα kk )) jj (( αα kk )) ii ++ (( αα kk )) jj (( (( uu kk )) jj -- (( uu kk )) ii )) xx ii jj WW ii jj )) (( ΣΣ jj mm jj ρρ jj 22 (( αα kk )) ii (( αα kk )) jj (( αα kk )) ii ++ (( αα kk )) jj xx ii jj xx ii jj WW ii jj )) -- 11 -- -- -- (( 88 ))

xij为粒子i和粒子j的距离向量，这样计算得到的不同组分的偏应力张量只适合独立组分的计算，在可混的情况下，计算混合偏应力张量的散度的时候，本方法同样需要根据体积分数进行平均化：x ij is the distance vector between particle i and particle j. The deviatoric stress tensors of different components calculated in this way are only suitable for the calculation of independent components. In the case of mixing, when calculating the divergence of the mixed deviatoric stress tensor , this method also needs to average according to the volume fraction:

Ff ii ττ SS mm == mm ii ΣΣ kk ΣΣ jj mm jj ρρ jj 22 (( αα kk )) ii (( αα kk )) jj (( αα kk )) ii ++ (( αα kk )) jj (( (( σσ kk )) ii ′′ ρρ ii 22 ++ (( σσ kk )) jj ′′ ρρ jj 22 )) ·&Center Dot; ▿▿ WW ii jj -- -- -- (( 99 ))

其中k遍历所有可能的组分，表示混合偏引力张量所对应的力，(σk)′表示第k个组分算出来的偏应力张量。这样就计算得到了混合的偏应力张量的散度(本质上就是这个偏应力张量的力形式：)。where k iterates over all possible components, denotes the force corresponding to the mixed deviatoric gravitational tensor, and (σ k )' denotes the deviatoric stress tensor calculated from the kth component. In this way, the divergence of the mixed deviatoric stress tensor is calculated (essentially the force form of this deviatoric stress tensor: ).

对于整体的黏滞力计算，本方法针对可混组分和不可混组分进行单独分析，如果两种物质是不可混的，那么需要加强黏滞力防止穿透，如果是可混的，利用平均的方法平均黏滞力：For the calculation of the overall viscous force, this method analyzes the miscible and immiscible components separately. If the two substances are immiscible, the viscous force needs to be strengthened to prevent penetration. If they are miscible, use Averaging method Average Viscosity:

Πij为i粒子和j粒子的相互黏滞力，γb为不可混情况的粘滞系数，γk为可混情况的组分k的粘滞系数，πij为黏滞变量，它和两个粒子之间的距离成反比，和两个粒子之间的相对速度成正比。Π ij is the mutual viscous force between particle i and particle j, γ b is the viscosity coefficient in the immiscible case, γ k is the viscosity coefficient of component k in the mixable case, and π ij is the viscosity variable, which is related to the two The distance between two particles is inversely proportional to the relative velocity between two particles.

对于溶解的模拟，本方法修改了扩散方程，使得扩散方程不会使得溶质的浓度超过饱和浓度，并且对于速度梯度的加权平均的方式，本方法也设置为当固体的体积分数小于饱和浓度的时候，其加权系数成为零：For the simulation of dissolution, this method modifies the diffusion equation so that the diffusion equation will not cause the concentration of solute to exceed the saturation concentration, and for the weighted average of the velocity gradient, this method is also set to when the volume fraction of solids is less than the saturation concentration , whose weighting coefficient becomes zero:

μμ kk (( ▿▿ ββ kk ββ kk -- ΣΣ kk ′′ ρρ kk ′′ ρρ mm ▿▿ ββ kk ′′ )) -- -- -- (( 1010 ))

表达式(10)是对扩散速度的修改，其中，μk为组分k的扩散系数，βk组分k的体积分数和饱和浓度的最小者，即βk＝(αs,αk)，其中αs为饱和浓度。Expression (10) is a modification of the diffusion velocity, where μ k is the diffusion coefficient of component k, and β k is the minimum of the volume fraction of component k and the saturation concentration, that is, β k = (α s ,α k ) , where α s is the saturation concentration.

对于用于速度梯度的体积分数的平均，本方法将公式(9)中的修改为：For the averaging of volume fractions used for velocity gradients, this method takes the change into:

22 (( (( αα kk )) ii -- αα sthe s )) (( (( αα kk )) jj -- αα sthe s )) (( (( αα kk )) ii -- αα sthe s )) ++ (( (( αα kk )) jj -- αα sthe s )) -- -- -- (( 1111 ))

如果这个表达式中(αk)i-αs或者(αk)j-αs小于零，那么就另表达式(11)变成零。If (α k ) i -α s or (α k ) j -α s in this expression is less than zero, then the other expression (11) becomes zero.

本发明根据物理学的守恒方程得到，可以在一个统一的框架下模拟弹塑性体，颗粒物质，多组分流体，固体和多组分流体的耦合，多孔介质流体以及溶解等现象。本方法可以扩展多个组分并且保证了鲁棒性并且不需要非常大的开销。The invention is obtained according to the conservation equation of physics, and can simulate elastic-plastic bodies, granular materials, multi-component fluids, coupling of solids and multi-component fluids, porous media fluids, and dissolution under a unified framework. This method can be extended to multiple components and guarantees robustness without very large overhead.

本发明实施例提供的基于SPH方法的多组分固体和流体模拟方法，基于物理模型中的质量守恒、动量守恒和固体的本构，可以模真实世界中的固体，多组分流体以及多组分流体和固体之间的相互作用，在模拟固体时只需要设置粒子中固体的体积分数为1即可，而多组分流体和固体之间的相互作用的模拟包括多组分流体的耦合和多孔介质流的模拟。对于溶解的模拟，只需要修改扩散方程中的基于体积分数的扩散部分，约束使体积分数不超过饱和浓度即可。The multi-component solid and fluid simulation method based on the SPH method provided by the embodiment of the present invention, based on the mass conservation, momentum conservation and solid constitutive in the physical model, can simulate solids in the real world, multi-component fluids and multi-group When simulating the interaction between the fluid and the solid, it is only necessary to set the volume fraction of the solid in the particle to 1, while the simulation of the interaction between the multi-component fluid and the solid includes the coupling of the multi-component fluid and Simulation of Flow in Porous Media. For the simulation of dissolution, it is only necessary to modify the volume fraction-based diffusion part of the diffusion equation, and constrain the volume fraction not to exceed the saturation concentration.

可选地，在本发明基于SPH方法的多组分固体和流体模拟方法的另一实施例中，所述S1，包括：Optionally, in another embodiment of the multi-component solid and fluid simulation method based on the SPH method of the present invention, the S1 includes:

通过标准SPH的状态方程计算每个粒子的整体压强。The bulk pressure of each particle is calculated by the equation of state of the standard SPH.

可选地，在本发明基于SPH方法的多组分固体和流体模拟方法的另一实施例中，所述S2，包括：Optionally, in another embodiment of the multi-component solid and fluid simulation method based on the SPH method of the present invention, the S2 includes:

根据SPH粒子中的每个组分相对于粒子总体速度的相对速度加上粒子的总体速度得到的每个组分的绝对速度，以及固体的体积分数利用移动最小二乘法计算组分的速度梯度。The absolute velocity of each component is obtained from the relative velocity of each component in the SPH particle relative to the overall velocity of the particle plus the overall velocity of the particle, and the volume fraction of the solid using the moving least squares method to calculate the velocity gradient of the component.

可选地，在本发明基于SPH方法的多组分固体和流体模拟方法的另一实施例中，所述S3，包括：Optionally, in another embodiment of the multi-component solid and fluid simulation method based on the SPH method of the present invention, the S3 includes:

根据每个组分的速度梯度，利用弹塑性体的本构方程计算得到应力张量，之后求解所述应力张量的偏量部分，得到弹塑性体每个组分的偏应力张量，最后使用von Mises准则对所述弹塑性体每个组分的偏应力张量进行修正；或者According to the velocity gradient of each component, the stress tensor is calculated by using the constitutive equation of the elastoplastic body, and then the deviator part of the stress tensor is solved to obtain the deviatoric stress tensor of each component of the elastoplastic body, and finally Correcting the deviatoric stress tensor for each component of the elastoplastic body using the von Mises criterion; or

根据每个组分的速度梯度，利用亚塑性体的本构方程计算得到应力张量，之后求解所述应力张量的偏量部分，得到颗粒物质每个组分的偏应力张量，最后使用Drucker-Prager准则对所述颗粒物质每个组分的偏应力张量进行修正。According to the velocity gradient of each component, the stress tensor is calculated by using the constitutive equation of the hypoplastic body, and then the deviator part of the stress tensor is solved to obtain the deviatoric stress tensor of each component of the granular matter, and finally use The Drucker-Prager criterion modifies the deviatoric stress tensor for each component of the particulate matter.

可选地，在本发明基于SPH方法的多组分固体和流体模拟方法的另一实施例中，所述S4，包括：Optionally, in another embodiment of the multi-component solid and fluid simulation method based on the SPH method of the present invention, the S4 includes:

对于步骤S3得到的固体独有的偏应力张量，通过体积分数的加权平均得到SPH粒子的偏应力张量。For the unique deviatoric stress tensor of the solid obtained in step S3, the deviatoric stress tensor of the SPH particle is obtained by weighted average of the volume fraction.

可选地，在本发明基于SPH方法的多组分固体和流体模拟方法的另一实施例中，所述S5，包括：Optionally, in another embodiment of the multi-component solid and fluid simulation method based on the SPH method of the present invention, said S5 includes:

对于可混粒子，通过体积分数加权计算粒子之间的摩擦力和黏滞力；或者For miscible particles, compute friction and viscous forces between particles weighted by volume fraction; or

对于不可混粒子，直接使用设定好的粘滞系数计算粒子之间的摩擦力和黏滞力；For immiscible particles, directly use the set viscosity coefficient to calculate the friction and viscosity between particles;

其中，所述S6，包括：Among them, the S6 includes:

根据所述偏应力张量、摩擦力和粘滞力大小进行多组分物理模拟。According to the deviatoric stress tensor, friction force and viscous force, multi-component physical simulation is carried out.

如图2所示，本实施例公开一种基于SPH方法的多组分固体和流体模拟系统，包括：As shown in Figure 2, this embodiment discloses a multi-component solid and fluid simulation system based on the SPH method, including:

第一计算单元1，用于将多组分流体的动量守恒公式中加入固体的偏应力张量的散度部分，保留压强部分，并计算每个粒子的整体压强；The first calculation unit 1 is used to add the divergence part of the deviatoric stress tensor of the solid to the momentum conservation formula of the multi-component fluid, retain the pressure part, and calculate the overall pressure of each particle;

第二计算单元2，用于利用标准的多组分流体计算方法更新体积分数，根据空间中固体和流体的体积分数以及相对速度计算出每个粒子所有组分的速度梯度；The second calculation unit 2 is used to update the volume fraction using a standard multi-component fluid calculation method, and calculate the velocity gradient of all components of each particle according to the volume fraction and relative velocity of solid and fluid in the space;

第三计算单元3，用于根据所述速度梯度利用固体的本构方程计算固体的偏应力张量，并根据屈服准则修正所述偏应力张量；The third calculation unit 3 is used to calculate the deviatoric stress tensor of the solid by using the constitutive equation of the solid according to the velocity gradient, and correct the deviatoric stress tensor according to the yield criterion;

第一修正单元4，用于利用固体的体积分数修正所述偏应力张量；The first correction unit 4 is used to correct the deviatoric stress tensor by using the volume fraction of solids;

第二修正单元5，用于计算不同粒子之间的粘滞力大小；The second correction unit 5 is used to calculate the size of the viscous force between different particles;

模拟单元6，用于根据所述偏应力张量和粘滞力大小，利用所述动量守恒公式进行多组分物理模拟。The simulation unit 6 is used to perform multi-component physical simulation by using the momentum conservation formula according to the magnitude of the deviatoric stress tensor and the viscous force.

虽然结合附图描述了本发明的实施方式，但是本领域技术人员可以在不脱离本发明的精神和范围的情况下做出各种修改和变型，这样的修改和变型均落入由所附权利要求所限定的范围之内。Although the embodiments of the present invention have been described in conjunction with the accompanying drawings, those skilled in the art can make various modifications and variations without departing from the spirit and scope of the present invention. within the bounds of the requirements.

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.
CN201610280000.0A 2016-04-28 2016-04-28 Multi-component solid and fluid simulation method and system based on SPH method Active CN105956262B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201610280000.0A CN105956262B (en) 2016-04-28 2016-04-28 Multi-component solid and fluid simulation method and system based on SPH method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201610280000.0A CN105956262B (en) 2016-04-28 2016-04-28 Multi-component solid and fluid simulation method and system based on SPH method

Publications (2)

Publication Number Publication Date
CN105956262A true CN105956262A (en) 2016-09-21
CN105956262B CN105956262B (en) 2019-08-09

Family

ID=56916885

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201610280000.0A Active CN105956262B (en) 2016-04-28 2016-04-28 Multi-component solid and fluid simulation method and system based on SPH method

Country Status (1)

Country Link
CN (1) CN105956262B (en)

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107133388A (en) * 2017-04-12 2017-09-05 傅德生 A kind of method to cast steel node finite element analysis
CN109271696A (en) * 2018-09-07 2019-01-25 中山大学 Blood clotting analogy method and system based on MPM
CN109726496A (en) * 2019-01-07 2019-05-07 北京航空航天大学 A Realization Method for Improving Incompressible Water Simulation Efficiency Based on IISPH
CN109783935A (en) * 2019-01-15 2019-05-21 北京航空航天大学 A kind of implementation method improving splash fluid stability based on ISPH
CN113486566A (en) * 2021-07-13 2021-10-08 南开大学 Method for performing porous medium and multiphase fluid mixed simulation by using virtual particle phase
CN120102266A (en) * 2025-05-08 2025-06-06 汉中天然谷生物科技股份有限公司 A component detection method suitable for plant extracts

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20100161298A1 (en) * 2008-12-22 2010-06-24 Electronics And Telecommunications Research Institute Method for calculating force acting on interface between immiscible fluids in fluid simulation
CN104200015A (en) * 2014-08-20 2014-12-10 清华大学 Fluid simulating method and device
CN104268943A (en) * 2014-09-28 2015-01-07 北京航空航天大学 Fluid simulation method based on Eulerian-Lagrangian coupling method
CN104573166A (en) * 2013-10-19 2015-04-29 利弗莫尔软件技术公司 Smoothed particle galerkin formulation for simulating physical behaviors in solids mechanics

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20100161298A1 (en) * 2008-12-22 2010-06-24 Electronics And Telecommunications Research Institute Method for calculating force acting on interface between immiscible fluids in fluid simulation
CN104573166A (en) * 2013-10-19 2015-04-29 利弗莫尔软件技术公司 Smoothed particle galerkin formulation for simulating physical behaviors in solids mechanics
CN104200015A (en) * 2014-08-20 2014-12-10 清华大学 Fluid simulating method and device
CN104268943A (en) * 2014-09-28 2015-01-07 北京航空航天大学 Fluid simulation method based on Eulerian-Lagrangian coupling method

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
BO REN 等: "Multiple-fluid SPH Simulation Using a Mixture Model", 《ACM TRANSATIONS ON GRAPHICS》 *

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107133388A (en) * 2017-04-12 2017-09-05 傅德生 A kind of method to cast steel node finite element analysis
CN109271696A (en) * 2018-09-07 2019-01-25 中山大学 Blood clotting analogy method and system based on MPM
CN109726496A (en) * 2019-01-07 2019-05-07 北京航空航天大学 A Realization Method for Improving Incompressible Water Simulation Efficiency Based on IISPH
CN109726496B (en) * 2019-01-07 2020-12-11 北京航空航天大学 A Realization Method for Improving Incompressible Water Simulation Efficiency Based on IISPH
CN109783935A (en) * 2019-01-15 2019-05-21 北京航空航天大学 A kind of implementation method improving splash fluid stability based on ISPH
CN109783935B (en) * 2019-01-15 2020-12-11 北京航空航天大学 An implementation method for improving the stability of splash fluid based on ISPH
CN113486566A (en) * 2021-07-13 2021-10-08 南开大学 Method for performing porous medium and multiphase fluid mixed simulation by using virtual particle phase
CN113486566B (en) * 2021-07-13 2022-05-27 南开大学 A Method for Mixing Simulation of Porous Media and Multiphase Fluids Using Virtual Particle Phases
CN120102266A (en) * 2025-05-08 2025-06-06 汉中天然谷生物科技股份有限公司 A component detection method suitable for plant extracts

Also Published As

Publication number Publication date
CN105956262B (en) 2019-08-09

Similar Documents

Publication Publication Date Title
CN105956262B (en) Multi-component solid and fluid simulation method and system based on SPH method
Bender et al. Position-based simulation of continuous materials
Boedec et al. 3D vesicle dynamics simulations with a linearly triangulated surface
Hwang et al. Direct simulation of particle suspensions in sliding bi-periodic frames
CN104318598B (en) A kind of realization method and system of the solid unidirectional couplings of three-dimensional flow
Liang et al. An efficient phase-field-based multiple-relaxation-time lattice Boltzmann model for three-dimensional multiphase flows
Bogner et al. Simulation of floating bodies with the lattice Boltzmann method
Thömmes et al. A lattice Boltzmann method for immiscible multiphase flow simulations using the level set method
O’Brien et al. A volume-of-fluid ghost-cell immersed boundary method for multiphase flows with contact line dynamics
Hutter et al. Dynamic response of granular and porous materials under large and catastrophic deformations
Kim A diffuse-interface model for axisymmetric immiscible two-phase flow
Fei et al. A multi-scale model for coupling strands with shear-dependent liquid
CN110532727A (en) It can be used for the method for numerical simulation of common non-newtonian fluid
Deng et al. Three-dimensional flow computation with Reynolds stress and algebraic stress models
De Rosis A lattice Boltzmann model for multiphase flows interacting with deformable bodies
Minaki et al. Multiscale modeling and simulation of dynamic wetting
CN109726496B (en) A Realization Method for Improving Incompressible Water Simulation Efficiency Based on IISPH
CN115310339A (en) Solid-liquid coupling simulation method with surface tension effect based on material point method
CN113158531B (en) Single-component and multi-component incompressible fluid simulation method utilizing deformation gradient
Caltagirone et al. A Kinematic Scalar Projection method (KSP) for incompressible flows with variable density
Nan et al. Volume of fluid-discrete element method based the simulation of floating object motion characteristics in complex terrain
Mai-Duy et al. An improved dissipative particle dynamics scheme
Uchiyama Numerical simulation of particle-laden gas flow by vortex in cell method
Fournié et al. A fictitious domain approach for Fluid-Structure Interactions based on the eXtended Finite Element Method.
Zhao et al. Stress recovery procedure for discontinuous deformation analysis

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant