基于PCISPH方法的海底滑坡仿真研究

吴蔚

摘  要：海底滑坡是一种复杂的地质现象，对海洋工程设施和人类生命财产安全有着巨大的潜在威胁。为了对海底滑坡研究提供可视化支持，提出基于PCISPH非可压缩流体模型的海底滑坡模拟方法，使用一种改进的算法进行边界处理，并通过移动立方体算法提取流体表面。实验表明该方法可以高效地实现海底滑坡仿真过程，具有较强的应用价值。

关键词：海底滑坡;PCISPH;表面重建

中图分类号：TU470         文献标志码：A         文章编号：2095-2945（2020）14-0018-02

Abstract： Submarine landslide is a complex geological phenomenon， which poses a great potential threat to marine engineering facilities and to the safety of life and property of the people. In order to provide visual support for submarine landslide research， a PCISPH-based incompressible fluid simulation model to simulate submarine landlide is proposed. A boundary handling method is improved， and the Marching Cubes Algorithmis used to extract the fluid surface. The experiment suggests that our method can realize the submarine landlide simulation effectively and has strong practical value.

Keywords： submarine landslide; PCISPH; surface reconstruction

1 概述

海底滑坡是一种复杂的深海地质现象，目前对于滑坡失稳后的滑块滑动过程很难用一个统一的物理模型进行准确描述。然而现阶段针对滑块滑动的物理机制和滑动特性的研究已经比较深入，比如大规模的滑块滑动和较长的滑移距离是海底滑坡相比陆上滑坡最典型的特征。

为了更深入的理解海底滑坡活动，研究者采用很多方法对滑块滑动过程进行了模拟，特别是基于粒子的SPH（Smoothed Particle Hydrodynamics）方法已經被广泛应用到海洋工程和地质灾害等复杂的流体问题的研究中[1]。然而当面对大规模的流体场景时，标准SPH方法模拟的流体会被极大地压缩，造成严重的失真现象。为了解决这个问题，很多不可压缩流体模拟方法被提出，典型的如PCISPH（Predictive-Corrective Incompressible SPH）[2]方法，不仅可以很好地增强流体的不可压缩性，而且通过减少时间步长可极大提升模拟效率。

本文提出了基于PCISPH的海底滑坡仿真方法，使用一种改进的算法进行边界处理，并采用移动立方体算法[3]对流体粒子进行表面重建。整个流体模拟算法完全基于CUDA架构进行并行化处理，场景的渲染通过OpenGL实现。

2 PCISPH

SPH是一种模拟流体运动的无网格拉格朗日方法，相比于有限元分析方法可以很好地解决流体大变形问题。该方法将连续的流体离散成有限数目的粒子，每个粒子携带多个物理量。任意位置的物理量可通过搜索邻域内一定范围的粒子进行插值得到。根据对粒子压强计算方法的不同，SPH方法模拟的流体可被分为可压缩流体和不可压缩流体。标准SPH模拟的流体具有可压缩性，通常使用增强的气体状态方程计算压强力从而实现不可压缩（或弱可压缩）效果[4]。这种方法实现简单，但相应的代价是需要设置很小的时间步长才能保证流体模拟的稳定性，这在很大程度上降低了算法整体运行效率，因此通常用于小规模的流体模拟中。

为了更好适应海底滑坡场景中大规模滑块滑动的特点，实现高效的不可压缩流体模拟，本文提出使用基于预测-矫正策略的PCISPH方法进行海底滑坡仿真。该方法允许设置较大的时间步长，每个时间步迭代地预测下一帧的速度和位置从而更新粒子密度值。若密度大于设定的参数，说明粒子之间过于紧凑，这就需要增加粒子的压强和压强力从而使粒子间变得稀疏。经过若干次迭代就可满足设定的流体压缩率的要求。这种迭代式的计算方法使得每一时间步的计算量比标准SPH更多，但由于允许设置更大的时间步长，因此总体上仍具有更高的流体模拟效率。

3 边界处理

