CN109344450B - 基于pbf的流体凝固模拟方法及系统 - Google Patents
基于pbf的流体凝固模拟方法及系统 Download PDFInfo
- Publication number
- CN109344450B CN109344450B CN201811043656.6A CN201811043656A CN109344450B CN 109344450 B CN109344450 B CN 109344450B CN 201811043656 A CN201811043656 A CN 201811043656A CN 109344450 B CN109344450 B CN 109344450B
- Authority
- CN
- China
- Prior art keywords
- particle
- curdled appearance
- calculates
- fluid
- grid
- 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.)
- Active
Links
- 239000012530 fluid Substances 0.000 title claims abstract description 104
- 238000000034 method Methods 0.000 title claims abstract description 79
- 239000002245 particle Substances 0.000 claims abstract description 378
- 238000004088 simulation Methods 0.000 claims abstract description 44
- 238000012937 correction Methods 0.000 claims description 29
- 238000007711 solidification Methods 0.000 claims description 27
- 230000008023 solidification Effects 0.000 claims description 27
- 230000033001 locomotion Effects 0.000 claims description 13
- 238000004422 calculation algorithm Methods 0.000 claims description 9
- 238000006073 displacement reaction Methods 0.000 claims description 9
- 230000000149 penetrating effect Effects 0.000 claims description 9
- 230000001133 acceleration Effects 0.000 claims description 7
- HPTJABJPZMULFH-UHFFFAOYSA-N 12-[(Cyclohexylcarbamoyl)amino]dodecanoic acid Chemical compound OC(=O)CCCCCCCCCCCNC(=O)NC1CCCCC1 HPTJABJPZMULFH-UHFFFAOYSA-N 0.000 claims 2
- 230000035515 penetration Effects 0.000 abstract description 9
- 230000006870 function Effects 0.000 description 32
- 239000007788 liquid Substances 0.000 description 12
- 239000007787 solid Substances 0.000 description 12
- 230000006872 improvement Effects 0.000 description 11
- 230000008569 process Effects 0.000 description 11
- 238000004364 calculation method Methods 0.000 description 10
- 206010053567 Coagulopathies Diseases 0.000 description 8
- 230000035602 clotting Effects 0.000 description 8
- 238000010586 diagram Methods 0.000 description 7
- 238000012545 processing Methods 0.000 description 6
- 230000000694 effects Effects 0.000 description 5
- 230000008901 benefit Effects 0.000 description 4
- 230000008859 change Effects 0.000 description 3
- 230000003068 static effect Effects 0.000 description 3
- YBJHBAHKTGYVGT-ZKWXMUAHSA-N (+)-Biotin Chemical compound N1C(=O)N[C@@H]2[C@H](CCCCC(=O)O)SC[C@@H]21 YBJHBAHKTGYVGT-ZKWXMUAHSA-N 0.000 description 2
- 230000005540 biological transmission Effects 0.000 description 2
- 230000015271 coagulation Effects 0.000 description 2
- 238000005345 coagulation Methods 0.000 description 2
- 239000003086 colorant Substances 0.000 description 2
- 238000001514 detection method Methods 0.000 description 2
- 230000008034 disappearance Effects 0.000 description 2
- 230000005484 gravity Effects 0.000 description 2
- 238000009499 grossing Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000009877 rendering Methods 0.000 description 2
- FEPMHVLSLDOMQC-UHFFFAOYSA-N virginiamycin-S1 Natural products CC1OC(=O)C(C=2C=CC=CC=2)NC(=O)C2CC(=O)CCN2C(=O)C(CC=2C=CC=CC=2)N(C)C(=O)C2CCCN2C(=O)C(CC)NC(=O)C1NC(=O)C1=NC=CC=C1O FEPMHVLSLDOMQC-UHFFFAOYSA-N 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 239000008280 blood Substances 0.000 description 1
- 210000004369 blood Anatomy 0.000 description 1
- 238000005094 computer simulation Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000009795 derivation Methods 0.000 description 1
- 235000012907 honey Nutrition 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 230000002452 interceptive effect Effects 0.000 description 1
- 238000012423 maintenance Methods 0.000 description 1
- 239000003973 paint Substances 0.000 description 1
- 239000006072 paste Substances 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Hardware Design (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
本发明公开了一种基于PBF的流体凝固模拟方法，包括：S1，初始化数据；S2，粒子固定半径邻域搜索；S3，更新粒子的凝固状态；S4，将速度衰减因子作用于粒子，计算粒子的位置中间量；S5，计算粒子密度及运动约束量；S6，计算粒子的位置变化量，并进行防穿透修正；S7，更新粒子速度及位置，并处理粒子与边界碰撞问题；S8，计算人造粘性和涡流约束。本发明还公开了一种基于PBF的流体凝固模拟系统。本发明，可通过不同的粒子速度衰减因子来模拟不同的流体模拟的快慢，还可通过位置变化量修正粒子穿透现象，在GPU的支持下实现更大规模更有效率的流体凝固模拟。
Description
技术领域
本发明涉及一种流体模拟技术领域，尤其涉及一种基于PBF的流体凝固模拟方法及一种基于PBF的流体凝固模拟系统。
背景技术
流体在现实世界中无处不在，如水、蜂蜜、血液、油漆或者浆糊等等。在出现流体的VR虚拟现实场景中，绘制出符合物理运动规律的流体图像是使得VR内容让人有真实感体验的主要环节之一。描述粘性不可压缩流体运动的物理方程是由两位科学家Claude-LouisNavier和George Gabriel Stokes分别归纳推导出的，并以他们的名字命名——Navier-Stokses方程(以下简称NS方程)。在直角坐标系下NS方程具有如下形式：
其中，D(·)/Dt表示物质导数，ρ表示流体的密度，v表示流体的速度，p表示流体的压强，μ表示流体的运动粘性系数，F表示流体所受外力；方程(1)表示的是流体的不可性方程，可以看到不可压缩性等价于流体处于无散度速度场中；方程(2)表示的是流体的动量守恒方程。
由于求解NS方程的精确解是十分困难的，因此现有图形学领域内的流体模拟技术主要是通过某种离散化数值方法来求解NS方程从而计算出流体的物理量，并通过某种可视化方式显示出来。通常有两种不同的方法来离散化流体这种连续介质，分别是欧拉法和拉格朗日法。如图1所示，欧拉法将流体的运动空间划分为不随流动运动所改变的固定均匀网格，在网格上使用有限差分法来计算该点处流体的物理量；如图2所示，拉格朗日法则将流体划分成大小相等的粒子，直接在每个粒子上计算流体的物理量，通过粒子的运动来表现出流体的整体运动。
SPH(Smoothed Particle Hydrodynamics，光滑粒子流体动力学)方法是一种基于拉格朗日法的流体模拟计算。在SPH方法中，流体粒子的物理量(包括密度、压力和粘力等)由其固定半径h内的邻域粒子平滑插值计算得到，插值计算公式如下：
其中，xij表示从粒子i到粒子j的位移，mj和ρj分别表示粒子j的质量和密度，W(xij,h)表示平滑插值函数。相应的，物理量的一次和二次空间导数可以分别写成和因此原先NS方程中难以求解的、非线性的带有哈密顿算子的项和带有拉普拉斯算子的项可以通过易于求解的、线性的加权求和计算来近似求解，在不明显失真的前提下极大地提高了计算效率。VR虚拟场景所要求的，模拟流体的技术要具备计算高效性、效果真实性和模拟实时性，因此传统的SPH方法的优势自然而然地满足了前两个条件。但是由于传统SPH方法存在一个明显的缺陷，就是为了保证数值的稳定性必须严格限制模拟的时间步长，太小的模拟时间步长导致模拟的过程很慢，通俗的理解就是真实流体的1秒钟的运动可能需要计算机模拟计算1分钟才能完成，显然无法满足VR动态场景所要求的实时性。
PBF(Position-Based Fluids，基于位置的流体)方法是传统的SPH方法的一种改进技术，它将PBD(Position-Based Dynamics，基于位置的动力学)框架应用于SPH中，继承了PBD本身较大模拟时间步长的优点，从而可以有效地解决传统SPH方法存在的模拟时间步长太小的问题，满足VR场景的实时性要求。相较于SPH方法，PBF方法不直接通过气体状态方程计算流体粒子的压强以实现流体的不可压缩性，而是通过密度不变这个约束条件来直接修正粒子的位置从而保证流体的不可压缩性。因此将PBF方法应用于有流体的运动的VR动态场景中，并使用GPU并行化来计算，可以实现实时的、具有真实感的流体模拟。但是，仅仅依靠PBF方法并不能实现流体凝固的模拟，PBF方法主要是将PBD框架应用于SPH方法上，是SPH方法在保证密度不变、计算压强的方法上的改进，因此PBF方法本质上也是对NS方程的一种数值求解近似解的方法；NS方程则是描述粘性不可压缩性流体在运动的过程中的密度保持和动量守恒，求解的是流体运动时由其内力和外力共同作用下所产生的加速度；然而当流体发生凝固时，流体分子间的相互作用发生改变，流体的气相状态发生改变，此时NS方程已经不适用于描述凝固的流体，因此仅仅使用PBF方法是无法进行流体凝固的模拟的。另外，使用PBF方法模拟的流体粒子容易与作为边界的固体粒子发生穿透；PBF方法实现流体的不可压缩性的方法是，基于流体粒子的密度构造一个流体粒子位置相关的运动约束函数，通过是运动约束函数为零来计算流体粒子的位置变化量，使得流体粒子的密度基本等于静态密度，从而实现流体的不可压缩性；流体粒子的位置变化量能否准确计算与其邻域粒子的数目有关，因此流体与固体边界相交处的流体粒子的位置变化量计算可能出现较大误差；并且在外力的作用下可能使得误差变得更大，最终使得流体粒子穿透边界的固体粒子；因此流体凝固后，需要考虑如何防止流动的粒子穿透入凝固立体。
发明内容
本发明所要解决的技术问题在于，提供一种基于PBF的流体凝固模拟方法及系统，可通过不同的粒子速度衰减因子来模拟不同的流体模拟的快慢，还可通过位置变化量修正粒子穿透现象，在GPU的支持下实现更大规模更有效率的流体凝固模拟。
为了解决上述技术问题，本发明提供了一种基于PBF的流体凝固模拟方法，包括：S1，初始化数据；S2，粒子固定半径邻域搜索；S3，更新粒子的凝固状态；S4，将速度衰减因子作用于粒子，计算粒子的位置中间量；S5，计算粒子密度及运动约束量；S6，计算粒子的位置变化量，并进行防穿透修正；S7，更新粒子速度及位置，并处理粒子与边界碰撞问题；S8，计算人造粘性和涡流约束。
作为上述方案的改进，所述步骤S3包括：若粒子的凝固状态为非凝固状态，计算凝固发生值，当凝固发生值达到指定阈值时，更新粒子的凝固状态为在凝固状态；若粒子的凝固状态为在凝固状态，计算速度衰减因子，当速度衰减因子为0时，更新粒子的状态为已凝固状态；若粒子的凝固状态为已凝固状态，保持粒子的状态为已凝固状态。
作为上述方案的改进，所述步骤S4包括：S41，计算粒子所受合力、加速度及当前速度；S42，根据粒子的凝固状态计算出粒子的速度衰减因子，并将速度衰减因子作用于粒子，计算出粒子的位置中间量。
作为上述方案的改进，所述步骤S6包括：S61，根据PBF方法，计算缩放因子、人造压强修正系数及粒子的位置变化量；S62，判断未凝固粒子是否会穿透已凝固粒子，若发生穿透，则修正位置变化量以纠正穿透。
作为上述方案的改进，所述步骤S7包括：S71，将位置中间量和位置变化量之和设置成为粒子的当前位置，根据当前位置及位置变化量计算出粒子的当前速度；S72，判断粒子的当前位置是否与边界发生碰撞，若发生碰撞则相应地更新粒子的当前速度和当前位置。
相应地，本发明还提供了一种基于PBF的流体凝固模拟系统，包括：初始化模块，用于初始化数据；搜索模块，用于搜索粒子固定半径邻域；状态更新模块，用于更新粒子的凝固状态；中间量计算模块，用于将速度衰减因子作用于粒子，计算粒子的位置中间量；密度和运动约束量计算模块，用于计算粒子密度及运动约束量；穿透修正模块，用于计算粒子的位置变化量，并进行防穿透修正；边界碰撞模块，用于更新粒子速度及位置，并处理粒子与边界碰撞问题；人造粘性和涡流计算模块，用于计算人造粘性和涡流约束。
作为上述方案的改进，所述状态更新模块包括：第一状态更新单元，用于若粒子的凝固状态为非凝固状态，计算凝固发生值，当凝固发生值达到指定阈值时，更新粒子的凝固状态为在凝固状态；第二状态更新单元，用于若粒子的凝固状态为在凝固状态，计算速度衰减因子，当速度衰减因子为0时，更新粒子的状态为已凝固状态；第三状态更新单元，用于若粒子的凝固状态为已凝固状态，保持粒子的状态为已凝固状态。
作为上述方案的改进，所述中间量计算模块包括：第一中间量计算单元，用于计算粒子所受合力、加速度及当前速度；第二中间量计算单元，根据粒子的凝固状态计算出粒子的速度衰减因子，并将速度衰减因子作用于粒子，计算出粒子的位置中间量。
作为上述方案的改进，所述穿透修正模块包括：位置变化量计算单元，用于根据PBF方法，计算缩放因子、人造压强修正系数及粒子的位置变化量；修正单元，用于判断未凝固粒子是否会穿透已凝固粒子，若发生穿透，则修正位置变化量以纠正穿透。
作为上述方案的改进，所述边界碰撞模块包括：第一速度和位置更新单元，将位置中间量和位置变化量之和设置成为粒子的当前位置，根据当前位置及位置变化量计算出粒子的当前速度；第二速度和位置更新单元，判断粒子的当前位置是否与边界发生碰撞，若发生碰撞则相应地更新粒子的当前速度和当前位置。
实施本发明，具有如下有益效果：
(1)本发明通过速度衰减的方法来模拟液体的凝固效果。液体的凝固从宏观上看是一个从局部到整体、流动性逐渐变差的液体到固体的变化过程，这个过程包括温度传递和变化以及其他复杂的描述方程。本发明忽略复杂的描述方程，避免了求解复杂的物理方程，通过直接对粒子施加一个速度衰减因子来模拟液体发生凝固时整个宏观的过程，在不明显的失去真实性的前提下，简单高效实现粒子的凝固，极大提高了计算效率；同时，本发明通过改变速度衰减因子的计算函数可以实现不同的凝固快慢效果。
(2)本发明通过粒子位置变化量修正来防止流动的粒子穿透入已经凝固的粒子。穿透的发生是由计算不准确的粒子位置变化量所导致的，本发明将这部分导致穿透发生的位置变化量减去，修正粒子总的位置变化量，从而防止流动粒子穿透凝固粒子，解决了流体已凝固的部分和非凝固部分之间存在的穿透问题。
(3)本发明将整个液体凝固模拟算法使用CUDA实现基于GPU的并行化，并在PBF方法保证较大模拟时间步长的前提下，实现了在VR动态场景中实时模拟流体凝固，能够根据不同的GPU并行化框架实现适用于不同设备的完全并行化版本，在GPU的支持下实现更大规模更有效率的流体凝固模拟。整个模拟算法较大计算量的步骤，包括邻域粒子搜索、计算粒子的各种物理量和边界以及防穿透处理都可以完全并行化，通过CUDA来使GPU一个线程处理一个粒子的模拟计算，并且每一步模拟中每个粒子的物理量的计算都只依赖于上一个计算步骤或者是前一步模拟计算的相关的物理量，数据之间的相关性很低，因此在GPU上液体凝固模拟算法可以高效地并行化计算，从而在VR场景中实时模拟成为可能。
附图说明
图1是欧拉法的原理图；
图2是拉格朗日法的原理图；
图3是本发明基于PBF的流体凝固模拟方法的流程图；
图4是本发明基于PBF的流体凝固模拟系统的示意图；
图5是图4中状态更新模块的示意图；
图6是图4中中间量计算模块的示意图；
图7是图4中穿透修正模块的示意图；
图8是图4中边界碰撞模块的示意图。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明作进一步地详细描述。仅此声明，本发明在文中出现或即将出现的上、下、左、右、前、后、内、外等方位用词，仅以本发明的附图为基准，其并不是对本发明的具体限定。
本发明提出了一种基于PBF的流体凝固模拟方法以实现在VR动态场景中进行具有真实感的、实时的流体凝固模拟。对于可交互式三维场景的应用，像三维游戏、场景漫游甚至虚拟手术等，或着是三维动态场景的还原和重建等，场景中的流体对于真实性、高效性和实时性的要求都是显而易见的。因此，本发明主要为这些需求提供了一种模拟流体凝固的思路。
参见图3，图3显示了本发明基于PBF的流体凝固模拟方法的流程图，包括：
S1，初始化数据；
所述初始化数据包括对模拟过程中保持不变的常量进行初始化计算及对初始数据进行初始化计算，具体地：
计算出模拟过程中保持不变的常量，包括粒子大小、粒子质量、粒子静态密度、粒子平滑半径长度、粒子平滑核函数常数项以及场景和流体的初始规模和位置等，然后通过CUDA(Compute Unified Device Architecture)将这些常量传输到GPU(GraphicsProcessing Unit，图形处理器)的常量内存上，当后面的计算步骤需要使用这些常量是，可直接从GPU的常量内存上取，提高了计算效率。
计算初始数据，如根据粒子的大小和流体的初始规模和位置计算出每个粒子的位置，还有粒子的初始速度和颜色等等，根据粒子的总个数通过CUDA在GPU上申请所需大小的全局内存，并将在CPU内存上计算得的初始数据复制到相应的GPU全局内存上，以便于接下来通过CUDA执行所有计算步骤的并行化运算。
接下来的步骤S2到S8，都可以通过设置合适的CUDA线程块和每块线程数，以使得每个粒子的每步计算都可以分配到GPU的一个线程上单独执行，从而实现整个模拟算法的并行化计算。
S2，粒子固定半径邻域搜索；
根据流体运动空间的包围盒大小以及平滑半径长度将包围盒划分成大小均匀的网长方体网格。根据粒子的位置，计算出粒子所在网格格子的索引，并为每个格子维护一个其所包含的粒子的索引表。从而对于每个粒子，查找器固定半径内的邻域粒子只需要至多在27个网格格子内寻找，即每个粒子本身所在格子以及与该格子紧密邻接的另外至多26个格子。在这27个格子内，对于每个粒子而言，与其距离小于等于平滑核半径的粒子即为其邻域粒子。为每个粒子维护其邻域粒子索引列表以便于步骤S3至S8中需要邻域粒子信息的步骤能够快速查找。
S3，更新粒子的凝固状态；
粒子的凝固状态分为三种，分别为：非凝固状态NONCLOT、在凝固状态CLOTTING和已凝固状态CLOTTED。本发明根据粒子的凝固状态进行相应的计算。
(1)若粒子的凝固状态为非凝固状态，计算凝固发生值Thi，当凝固发生值达到指定阈值时，更新粒子的凝固状态为在凝固状态；
若粒子的凝固状态为“非凝固状态NONCLOT”，计算其凝固发生值Thi，计算公式如下：
其中，Wij＝W(xij,h)，xij和h分别表示粒子i到粒子j的位移和固定的核平滑半径，m、ρ0、Thj分别表示粒子的质量、初始密度和凝固发生值；当Thi的值达到指定阈值时，粒子开始凝固，此时更新粒子的凝固状态为“在凝固状态CLOTTING”；
(2)若粒子的凝固状态为在凝固状态，计算速度衰减因子，当速度衰减因子为0时，更新粒子的状态为已凝固状态；
若粒子的凝固状态为“在凝固状态CLOTTING”，计算其速度衰减因子di＝max(Decay(x),0)，其中函数Decay(x)为速度衰减函数，x随着模拟的推进呈线性增长，取为x＝x*+Δx，其中x*和Δx分别表示粒子从上一步模拟计算得到的位置和粒子在该步模拟计算中的位移量，函数Decay(x)给出如下：
其中，thresholdx应为接近1但是不大于1的数，目的是当粒子的速度足够缓慢时即可视为凝固，从而避免了一些不必要的计算，提高计算效率。本发明中，取x＝0.95，当di＝0时，更新粒子的状态为“已凝固状态CLOTTED”
(3)若粒子的凝固状态为已凝固状态，保持粒子的状态为已凝固状态。
若粒子的凝固状态为“已凝固状态CLOTTED”，粒子的速度为0，粒子不再发生位移，可以避免其他流动或未凝固粒子所需要的计算。
S4，将速度衰减因子作用于粒子，计算粒子的位置中间量；
具体地，所述步骤S4包括：
S41，计算粒子所受合力、加速度及当前速度；
粒子所受合力F主要包括重力G和涡流约束产生的力Fvor，则粒子的当前速度ai计算如下：
根据向前欧拉法，计算出粒子的当前速度其中vi表示粒子从上一步模拟计算得到的速度，Δt表示模拟的时间步长
S42，根据粒子的凝固状态计算出粒子的速度衰减因子，并将速度衰减因子作用于粒子，计算出粒子的位置中间量。
判断粒子的凝固状态，若凝固状态为“非凝固状态NONCLOT”，则粒子的速度衰减因子di＝1；若凝固状态为“在凝固状态CLOTTING”，则di由步骤S3计算出，且值域为(0，1)；若凝固状态为“已凝固状态CLOTTED”，则di＝0。将速度衰减因子作用于粒子，则有当前粒子的速度
进而使用向前欧拉法计算出粒子的位置中间量
其中xi表示粒子从上一步模拟计算得到的位置。
S5，计算粒子密度及运动约束量；
粒子的密度ρi的计算公式给出如下，
其中，m为粒子质量，Wij为Poly6平滑核函数，
本发明中需要用到核函数进行平滑插值的计算本均使用Poly6平滑核函数。根据PBF方法，运动约束量Ci计算公式如下：
S6，计算粒子的位置变化量，并进行防穿透修正；
具体地，所述步骤S6包括：
S61，根据PBF方法，计算缩放因子、人造压强修正系数及粒子的位置变化量；
由于流体的不可压缩性的要求，粒子的密度应该保持不变，即ρi＝ρ0，因此由泰勒展开公式可得，其中表示Ci的梯度，表示的转置向量，λi表示缩放因子，ε为不小于0的常数，用于防止梯度消失造成的数值不稳定的现象发生。因此缩放因子的计算公式如下：
其中为Spiky核平滑函数的梯度函数，
表示粒子k的Spiky核平滑函数的梯度函数，k为粒子的索引，当k＝i表示粒子本身，否则当k＝j表示粒子的固定半径内的邻域粒子。
本发明中需要用到核函数的梯度进行平滑插值的计算本均使用Spiky平滑函数的梯度函数。人造压强用于防止粒子出现不自然的聚集现象，人造压强修正系数Scorr的计算公式如下：
其中|Δq|＝0.3h，根据缩放因子和人造压强修正系数，粒子的位置变化量Δxi按如下公式计算
S62，判断未凝固粒子是否会穿透已凝固粒子，若发生穿透，则修正位置变化量以纠正穿透。
假设将位置变化量加入到当前位置中，未凝固粒子i与已凝固粒子j之间的距离为dist＝|xij|，粒子的半径为r，若dist<1.8r，则认为穿透发生。这是由于过多的位置变化量所导致的，因此需要将这部分变化量消除。计算出从位置xj指向位置xi的向量xji的法向量修正位置变化量Δxi
Δxi＝Δxi+(1.8r-dist)*nji (18)
通过修正位置变化量使得该防粒子穿透方法可以直接应用于PBF中。
S7，更新粒子速度及位置，并处理粒子与边界碰撞问题；
具体地，所述步骤S7包括：
S71，将位置中间量和位置变化量之和设置成为粒子的当前位置，根据当前位置及位置变化量计算出粒子的当前速度；
根据向前欧拉法，粒子的当前速度更新粒子的当前位置
S72，判断粒子的当前位置是否与边界发生碰撞，若发生碰撞则相应地更新粒子的当前速度和当前位置。
为了方便流体与边界的碰撞检测处理，本发明将流体所能到达的运动空间设置为轴对齐长方体包围盒。长方体包围盒由最小点N(xn,yn,zn)和最大点M(xm,ym,zm)，则包围盒内的流体粒子的当前位置xi(xi,yi,zi)更新为
xi＝min(max(xi,xn),xm) (19)
yi和zi同理。当发生碰撞时，粒子的入射速度方向关于碰撞边界面的法向np发生反射，因此发生碰撞时，粒子的速度更新公式为
vi＝vi-2*vi*(vi·np) (20)
S8，计算人造粘性和涡流约束。
本发明直接应用XSPH方法的来计算人造粘性对粒子速度产生的影响，该方法的计算公式为
涡流约束用于防止模拟的流体出现不自然的能量损失现象。粒子的涡流ωi计算公式如下：
其中vij表示粒子i相对于其某一固定半径内的邻域粒子j的速度差。从而可以计算出涡流力
Fi vor＝ε(N×ωi) (23)
其中
最后，将通过CUDA让OpenGL(Open Graphics Library)或者是DirectX 3D直接从GPU中读取渲染流体粒子所需要的数据，包括位置、速度或颜色等，其中位置是必须的，从而避免了在CPU和GPU之间来回传输大量数据所占用的大量时间。
参见图4，图4显示了本发明基于PBF的流体凝固模拟系统100的具体结构，其包括初始化模块1、搜索模块2、状态更新模块3、中间量计算模块4、密度和运动约束量计算模块5、穿透修正模块6、边界碰撞模块7及人造粘性和涡流计算模块8，具体地：
初始化模块：
初始化模块1，用于初始化数据；所述初始化数据包括对模拟过程中保持不变的常量进行初始化计算及对初始数据进行初始化计算，具体地：计算出模拟过程中保持不变的常量，包括粒子大小、粒子质量、粒子静态密度、粒子平滑半径长度、粒子平滑核函数常数项以及场景和流体的初始规模和位置等，然后通过CUDA(Compute Unified DeviceArchitecture)将这些常量传输到GPU(Graphics Processing Unit，图形处理器)的常量内存上，当后面的计算步骤需要使用这些常量是，可直接从GPU的常量内存上取，提高了计算效率。计算初始数据，如根据粒子的大小和流体的初始规模和位置计算出每个粒子的位置，还有粒子的初始速度和颜色等等，根据粒子的总个数通过CUDA在GPU上申请所需大小的全局内存，并将在CPU内存上计算得的初始数据复制到相应的GPU全局内存上，以便于接下来通过CUDA执行所有计算步骤的并行化运算。
搜索模块：
搜索模块2，用于搜索粒子固定半径邻域；搜索模块2根据流体运动空间的包围盒大小以及平滑半径长度将包围盒划分成大小均匀的网长方体网格。根据粒子的位置，计算出粒子所在网格格子的索引，并为每个格子维护一个其所包含的粒子的索引表。从而对于每个粒子，查找器固定半径内的邻域粒子只需要至多在27个网格格子内寻找，即每个粒子本身所在格子以及与该格子紧密邻接的另外至多26个格子。在这27个格子内，对于每个粒子而言，与其距离小于等于平滑核半径的粒子即为其邻域粒子。为每个粒子维护其邻域粒子索引列表以便于步骤S3至S8中需要邻域粒子信息的步骤能够快速查找。
状态更新模块：
状态更新模块3，用于更新粒子的凝固状态；粒子的凝固状态分为三种，分别为：非凝固状态NONCLOT、在凝固状态CLOTTING和已凝固状态CLOTTED。状态更新模块3根据粒子的凝固状态进行相应的计算。
如图5所示，所述状态更新模块3包括第一状态更新单元31、第二状态更新单元32及第三状态更新单元33，其中：
第一状态更新单元31，用于若粒子的凝固状态为非凝固状态，计算凝固发生值，当凝固发生值达到指定阈值时，更新粒子的凝固状态为在凝固状态；
若粒子的凝固状态为“非凝固状态NONCLOT”，计算其凝固发生值Thi，计算公式如下：
其中，Wij＝W(xij,h)，xij和h分别表示粒子i到粒子j的位移和固定的核平滑半径，m、ρ0、Thj分别表示粒子的质量、初始密度和凝固发生值；当Thi的值达到指定阈值时，粒子开始凝固，此时更新粒子的凝固状态为“在凝固状态CLOTTING”；
第二状态更新单元32，用于若粒子的凝固状态为在凝固状态，计算速度衰减因子，当速度衰减因子为0时，更新粒子的状态为已凝固状态；
若粒子的凝固状态为“在凝固状态CLOTTING”，计算其速度衰减因子di＝max(Decay(x),0)，其中函数Decay(x)为速度衰减函数，x随着模拟的推进呈线性增长，取为x＝x*+Δx，其中x*和Δx分别表示粒子从上一步模拟计算得到的位置和粒子在该步模拟计算中的位移量，函数Decay(x)给出如下：
其中，thresholdx应为接近1但是不大于1的数，目的是当粒子的速度足够缓慢时即可视为凝固，从而避免了一些不必要的计算，提高计算效率。本发明中，取x＝0.95，当di＝0时，更新粒子的状态为“已凝固状态CLOTTED”
第三状态更新单元33，用于若粒子的凝固状态为已凝固状态，保持粒子的状态为已凝固状态。
若粒子的凝固状态为“已凝固状态CLOTTED”，粒子的速度为0，粒子不再发生位移，可以避免其他流动或未凝固粒子所需要的计算。
中间量计算模块：
中间量计算模块4，用于将速度衰减因子作用于粒子，计算粒子的位置中间量；
如图6所示，所述中间量计算模块4包括第一中间量计算单元41及第二中间量计算单元42：
第一中间量计算单元41，用于计算粒子所受合力、加速度及当前速度；
粒子所受合力F主要包括重力G和涡流约束产生的力Fvor，则粒子的当前速度ai计算如下：
根据向前欧拉法，计算出粒子的当前速度其中vi表示粒子从上一步模拟计算得到的速度，Δt表示模拟的时间步长
第二中间量计算单元42，根据粒子的凝固状态计算出粒子的速度衰减因子，并将速度衰减因子作用于粒子，计算出粒子的位置中间量。
判断粒子的凝固状态，若凝固状态为“非凝固状态NONCLOT”，则粒子的速度衰减因子di＝1；若凝固状态为“在凝固状态CLOTTING”，则di由步骤S3计算出，且值域为(0，1)；若凝固状态为“已凝固状态CLOTTED”，则di＝0。将速度衰减因子作用于粒子，则有当前粒子的速度
进而使用向前欧拉法计算出粒子的位置中间量
其中xi表示粒子从上一步模拟计算得到的位置。
密度和运动约束量计算模块5，用于计算粒子密度及运动约束量；
粒子的密度ρi的计算公式给出如下，
其中，m为粒子质量，Wij为Poly6平滑核函数，
本发明中需要用到核函数进行平滑插值的计算本均使用Poly6平滑核函数。根据PBF方法，运动约束量Ci计算公式如下：
穿透修正模块：
穿透修正模块6，用于计算粒子的位置变化量，并进行防穿透修正；
如图7所示，所述穿透修正模块6包括位置变化量计算单元61及修正单元62：
位置变化量计算单元61，用于根据PBF方法，计算缩放因子、人造压强修正系数及粒子的位置变化量；
由于流体的不可压缩性的要求，粒子的密度应该保持不变，即ρi＝ρ0，因此由泰勒展开公式可得，其中表示Ci的梯度，表示的转置向量，λi表示缩放因子，ε为不小于0的常数，用于防止梯度消失造成的数值不稳定的现象发生。因此缩放因子的计算公式如下：
其中为Spiky核平滑函数的梯度函数
表示粒子k的Spiky核平滑函数的梯度函数，k为粒子的索引，当k＝i表示粒子本身，否则当k＝j表示粒子的固定半径内的邻域粒子。
本发明中需要用到核函数的梯度进行平滑插值的计算本均使用Spiky平滑函数的梯度函数。人造压强用于防止粒子出现不自然的聚集现象，人造压强修正系数Scorr的计算公式如下：
其中|Δq|＝0.3h，根据缩放因子和人造压强修正系数，粒子的位置变化量Δxi按如下公式计算
修正单元62，用于判断未凝固粒子是否会穿透已凝固粒子，若发生穿透，则修正位置变化量以纠正穿透。
假设将位置变化量加入到当前位置中，未凝固粒子i与已凝固粒子j之间的距离为dist＝|xij|，粒子的半径为r，若dist<1.8r，则认为穿透发生。这是由于过多的位置变化量所导致的，因此需要将这部分变化量消除。计算出从位置xj指向位置xi的向量xji的法向量修正位置变化量Δxi
Δxi＝Δxi+(1.8r-dist)*nji (18)
通过修正位置变化量使得该防粒子穿透方法可以直接应用于PBF中。
边界碰撞模块7，用于更新粒子速度及位置，并处理粒子与边界碰撞问题；
边界碰撞模块：
如图8所示，所述边界碰撞模块7包括：
第一速度和位置更新单元71，将位置中间量和位置变化量之和设置成为粒子的当前位置，根据当前位置及位置变化量计算出粒子的当前速度；
根据向前欧拉法，粒子的当前速度更新粒子的当前位置
第二速度和位置更新单元72，判断粒子的当前位置是否与边界发生碰撞，若发生碰撞则相应地更新粒子的当前速度和当前位置。
为了方便流体与边界的碰撞检测处理，本发明将流体所能到达的运动空间设置为轴对齐长方体包围盒。长方体包围盒由最小点N(xn,yn,zn)和最大点M(xm,ym,zm)，则包围盒内的流体粒子的当前位置xi(xi,yi,zi)更新为
xi＝min(max(xi,xn),xm) (19)
yi和zi同理。当发生碰撞时，粒子的入射速度方向关于碰撞边界面的法向np发生反射，因此发生碰撞时，粒子的速度更新公式为
vi＝vi-2*vi*(vi·np) (20)
人造粘性和涡流计算模块：
人造粘性和涡流计算模块8，用于计算人造粘性和涡流约束。
本发明直接应用XSPH方法的来计算人造粘性对粒子速度产生的影响，该方法的计算公式为
涡流约束用于防止模拟的流体出现不自然的能量损失现象。粒子的涡流ωi计算公式如下：
其中vij表示粒子i相对于其某一固定半径内的邻域粒子j的速度差。从而可以计算出涡流力
Fi vor＝ε(N×ωi) (23)
其中
最后，将通过CUDA让OpenGL(Open Graphics Library)或者是DirectX 3D直接从GPU中读取渲染流体粒子所需要的数据，包括位置、速度或颜色等，其中位置是必须的，从而避免了在CPU和GPU之间来回传输大量数据所占用的大量时间。
由上可知，本发明具有以下有益效果：
(1)本发明通过速度衰减的方法来模拟液体的凝固效果。液体的凝固从宏观上看是一个从局部到整体、流动性逐渐变差的液体到固体的变化过程，这个过程包括温度传递和变化以及其他复杂的描述方程。本发明忽略复杂的描述方程，避免了求解复杂的物理方程，通过直接对粒子施加一个速度衰减因子来模拟液体发生凝固时整个宏观的过程，在不明显的失去真实性的前提下，简单高效实现粒子的凝固，极大提高了计算效率；同时，本发明通过改变速度衰减因子的计算函数可以实现不同的凝固快慢效果，并在PBF方法保证较大模拟时间步长的前提下，实现了在VR动态场景中实时模拟流体凝固。
(2)本发明通过粒子位置变化量修正来防止流动的粒子穿透入已经凝固的粒子。穿透的发生是由计算不准确的粒子位置变化量所导致的，本发明将这部分导致穿透发生的位置变化量减去，修正粒子总的位置变化量，从而防止流动粒子穿透凝固粒子，解决了流体已凝固的部分和非凝固部分之间存在的穿透问题。
(3)本发明将整个液体凝固模拟算法使用CUDA实现基于GPU的并行化，能够根据不同的GPU并行化框架实现适用于不同设备的完全并行化版本，在GPU的支持下实现更大规模更有效率的流体凝固模拟。整个模拟算法较大计算量的步骤，包括邻域粒子搜索、计算粒子的各种物理量和边界以及防穿透处理都可以完全并行化，通过CUDA来使GPU一个线程处理一个粒子的模拟计算，并且每一步模拟中每个粒子的物理量的计算都只依赖于上一个计算步骤或者是前一步模拟计算的相关的物理量，数据之间的相关性很低，因此在GPU上液体凝固模拟算法可以高效地并行化计算，从而在VR场景中实时模拟成为可能。
以上所述是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也视为本发明的保护范围。
Claims (6)
1.一种基于PBF的流体凝固模拟方法，其特征在于，包括：
S1，初始化数据；
S2，粒子固定半径邻域搜索；具体地，根据流体运动空间的包围盒大小以及平滑半径长度将包围盒划分成大小均匀的网长方体网格，根据粒子的位置，计算出粒子所在网格格子的索引，并为每个格子维护一个其所包含的粒子的索引表，从而对于每个粒子，查找器固定半径内的邻域粒子只需要至多在27个网格格子内寻找，即每个粒子本身所在格子以及与该格子紧密邻接的另外至多26个格子，在这27个格子内，对于每个粒子而言，与其距离小于等于平滑核半径的粒子即为其邻域粒子；
S3，更新粒子的凝固状态；粒子的凝固状态分为三种，分别为：非凝固状态、在凝固状态和已凝固状态；若粒子的凝固状态为“非凝固状态”，计算凝固发生值Thi，当Thi的值达到指定阈值时，更新粒子的凝固状态为“在凝固状态”，其中，Wij＝W(xij,h)，xij和h分别表示粒子i到粒子j的位移和固定的核平滑半径，m、ρ0、Thj分别表示粒子的质量、初始密度和凝固发生值；若粒子的凝固状态为在凝固状态，计算速度衰减因子，当速度衰减因子为0时，更新粒子的状态为已凝固状态；
S4，将速度衰减因子作用于粒子，计算粒子的位置中间量；具体地，所述步骤S4包括：S41，计算粒子所受合力、加速度及当前速度；S42，根据粒子的凝固状态计算出粒子的速度衰减因子，并将速度衰减因子作用于粒子，计算出粒子的位置中间量；
S5，计算粒子密度及运动约束量；
S6，计算粒子的位置变化量，并进行防穿透修正；
S7，更新粒子速度及位置，并处理粒子与边界碰撞问题；
S8，计算人造粘性和涡流约束；
所述步骤S2到S8，通过设置合适的CUDA线程块和每块线程数，以使得每个粒子的每步计算都可以分配到GPU的一个线程上单独执行，从而实现整个模拟算法的并行化计算。
2.如权利要求1所述的基于PBF的流体凝固模拟方法，其特征在于，所述步骤S6包括：
S61，根据PBF方法，计算缩放因子、人造压强修正系数及粒子的位置变化量；
S62，判断未凝固粒子是否会穿透已凝固粒子，若发生穿透，则修正位置变化量以纠正穿透。
3.如权利要求1所述的基于PBF的流体凝固模拟方法，其特征在于，所述步骤S7包括：
S71，将位置中间量和位置变化量之和设置成为粒子的当前位置，根据当前位置及位置变化量计算出粒子的当前速度；
S72，判断粒子的当前位置是否与边界发生碰撞，若发生碰撞则相应地更新粒子的当前速度和当前位置。
4.一种基于PBF的流体凝固模拟系统，其特征在于，包括：
初始化模块，用于初始化数据；
搜索模块，用于搜索粒子固定半径邻域；具体地，根据流体运动空间的包围盒大小以及平滑半径长度将包围盒划分成大小均匀的网长方体网格，根据粒子的位置，计算出粒子所在网格格子的索引，并为每个格子维护一个其所包含的粒子的索引表，从而对于每个粒子，查找器固定半径内的邻域粒子只需要至多在27个网格格子内寻找，即每个粒子本身所在格子以及与该格子紧密邻接的另外至多26个格子，在这27个格子内，对于每个粒子而言，与其距离小于等于平滑核半径的粒子即为其邻域粒子；
状态更新模块，用于更新粒子的凝固状态；粒子的凝固状态分为三种，分别为：非凝固状态、在凝固状态和已凝固状态；若粒子的凝固状态为“非凝固状态”，计算凝固发生值Thi，当Thi的值达到指定阈值时，更新粒子的凝固状态为“在凝固状态”，其中，Wij＝W(xij,h)，xij和h分别表示粒子i到粒子j的位移和固定的核平滑半径，m、ρ0、Thj分别表示粒子的质量、初始密度和凝固发生值；若粒子的凝固状态为在凝固状态，计算速度衰减因子，当速度衰减因子为0时，更新粒子的状态为已凝固状态；
中间量计算模块，用于将速度衰减因子作用于粒子，计算粒子的位置中间量；具体地，所述中间量计算模块包括：第一中间量计算单元，用于计算粒子所受合力、加速度及当前速度；第二中间量计算单元，根据粒子的凝固状态计算出粒子的速度衰减因子，并将速度衰减因子作用于粒子，计算出粒子的位置中间量；
密度和运动约束量计算模块，用于计算粒子密度及运动约束量；
穿透修正模块，用于计算粒子的位置变化量，并进行防穿透修正；
边界碰撞模块，用于更新粒子速度及位置，并处理粒子与边界碰撞问题；
人造粘性和涡流计算模块，用于计算人造粘性和涡流约束；
系统通过设置合适的CUDA线程块和每块线程数，以使得每个粒子的每步计算都可以分配到GPU的一个线程上单独执行，从而实现整个模拟算法的并行化计算。
5.如权利要求4所述的基于PBF的流体凝固模拟系统，其特征在于，所述穿透修正模块包括：
位置变化量计算单元，用于根据PBF方法，计算缩放因子、人造压强修正系数及粒子的位置变化量；
修正单元，用于判断未凝固粒子是否会穿透已凝固粒子，若发生穿透，则修正位置变化量以纠正穿透。
6.如权利要求4所述的基于PBF的流体凝固模拟系统，其特征在于，所述边界碰撞模块包括：
第一速度和位置更新单元，将位置中间量和位置变化量之和设置成为粒子的当前位置，根据当前位置及位置变化量计算出粒子的当前速度；
第二速度和位置更新单元，判断粒子的当前位置是否与边界发生碰撞，若发生碰撞则相应地更新粒子的当前速度和当前位置。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201811043656.6A CN109344450B (zh) | 2018-09-07 | 2018-09-07 | 基于pbf的流体凝固模拟方法及系统 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201811043656.6A CN109344450B (zh) | 2018-09-07 | 2018-09-07 | 基于pbf的流体凝固模拟方法及系统 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN109344450A CN109344450A (zh) | 2019-02-15 |
| CN109344450B true CN109344450B (zh) | 2019-07-23 |
Family
ID=65304964
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201811043656.6A Active CN109344450B (zh) | 2018-09-07 | 2018-09-07 | 基于pbf的流体凝固模拟方法及系统 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109344450B (zh) |
Families Citing this family (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111047707B (zh) * | 2019-11-11 | 2021-09-28 | 南昌大学 | 一种基于sph的混合粒子血液模型的流血仿真方法 |
| CN111581875A (zh) * | 2020-04-16 | 2020-08-25 | 天津大学 | 一种自适应粒子流体的计算机模拟方法 |
| CN112561775B (zh) * | 2021-02-22 | 2021-05-04 | 欢乐互娱(上海)科技股份有限公司 | 基于gpu实现粒子阻尼效果的方法、装置及存储介质 |
| CN113033068B (zh) * | 2021-04-23 | 2022-07-22 | 上海交通大学 | 容器内流体粒子沸腾时的视觉仿真方法及电子设备 |
| CN113824990A (zh) * | 2021-08-18 | 2021-12-21 | 北京达佳互联信息技术有限公司 | 视频生成方法、装置及存储介质 |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103729555A (zh) * | 2013-12-20 | 2014-04-16 | 深圳先进技术研究院 | 一种模拟血流与血管壁作用的方法和装置 |
| CN104200015A (zh) * | 2014-08-20 | 2014-12-10 | 清华大学 | 一种流体模拟方法及装置 |
| CN104268943A (zh) * | 2014-09-28 | 2015-01-07 | 北京航空航天大学 | 一种基于欧拉-拉格朗日耦合方法的流体仿真方法 |
Family Cites Families (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20150161810A1 (en) * | 2013-12-10 | 2015-06-11 | Nvidia Corporation | Position based fluid dynamics simulation |
| CN104360896B (zh) * | 2014-12-04 | 2017-12-15 | 北京航空航天大学 | 一种基于gpu集群的并行流体仿真加速方法 |
| CN104991999A (zh) * | 2015-06-17 | 2015-10-21 | 大连理工大学 | 一种基于二维sph的溃坝洪水演进模拟方法 |
-
2018
- 2018-09-07 CN CN201811043656.6A patent/CN109344450B/zh active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103729555A (zh) * | 2013-12-20 | 2014-04-16 | 深圳先进技术研究院 | 一种模拟血流与血管壁作用的方法和装置 |
| CN104200015A (zh) * | 2014-08-20 | 2014-12-10 | 清华大学 | 一种流体模拟方法及装置 |
| CN104268943A (zh) * | 2014-09-28 | 2015-01-07 | 北京航空航天大学 | 一种基于欧拉-拉格朗日耦合方法的流体仿真方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN109344450A (zh) | 2019-02-15 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109344450B (zh) | 基于pbf的流体凝固模拟方法及系统 | |
| Rasmussen et al. | Directable photorealistic liquids | |
| Foster et al. | Practical animation of liquids | |
| Schechter et al. | Evolving sub-grid turbulence for smoke animation | |
| US7647214B2 (en) | Method for simulating stable but non-dissipative water | |
| Diziol et al. | Robust real-time deformation of incompressible surface meshes | |
| Ando et al. | A particle-based method for preserving fluid sheets | |
| CN110717269B (zh) | 一种基于网格和粒子耦合的流体表面细节保护方法 | |
| Band et al. | Moving Least Squares Boundaries for SPH Fluids. | |
| Li et al. | Fast simulation of deformable characters with articulated skeletons in projective dynamics | |
| Huang et al. | Ships, splashes, and waves on a vast ocean | |
| Shah et al. | Extended galilean invariance for adaptive fluid simulation | |
| Keiser et al. | Multiresolution particle-based fluids | |
| Liu et al. | Turbulent details simulation for SPH fluids via vorticity refinement | |
| Auer et al. | A semi‐Lagrangian closest point method for deforming surfaces | |
| Feng et al. | Detail‐preserving SPH fluid control with deformation constraints | |
| Gao et al. | Accelerating liquid simulation with an improved data‐driven method | |
| Angst et al. | Robust and efficient wave simulations on deforming meshes | |
| Wu et al. | Improved divergence‐free smoothed particle hydrodynamics via priority of divergence‐free solver and SOR | |
| Jang et al. | A geometric approach to animating thin surface features in smoothed particle hydrodynamics water | |
| Akinci | Interface handling in smoothed particle hydrodynamics | |
| Shi et al. | Controllable motion synthesis in a gaseous medium | |
| Abdelnaim et al. | Fluid-structure interactions simulation and visualization using ISPH approach | |
| Shao | Fluids, threads and fibers: towards high performance physics-based modeling and simulation | |
| Yang et al. | A semi-explicit surface tracking mechanism for multi-phase immiscible liquids |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant | ||
| OL01 | Intention to license declared | ||
| OL01 | Intention to license declared |