一种用于光滑粒子流体动力学模拟出血及处理加速的方法技术

一种用于光滑粒子流体动力学模拟出血及处理加速的方法，利用基于GPU加速的模拟方法来模拟虚拟手术中的出血及利用虹吸原理吸出渗血的效果；同时，利用网格法实时划分问题区域，创建以支持域为边长的空间网格，通过临近网格搜索最近相邻粒子，并且通过并行计算架构(CUDA)多线程并行加速技术完成粒子控制方程的求解以及血液与固体交互的计算，大大提高了运算效率，从而提高了手术训练的实时性。此外；改进的移动立方体算法(marching cube)被用于流体表面的渲染，大大提高了手术训练的真实性。

A method of haemorrhage and acceleration for fluid dynamics simulation of smooth particles

【技术实现步骤摘要】
一种用于光滑粒子流体动力学模拟出血及处理加速的方法
本专利技术涉及计算机图形学下的医学成像、流体动力学领域，尤其是一种基于GPU加速的流体生成及在需要擦除时利用虹吸原理对渗血吸除的模拟分析方法。
技术介绍
随着计算机硬件的发展，医学成像技术的进步，虚拟现实在医学领域中占有越来越大的影响力，拥有先进的虚拟现实技术是这个时代不可或缺的标志。虚拟现实技术是一种通过计算机图形学等学科高度仿真模拟现实的技术，在医学领域中，目前正被广泛应用于医学仿真训练的研究。虚拟现实技术的模拟真实度主要取决于模型的仿真效果以及运行的流畅度。虚拟训练系统中的对血液的模拟在传统上采用光滑粒子流体动力学(SPH)方法。光滑粒子流体动力学将流体离散成带有属性和体积的粒子，再结合有限元的思想为每个粒子属性求解流体力学方程得到下一时刻的状态，粒子的属性不是单独存在的，是通过对周围光滑核函数半径范围内的粒子属性值使用核函数加权求和得到。然而，传统血液模拟的质量一直以来都不是很理想，且由于算法计算量偏大使得成像的反应时间较长，导致模拟的效果难以满足实时性的要求。
技术实现思路
为了解决现有系统中血液的计算量大而导致的实时性和仿真效果的不理想的问题，本专利技术使用CUDA指令集来作为GPU加速的技术来提高计算的速度。同时采用一种改进的移动立方体方法对血液表面进行实时可视化渲染，增强了血液模拟的真实性。为了解决上述技术问题提供了如下技术方案：一种用于对光滑粒子流体动力学血液模拟过程加速的方法，所述方法包括以下步骤：步骤一，使用PCISPH方法对粒子建模，过程如下：1.1使用平坦和归一化内核得到的二阶精度插值设计核函数；1.2由核函数计算出单个粒子所受的压力、粘性力、外部作用力，再根据牛顿第二定律计算出合力；1.3通过各粒子的内部和外部的受力分析对出血做出模拟；1.4根据对当前压力计算下一步速度和位置；步骤二，基于CUDA对出血模拟的加速；步骤三，对流体表面渲染，过程如下：每进行一步的SPH模拟都要将点云转化为可着色的表面，将密度场均匀划分为大小相等的格子，使用MarchingCubes算法对格子的胞元的每个角点的密度场交点进行评估，根据其是否大于某个常量阈值来决定生成三角形的形式；使用PN三角形方法来平滑三角形网格边缘并为顶点着色操作提供更多的采样点，根据三角形的三个顶点及其法线就生成平滑的细分三角形。进一步，在将已经生成的血液进行清除操作的模拟时，需要将吸引器管口附近的血液粒子按照距离管口距离的大小施加大小不同指向管口的吸引力，使得粒子以不同的速度向管口运动，且在运动到管口处与管口相碰撞后使该粒子“消失”。由于该过程中设计到判断有效范围，并在有效范围内分多个区域，每个区域会有不同的力施加给粒子等相关运算，计算量非常大，因此也采用CUDA进行加速运算；所述方法还包括以下步骤：步骤四确定粒子系统中在该吸管吸引范围中的粒子，为吸引范围内不同的区域设置不同的吸引力以吸引器的管口为中心，根据设定的吸管吸引力的大小确定吸引半径，做出的吸引范围，在该范围内吸管才能做出有效的吸引，超出该范围则是无效吸引；按照吸管吸引的反方向为法线设置不同的角度设置不同的吸引力大小；步骤五将吸到吸引器管口位置的粒子“吸掉”，对吸血过程基于CUDA加速将运行到吸引器管口位置并与该管口相碰撞的粒子的生命周期设为0，删除该粒子，注销该粒子的内存空间；由于吸引范围内有不同区域内吸引力大小不同使得粒子在不同的区域会有不同的指向吸引器管口的加速度。再进一步，所述CUDA的计算包含以下参数：a.哈密顿算子▽：对粒子的受力情况分析以及粒子速度加速度的计算时需对其所在的矢量场进行分析，因此在设计核函数的时候得计算出矢量场，即粒子所在位置的“梯度”，用“梯度”来表示标量场中粒子的变化快慢和方向；在计算加速度时还需使用拉普拉幸算子，即▽·▽或▽2；b.Helmholtz-Hodge定理：任意的矢量场都可以唯一的分解成一个零散度的矢量场和一个标量场的梯度之和，如下所示：其中u是零散度矢量场，p是标量场。对上述公式求散度得：因为u的散度为零，所以得压强泊松方程：c.扩散项的CUDA映射内核函数对计算部分加速迭代：该内核函数对应的Grid线程结构被设计为二维成结构，其在x维度上的block数量为：(M+blockDim.x-1)/blockDim.x在y维度上的block数量为：(M+blockDim.y-1)/blockDim.yd.对流项的计算：本专利技术采用的是JosStam提出的隐式方法求解对流项，即从当前时刻起，根据粒子的运动轨迹，计算出上一时刻粒子的位置，并将该位置粒子的q量赋值给现在的这个粒子，其计算公式为:q(x,t+δt)＝q(x-u(x,t)δt,t)其中，q(x,t+δt)是t+δt时刻粒子的位置，u(x,t)表示x位置处t时刻的速度；e.Navie-Stokes方程组：用于描述流体的运动，当描述牛顿流体的不可压缩流体时，方程组定义的方程可被简化为不可压缩的质量守恒方程和动量守恒方程：其中ρ是密度，p是压强，v是速度，μ是粘性系数，fext是所有外力之和；f.邻近区域内粒子插值计算属性：流体的模拟是一种基于拉格朗日方法的流体建模方法，将流体离散为大量的粒子，是一种基于粒子的插值方法，在位置ri处的粒子i的一个标量A(ri)由该粒子邻近区域内粒子的属性通过插值得到：其中，n是领域内粒子的数量，j是粒子的索引号，Aj是粒子粒子j的属性值，mj是粒子j的质量，ρj是粒子j的吗，密度，ri是粒子j的位置，W(ri-rj,h)是插值权重函数在SPH方法中称为光滑核函数，h是光滑核函数W影响区域的半径，在支持域以外光滑核函数W被定义为零；依据该插值计算可依次求得标量A(ri)的梯度▽A(ri)、拉普拉新算子▽2A(ri)，进而求得粒子的各属性值；h.圆管流雷诺数：表针粘性影响的相似准数，表征流体流动情况的无量钢数，雷诺数大小决定了粘性流体的流动特性，u为流速，v为运动粘度，d为吸管直径。更进一步，所述CUDA是由CPU和GPU协同完成的，CPU作为主机端，主要处理逻辑性较强的事务及串行程序的计算，GPU扮演设备端通过kernel函数并行处理负责密集型大数据的并行运算；所需要的出血模拟是一个由粒子系统模拟的模型，在该模拟中将大量的重复的运算放置到GPU端运算处理，CUDA的每一个线程都有一个唯一的线程ID，所以CUDA线程可使用一维、二维、三维索引标识，从而形成一维、二维、三维线程块，线程块包含在线程格(grid)中，线程格的尺寸根据粒子系统的规模来指定。所述步骤二中，基于CUDA对出血模拟加速的过程如下：2.1初始化数据，获取流体的外力大小，创建密度场和速度场的存储空间及其相关的初始化；2.2数据传输，将相关粒子属性等数据拷贝到显存中，为GPU上流体的计算提供数据；2.3在GPU上实现模拟，包括了扩散项、外力项、对流项、投影过程中各个步骤的计算，从而得到每一个时间步长后流体的密度场和密度场，计算出下一个时间步长后的速度和位置；2.4将GPU计算的速度场和密度场渲染到屏幕的窗口中去。2.5边界的处理速度在边界处遵循无滑条件，且密度保持不变，在水平方向上的流体边界条件为：水平边界1：Dens(i,0)本文档来自技高网...

