专利名称：半隐式类粒子法的大规模并行处理方法
技术领域：
本发明属于流体流动计算机仿真领域，涉及对不可压缩流体的自由表面流动的模拟，具体涉及一种半隐式类粒子法的大规模并行处理方法。
背景技术：
自由表面流动是一种常见的流动，在溃坝流动、河流的明渠流动，以及蒸发沸腾过程的传热传质过程中都涉及到自由表面流动，采用计算机准确的模拟自由表面流动对工程应用有重大意义；粒子类方法是在拉格朗日框架下追踪粒子的运动，因此易于仿真自由表面流动。随着计算机技术的发展，采用计算机模拟流体流动已经成为研究流体运动的一种 流动的时候需要引入专门的处理界面运动的方法，如VOF方法和Level Set方法。对于剧烈波动的自由表面流动问题，大变形自由表面的处理是很复杂且很耗时的工作，并且有时会导致求解的错误；另外，网格类方法中对流项的数值耗散也是一个很严重的问题。粒子类方法是近年来快速发展的一种算法，它们采用的是在拉格朗日框架下追踪粒子的运动，可以追踪自由表面的位置，因此易于仿真自由表面流动。光滑粒子动力学(SPH)和移动粒子半隐式算法(MPS)是两种常见的粒子类方法。前者采用显式算法求解压力，主要用于求解固体力学问题以及可压缩的流体力学问题；而后者采用半隐式算法通过求解压力泊松方程获得压力，主要用于不可压缩性流动的模拟。通常来说，粒子类方法的计算量是大于网格类方法的计算量，同时工程应用领域也要求大规模粒子数的模拟，因此发展高效大规模并行的粒子类方法具有重大的现实意义。显式的粒子类方法(如SPH)已经得到了较好的并行。本发明主要是针对半隐式类粒子法(如MPS)的并行化。
发明内容
粒子类方法的算法可以分为显式算法和半隐式算法。在显式类算法(如SPH算法)中，粒子的压力是通过状态方程由该粒子所代表的流体微团的密度计算得到，每一个粒子压力的计算是独立于其他粒子的，因此显式类粒子方法是很容易并行的；而在半隐式类粒子算法(如MPS)中，压力是通过全局求解压力泊松方程而获得，因此每一个粒子压力的计算都是依赖于其他粒子的信息的，并且全局求解压力泊松方程是半隐式类粒子法中最为耗时的部分，因此半隐式类粒子法并行的关键在于找到一个高效且易于并行的用于求解压力泊松方程的求解器。本发明的目的在于提供一种半隐式类粒子法的大规模并行处理方法，该方法先给出一个高效且易于并行的求解器，然后通过并行实现大规模粒子数的模拟。为达到上述目的，本发明采用了以下技术方案该处理方法在并行服务器上采用并行模型并行加速半隐式类粒子法，所述半隐式类粒子法是指黏性项和外力源项显式计算，而压力项则通过求解压力泊松方程隐式计算的粒子类算法，半隐式类粒子法中的压力泊松方程采用对称兰乔斯算法求解。所述半隐式类粒子法包括移动粒子半隐式算法(Moving ParticleSemi-implicit method, MPS)和不可压缩的光滑粒子动力学方法(IncompressibleSmoothed Particle Hydrodynamics, ISPH)。所述并行服务器为对称多处理器、GPU服务器、由多CPU核的服务器组成的多节点的服务器、由多GPU核的服务器组成的多节点的服务器或由多CPU核和多GPU核服务器组成的混合多节点服务器。所述并行模型为共享内存的并行模型、消息传递的并行模型或二者的混合模型。所述并行加速具体包括以下几种优选情况
a)如果并行服务器是对称多处理器，采用OpenMP模型(http://www. openmp.org)、MPI模型(张武生，薛巍，李建江，郑伟民《MPI并行程序设计实例教程》清华大学出版社，北京，2009)或采用OpenMP和MPI的混合模型来并行加速；0penMP模型是基于共享内存的并行模型，由开源网站维护，网站上给出了具体的编程方法；MPI模型是基于消息传递的并行模型；b)如果并行服务器是GPU服务器,通过OpenCL模型(http://www. khronos. org/opencl/)或 CUDA 模型(http://developer, nvidia. com/category/zone/cuda-zone)并行加速；0penCL模型适用于GPU的开源并行模型，由开源网站维护；CUDA是由Nvidia公司开发使用于GPU的并行模型，在公司官网上给出了具体的编程方法；c)如果并行服务器是多节点的服务器，每个节点是多CPU核的服务器，在每一个节点上，采用OpenMP模型并行加速，而在节点之间采用MPI模型来并行加速；d)如果并行服务器是多节点的服务器，每个节点是多核的GPU服务器，在每一个节点上，采用OpenCL模型或CUDA模型并行加速，而在节点之间采用MPI模型来并行加速；如果节点中既有多核的CPU节点，又有多核GPU节点，对CPU节点采用OpenMP模型并行加速，对GPU节点采用OpenCL模型或CUDA模型并行加速，在节点之间仍采用MPI模型并行加速。所述处理方法的具体步骤如下I)设置计算区域内粒子的初始分布，将计算区域划分为多个子区域，将每一个子区域分配给一个进程负责计算；2)将进程上靠近划分边界的粒子信息发送给邻居进程，并且邻居进程上要完成相应的接受操作；3)在每一个进程上显式计算外力源项和黏性项，然后计算临时速度场和临时位置场，从而计算压力泊松方程的源项；步骤3)的操作是针对子区域上的每一个粒子在进程内并行加速；4)采用并行对称兰乔斯算法全局求解压力泊松方程，获得压力值；5)将获得的压力值发送到邻居进程上；6)传送压力值的通信完成之后，并行计算压力梯度项，用于修正临时速度场，从而获得最终的速度场；7)计算过程中要动态调节区域划分边界，以保证每个进程上的负载都近似相同；
8)判断本进程上由于粒子移动产生的粒子越界的情况，将已经越过本进程边界的粒子信息发送给相应的进程；9)判断是否已达到模拟设定时间，如果已达到模拟设定时间则退出；如果没有达到模拟设定时间，返回步骤2)，继续下一个时间步长的计算。粒子法的压力泊松方程离散后通常得到一个主对角占优的大型对称稀疏矩阵，不完全乔莱斯基分解算法(ICCG)是求解大型对称稀疏矩阵的理想迭代方法，但是难以并行。为了实现粒子法的大规模并行运算，需要引入一种高效且易于并行的求解线性方程组的方法。通过多次尝试，本发明最终采用对称兰乔斯算法(SLA)求解压力泊松方程，这种算法易于并行，并且其求解效率优于ICCG算法，并行加速效果显著。
图I为背景网格的划分以及中心粒子的作用范围。 图2为子区域划分。图3为不同子区域(进程)间相互通信。图4为溃坝问题的模型。图5为在求解溃坝问题中压力泊松方程时SLA算法和ICCG算法单次求解执行时间的比较。图6为SLA算法和ICCG算法中易于并行的操作所占的百分比。图7为OpenMP模型在不同CPU核数下的并行效率。图8为MPI模型在不同节点数下的并行效率。
具体实施例方式下面结合附图和实施例对本发明作进一步说明。SLA算法求解压力泊松方程的步骤如下SLA算法实际上是一种Krylov子空间方法,设m是Krylov子空间的维度,对称兰乔斯算法求解对称矩阵方程组Ax=b的具体算法如下I)给定误差极限ε >0,初始解向量x° e Rn;2)计算 r°=b-Ax°, β =I r0 , N1=V0/ β ;3) for k = Ito n doβ!=0for i=lto m doa ^(Avi, Vi)Vt=Ayi- a ^i- β JV^1if (i〈m) β i+1= | vt , vi+1=vt/ β i+1End for4)用追赶法求解如下的m维的三对角方程组
权利要求
1.一种半隐式类粒子法的大规模并行处理方法，其特征在于该处理方法在并行服务器上采用并行模型并行加速半隐式类粒子法，半隐式类粒子法中的压力泊松方程采用对称兰乔斯算法求解。
2.根据权利要求I所述一种半隐式类粒子法的大规模并行处理方法，其特征在于所述半隐式类粒子法包括移动粒子半隐式算法和不可压缩的光滑粒子动力学方法。
3.根据权利要求I所述一种半隐式类粒子法的大规模并行处理方法，其特征在于所述并行服务器为对称多处理器、GPU服务器、由多CPU核的服务器组成的多节点的服务器、由多GPU核的服务器组成的多节点的服务器或由多CPU核和多GPU核服务器组成的混合多节点服务器。
4.根据权利要求I所述一种半隐式类粒子法的大规模并行处理方法，其特征在于所述并行模型为共享内存的并行模型、消息传递的并行模型或二者的混合模型。
5.根据权利要求1、2、3或4所述一种半隐式类粒子法的大规模并行处理方法，其特征在于所述并行加速具体包括以下几种优选情况 a)如果并行服务器是对称多处理器，采用OpenMP模型、MPI模型或采用OpenMP和MPI的混合模型来并行加速； b)如果并行服务器是GPU服务器，通过OpenCL模型或⑶DA模型并行加速； c)如果并行服务器是多节点的服务器，每个节点是多CPU核的服务器，在每一个节点上，采用OpenMP模型并行加速，而在节点之间采用MPI模型来并行加速； d)如果并行服务器是多节点的服务器，每个节点是多核的GPU服务器，在每一个节点上，采用OpenCL模型或CUDA模型并行加速，而在节点之间采用MPI模型来并行加速；如果节点中既有多核的CPU节点，又有多核GPU节点，对CPU节点采用OpenMP模型并行加速，对GPU节点采用OpenCL模型或CUDA模型并行加速，在节点之间仍采用MPI模型并行加速。
6.根据权利要求I所述一种半隐式类粒子法的大规模并行处理方法，其特征在于所述处理方法的具体步骤如下 1)设置计算区域内粒子的初始分布，将计算区域划分为多个子区域，将每一个子区域分配给一个进程负责计算； 2)将进程上靠近划分边界的粒子信息发送给邻居进程，并且邻居进程上要完成相应的接受操作； 3)在每一个进程上显式计算外力源项和黏性项，然后计算临时速度场和临时位置场，从而计算压力泊松方程的源项；步骤3)的操作是针对子区域上的每一个粒子在进程内并行加速； 4)采用并行对称兰乔斯算法全局求解压力泊松方程，获得压力值； 5)将获得的压力值发送到邻居进程上； 6)传送压力值的通信完成之后，并行计算压力梯度项，用于修正临时速度场，从而获得最终的速度场； 7)计算过程中要动态调节区域划分边界，以保证每个进程上的负载都近似相同； 8)判断本进程上由于粒子移动产生的粒子越界的情况，将已经越过本进程边界的粒子信息发送给相应的进程； 9)判断是否已达到模拟设定时间，如果已达到模拟设定时间则退出；如果没有达到模拟设定 时间，返回步骤2)，继续下一个时间步长的计算。
全文摘要
本发明提供一种半隐式类粒子法的大规模并行处理方法1)采用高效且易并行的对称兰乔斯算法(Symmetric Lanczos Algorithm,SLA)求解半隐式粒子法中压力泊松方程；2)在服务器上可采用共享内存的并行模型(如OpenMP、OpenCL)和/或消息传递的并行模型(MPI)以及它们的混合模型并行加速算法，本发明可以大规模并行加速需要求解压力泊松方程的半隐式类粒子法，并行加速效果显著。
文档编号G06F17/50GK102902514SQ20121033105
公开日2013年1月30日 申请日期2012年9月7日 优先权日2012年9月7日
发明者陈斌, 段广涛 申请人:西安交通大学