滑块自身的流体式滑动基于PCISPH方法建模，滑坡表面采用三角形网格进行建模。当滑块沿着滑坡表面滑动时，需要正确处理流体粒子和边界面之间的交互过程。该过程可分为两个步骤：碰撞检测和碰撞响应。所谓碰撞检测即预测粒子前进的过程中可能与哪些滑坡表面三角形发生碰撞。这需要得到粒子当前位置和预测位置构成的线段，然后与三角面进行相交测试。碰撞响应则指对于可能发生碰撞的粒子通过调整速度、位置等属性从而避免潜在碰撞的发生。

传统的连续碰撞检测方法在检测到碰撞时，下一帧的预测位置被矫正到三角形表面。如图1所示，Pi为当前粒子位置，Pj为预测的位置，检测到线段PiPj与三角形发生碰撞。通过减去导致粒子位置由Pa到Pj的速度分量，下一帧的粒子实际位置被矫正到表面Pa。然而这种方法使得粒子在凹陷区域可能会被卡住从而产生抖动现象[5]，耗费大量的计算资源。为了解决这个问题，本文改进碰撞响应的方式，通过只减去造成粒子垂直方向穿透的速度分量，粒子被矫正到了表面位置Pb。这种改进可使被卡住的粒子有效脱离凹陷区域，从而矫正模拟中出现的错误。

4 流體表面重建

当每一帧流体粒子的位置被计算出以后，还需要基于这些位置信息重建出光滑的流体表面，从而提升海底滑坡仿真的效果。可通过两步实现：标量域计算和等值面提取。首先整个模拟区域被划分为网格，针对每个网格顶点计算出一个标量值。然后基于得到的离散标量域，采用移动立方体算法在每个网格单元内部提取出多边形构成的等值面，这样就重建出流体的表面区域。这里针对每个网格顶点和网格单元的计算是独立的，因此可以采用完全并行的方式实现。

表面重建的运算效率主要受标量域计算的影响。特别是在高分辨率网格中，每个网格顶点的标量值计算都需要对邻域内所有粒子进行搜索以完成插值。为了提高算法的运行效率，很多研究者使用局部域技术，即只针对流体表面附近区域的网格顶点的标量值进行计算。本文使用Akinci等人[6]的局部域算法，通过高效的并行化技术，提升了流体表面重建的运算和存储效率。

5 结论

为了对本文算法进行验证，设计了海底滑坡仿真实验，使用63K流体粒子模拟滑块滑动过程，如图2所示。初始滑块设置为长方体，在重力作用及与滑坡表面的边界交互作用下，滑块沿着复杂的沟壑滑动并逐渐分离，较完整地实现了滑坡失稳后滑块滑动的过程。

本文基于PCISPH方法高效地实现了海底滑坡仿真过程，对于深入理解海底滑坡地质现象提供了可视化支持，具有较强的应用价值。今后的研究将在本文算法的基础上，不仅实现滑块与滑坡体之间的交互，也要考虑更复杂的影响因素，比如流体颗粒之间的摩擦力、海水介质对颗粒流的影响等，以期进一步提升海底滑坡的仿真效果。

参考文献：

[1]Cleary， P. W.， Prakash， M.， Sinnott， M. D.， Rudman， M.， & Das， R. （2011）. Large scale simulation of industrial， engineering and geophysical flows using particle methods. In Particle-Based Methods （pp. 89-111）. Springer Netherlands.

[2]Solenthaler， B.， & Pajarola， R. （2009， July）. Predictive-corrective incompressible SPH. In ACM transactions on graphics （TOG） （Vol. 28， No. 3， p. 40）. ACM.

[3]Lorensen， W. E.， & Cline， H. E. （1987， August）. Marching cubes： A high resolution 3D surface construction algorithm. In ACM siggraph computer graphics（Vol. 21， No. 4， pp. 163-169）. ACM.

[4]Becker， M.， & Teschner， M. （2007， August）. Weakly compressible SPH for free surface flows. In Proceedings of the 2007 ACM SIGGRAPH/Eurographics symposium on Computer animation （pp. 209-217）. Eurographics Association.

[5]Grahn， A. （2008）. Interactive simulation of contrast fluid using smoothed particle hydrodynamics （Doctoral dissertation， Master's thesis， Ume■ University）.

[6]Akinci G， Ihmsen M， Akinci N， et al. Parallel Surface Reconstruction for Particle-Based Fluids[J]. Computer Graphics Forum， 2012，31（6）：1797-1809.