基于三角剖分算法的SPH粒子封闭曲面可视化方法与流程

本发明属于sph粒子封闭曲面可视化技术领域，具体涉及基于三角剖分算法的sph粒子封闭曲面可视化方法。

背景技术：

近年来随着人们对有限元法的成熟应用，网格数值方法中的一些不足之处也暴露出来。欧拉网格数值模拟的先决条件是要在问题域生成网格，但是在不规则或者很复杂的几何形状上构造规则的网格非常困难。而且还需要一些复杂的数学转换，这些转换有时甚至比问题本身更为繁琐。基于拉格朗日网格的方法在大变形求解中同样需要繁琐费时的网格再分法。

当需要研究大量离散物质点时候，近年来兴起的新一代计算方法——无网格法更加适用。其中光滑粒子法sph通过任意分布的粒子或节点求解具有各种边界条件的方程组得到精确数值解，这些粒子或节点不需要网格连接。sph是一种拉格朗日类型的无网格粒子法，具有非常好的自适应性，可以处理一些具有极大变形的问题，受到众多学者的追捧。

但是，sph方法由于没有固定关联的粒子，在后处理方面没有网格法处理简便直接。现在市场上主流的后处理软件都是基于网格的，不能直接由粒子进行模拟渲染，较好可视化更是遥遥无期。对于sph后处理的可视化来说，一条途径就是基于这些粒子建立合适的网格，将粒子信息转换为网格信息。这样的话就可以通过主流后处理软件得到较好的可视化效果。

delaunay三角剖分算法是一种简单直接的高效算法，离散粒子经过delaunay三角剖分后具有最接近、唯一性、最优性、最规则、区域性等优异特性。但是通过对近年来对各种sph后处理方法的研究，发现主要存在三个问题：1)空白单元过滤还主要停留在二维层面，没有发展到三维曲面上；2)三维粒子后处理主要停留在未封闭曲面剖分，对于封闭曲面主要靠凸包方法解决；3)对于粒子形成凹包的三角剖分，尤其是封闭凹包的曲面，没有较好的方法解决如何处理。

技术实现要素：

为克服上述现有技术的不足，本发明的目的是提供基于三角剖分算法的sph粒子封闭曲面可视化方法，具有

为实现上述目的，本发明采用的技术方案是：基于三角剖分算法的sph粒子封闭曲面可视化方法，包括以下步骤：

步骤1，离散粒子筛选表面粒子以及插值、拟合，具体做法是：利用sph粒子特有的光滑长度作为等值面的选取手段，由matlab的函数meshgrid对粒子生成网格采样点，并对于其中光滑长度hij在0.45～0.55的粒子选取为表面粒子；用delaunay三角剖分将这些表面粒子划分网格；对于网格化后的三角面片的法线坐标使用每个粒子点所连接的所有三角面的面法线平均值；发现筛选后的粒子数量较少，不足以划分出较为满意的网格，可在筛选粒子之前对粒子进行插值；

步骤2，粒子网格化后的可视化，具体做法是：

将处理后的粒子网格按照通用的obj文件格式保存，利用现有的网格法后处理软件实现可视化、渲染，可得。

所述的插值包括以下步骤：

首先，通过对粒子数据结果建立离散函数：

f(x,y,z)＝hij

其中x,y,z为粒子几何位置；

其次，在已有粒子的坐标范围按比例的增加坐标点，通过函数f插值出增加坐标点的hij值；然后再通过光滑长度对这些粒子筛选，就可以得到符合要求且数量可观的表面粒子；

最后，可以采用最近距离划分网格，利用粒子之间的距离来判断是否划分网格，得到液滴碰撞初始表面粒子网格化。

本发明的有益效果是：

本发明基于delaunay三角剖分，利用sph粒子的光滑长度筛选出有效表面粒子，应用插值、拟合以及求取等值面等办法，通过matlab、c++等软件解决了sph粒子在三维空间中的三角网格化而无需考虑sph粒子是否凸包、是否封闭的难点，使得sph粒子可以在较少的数据结果上剖分较为细腻的网格。这样就可以使用网格法的渲染手段来达到较好的可视化效果。

附图说明

图1为本发明实施例选取的液滴碰撞初始粒子数据图。

图2为本发明实施例的液滴碰撞初始粒子筛选后表面粒子数据图。

图3为本发明实施例的液滴碰撞初始表面粒子未插值的网格化图。

图4为本发明对比实施例液滴碰撞初始表面粒子插值后的图。

图5为本发明对比实施例采用凸包给粒子划分网格的图。

图6为本发明实施例液滴碰撞初始表面粒子插值后网格化的图。

图7为本发明实施例液滴碰撞初始粒子网格化后可视化图。