【技术保护点】
一种用于对光滑粒子流体动力学出血模拟过程加速的方法，其特征在于：所述方法包括以下步骤：步骤一，使用PCISPH方法对粒子建模，过程如下：1.1使用平坦和归一化的内核设计核函数；1.2由核函数计算出单个粒子所受的压力、粘性力、外部作用力，再根据牛顿第二定律计算出合力；1.3通过各粒子的内部和外部的受力分析对出血做出模拟；1.4根据当前压力计算下一步速度和位置；步骤二，基于CUDA对出血模拟的加速；步骤三，对流体表面渲染，过程如下：每进行一步的SPH模拟都要将点云转化为可着色的表面，将密度场均匀划分为大小相等的格子，使用Marching Cubes算法对格子的胞元的每个角点的密度场交点进行评估，根据其是否大于某个常量阈值来决定生成三角形的形式；使用PN三角形方法来平滑三角形网格边缘并为顶点着色操作提供更多的采样点，根据三角形的三个顶点及其法线就生成平滑的细分三角形。

【技术特征摘要】
1.一种用于对光滑粒子流体动力学出血模拟过程加速的方法，其特征在于：所述方法包括以下步骤：步骤一，使用PCISPH方法对粒子建模，过程如下：1.1使用平坦和归一化的内核设计核函数；1.2由核函数计算出单个粒子所受的压力、粘性力、外部作用力，再根据牛顿第二定律计算出合力；1.3通过各粒子的内部和外部的受力分析对出血做出模拟；1.4根据当前压力计算下一步速度和位置；步骤二，基于CUDA对出血模拟的加速；步骤三，对流体表面渲染，过程如下：每进行一步的SPH模拟都要将点云转化为可着色的表面，将密度场均匀划分为大小相等的格子，使用MarchingCubes算法对格子的胞元的每个角点的密度场交点进行评估，根据其是否大于某个常量阈值来决定生成三角形的形式；使用PN三角形方法来平滑三角形网格边缘并为顶点着色操作提供更多的采样点，根据三角形的三个顶点及其法线就生成平滑的细分三角形。2.如权利要求1所述的一种用于对光滑粒子流体动力学模拟血液清除过程加速的方法，其特征在于：在将已经生成的血液进行清除操作的模拟时，需要将吸引器管口附近的血液粒子按照距离管口距离的大小用不同的速度向管口运动，且在运动到管口使该粒子“消失”，计算过程也采用CUDA进行加速；所述方法还包括以下步骤：步骤四确定粒子系统中在该吸管吸引范围中的粒子，为吸引范围内不同的区域设置不同的吸引力以吸引器的管口为中心，根据设定的吸管吸引力的大小确定吸引半径，做出的吸引范围，在该范围内吸管才能做出有效的吸引，超出该范围则是无效吸引；按照吸管吸引的反方向为法线设置不同的角度设置不同的吸引力大小；步骤五将吸到吸引器管口位置的粒子“吸掉”，对渗血吸除过程基于CUDA加速将运行到吸引器管口位置并与该管口相碰撞的粒子的生命周期设为0，删除该粒子，注销该粒子的内存空间；由于吸引范围内有不同区域内吸引力大小不同使得粒子在不同的区域会有不同的加速度。3.如权利要求2所述的一种用于光滑粒子流体动力学模拟的出血及处理加速的方法，其特征在于：所述CUDA的计算包含以下参数：a.哈密顿算子▽：对粒子的受力情况分析以及粒子速度加速度的计算时需对其所在的矢量场进行分析，因此在设计核函数的时候得计算出矢量场，即粒子所在位置的“梯度”，用“梯度”来表示标量场中粒子的变化快慢和方向；在计算加速度时还需使用拉普拉幸算子，即▽·▽或▽2；b.Helmholtz-Hodge定理：任意的矢量场都可以唯一的分解成一个零散度的矢量场和一个标量场的梯度之和，如下所示：w＝u+▽p其中u是零散度矢量场，p是标量场。对上述公式求散度得：▽·w＝▽·u+▽2p因为u的散度为零，所以得压强泊松方程：▽2p＝▽·wc.扩散项的CUDA映射内核函数对计算部分加速迭代：该内核函数对应的Grid线程结构被设计为二维成结构，其在x维度上的block数量为：(M+blockDim.x-1)/blockDim.x在y维度上的block数量为：(M+blockDim.y-1)/blockDim.yd.对流项的计算：本发明采用的是JosStam提出的隐式方法求解对流项，即从当前时刻起，根据粒子的运动轨迹，计算出上一时刻粒子的位置，并将该位置粒子的q量赋值给现在的这个粒子，其计算公式为:q(x,t+δt)＝q(x-u(x,t)δt,t)其中，q(x,t+δt)是t+δt时刻粒子的位置，u(x,t)表示x位置处t时刻的速度；e.Navie-Stokes方程组：用于描述流体的运动，当描述牛顿流体的不可压缩流体时，方程组定义的方程可被简化为不可压缩的质量守恒方程和动量守恒方程：其中ρ是密度，p是压强，v是速度，μ是粘性系数，fext是所有外力之和；f.邻近区域内粒子插值计算属性：流体的模拟是一种基于拉格朗日方法的流体建模方法，将流体离散为大量的粒子，是一种基于粒子的插值方法，在位置ri处的粒子i的一个标量A(ri)由该粒子邻近区域内粒子的属性通过插值得到：其中，n是领域内粒子的数量，j是粒子的索...

【专利技术属性】
技术研发人员：冯远静，唐良甫，张弛，肖昌盛，朱珂权，曾庆润，曹喆文，
申请(专利权)人：浙江工业大学，
类型：发明
国别省市：浙江,33