具体实施方式

下面结合附图和实施例对本发明作进一步详细说明。

参见图1，选取两液滴相撞初始状态的未处理粒子作为例子；液滴数据存储在“液滴碰撞.dat”中，以下为本实施例的液滴碰撞的获取步骤：

1)dat文件中可以直接获取粒子的坐标值和光滑长度；

2)已知n个粒子的坐标值x、y、z和光滑长度h；

3)分别求取x、y、z的最大值xmax、ymax、zmax和最小值xmin、ymin、zmin等距线性插值得到向量xvector、yvector、zvector；

xvector＝[xmin,x2,x3,x4……xn……xmax]；

yvector＝[ymin,y2,y3,y4……yn……ymax]；

zvector＝[zmin,z2,z3,z4……zn……zmax]；

4)利用这三个递增等距向量可以构建出一个n×n×n的矩阵，即vx,y,z＝[sx,y,1…sx,y,s…sx,y,z]；

5)对应的每一个我们需要对应的光滑长度h；

6)通过已有的n个粒子的坐标值x、y、z和光滑长度h进行插值，得到一个近似函数f(x,y,z)＝h；

7)通过函数f(x,y,z)＝h可以求得矩阵s中每个元素对应的光滑长度h；

8)通过筛选滑长度h在0.5附近(0.45-0.55)的元素值，这些粒子可以组成等值面；

9)由于这些粒子数量足够并且位置x、y、z都是单调递增且等距，可以通过让每个粒子连接距离自己最近的两个点组成三角形，划分三角网格，避免了delaunay算法中的只能凸包的障碍，也避免里alphashapes算法中的需要设置滚边长度的不确定性；

10)获取的划分好的网格，对于每个三角面片可以通过三点坐标求出该面片的法线，将法线归一化。通过三点距离面片中心距离来取权重，获得点的法线值。将点所连接的面的法线值求和平均即使该点的法线值。

基于三角剖分算法的sph粒子封闭曲面可视化方法，包括以下步骤：

步骤1，离散粒子筛选表面粒子以及插值、拟合，具体做法是：利用sph粒子特有的光滑长度作为等值面的选取手段，由matlab的函数meshgrid对粒子生成网格采样点，并对于其中光滑长度hij在0.45～0.55的粒子选取为表面粒子，区间越小精度越高，得到如图2所示的液滴碰撞初始粒子筛选后表面粒子数据，图2中所示粒子是将图1中光滑长度在0.45～0.55之外的粒子剔除后所剩例子；用delaunay三角剖分将这些表面粒子划分网格，参见图3，图3是液滴碰撞初始表面粒子(未插值)网格化。将图2中的粒子用alphashapes的方法联结成三角网格；对于网格化后的三角面片的法线坐标使用每个粒子点所连接的所有三角面的面法线平均值；发现筛选后的粒子数量较少，不足以划分出较为满意的网格，可在筛选粒子之前对粒子进行插值；

步骤2，粒子网格化后的可视化，具体做法是：

将处理后的粒子网格按照通用的obj文件格式保存，利用现有的网格法后处理软件实现可视化、渲染，得到如图7所示的液滴碰撞初始粒子网格化后可视化图。

所述的插值包括以下步骤：

首先，通过对粒子数据结果建立离散函数：

f(x,y,z)＝hij

其中x,y,z为粒子几何位置；

其次，在已有粒子的坐标范围按比例的增加坐标点，通过函数f插值出增加坐标点的hij值；然后再通过光滑长度对这些粒子筛选，就可以得到符合要求且数量可观的表面粒子；

在三维空间中对离散点划分曲面网格较为复杂，如果采用传统delaunay三角剖分就会立体的划分为有限数量个三棱锥，不符合要求。如果采用离散点拟合曲面会发现由于是封闭曲面，会在同一个x,y坐标值时存在两个z值；采用球坐标系的话会存在类似的问题；

如果采用凸包来求取边界的网格划分的话，发现如图4中粒子在中间部位存在凹陷，不符合凸包，划分出的网格不正确，如图5所示。图5就是采用凸包的方法求取网格，没有在应该凹陷的地方凹陷，得不到正确的网格。

如果对粒子我们采用滚边法来划分粒子的话，可以基本达到要求。但是滚边法要求设置一个最小半径来确定滚边的大小，不具有一般性。对于不同的粒子要求采取不同的滚边，不值得采纳；

最后，可以采用最近距离划分网格，利用粒子之间的距离来判断是否划分网格，得到液滴碰撞初始表面粒子网格化，得到如图6所示的视图。图6所示的网格与实际数据(dat文件)基本吻合，图中坐标单位皆为毫米。