[发明专利]基于完全矩阵计算的离散元邻居搜索及求解方法和系统

说明书

本发明公开了一种基于完全矩阵计算的离散元邻居搜索及求解方法和系统，将空间划分为正交网格元胞；通过完全矩阵运算，基于网格元胞的邻域元胞建立初步接触搜索矩阵，并通过距离判别筛选接触信息及虚单元填充方法，得到接触搜索的邻居矩阵表示；基于邻居矩阵及相应布尔矩阵实现离散元程序的求解计算。本发明摒弃了循环语句结构，实现了离散元邻居搜索与求解计算的纯矩阵运算，进而通过GPU并行加速技术大幅提高三维离散元计算效率。

技术领域

本发明涉及一种邻居搜索及求解计算的纯矩阵实现方法，尤其适用于颗粒离散元GPU并行加速的高效计算方法。

背景技术

自1979年Cundall提出离散元法以来，离散元法因其天然的非连续特性被广泛应用到颗粒物质材料及大变形问题。离散元法基于牛顿第二定律，通过时间步迭代算法显示求解，因此随着计算的空间尺度和时间尺度的增加，计算单元数和时间步迭代次数逐渐增加，离散元法的计算量也显著增加。其中三维接触搜索与判别是颗粒离散元计算中最关键的问题之一，其计算量占比可达到总计算量的50％到80％。特别是对于大变形问题(如滑坡)、动力过程，如何在数百万单元运动过程中，实时高效地实现三维接触搜索判别一直是制约离散元法在工程尺度推广与应用的瓶颈。

近年来基于GPU的高性能计算逐渐兴起，而传统的离散元程序往往以串行代码为基础，以链表的形式存储接触信息，不利于代码的GPU并行计算。

发明内容

发明目的：为了克服现有技术中存在的不足，本发明提供一种基于完全矩阵计算的离散元邻居搜索及求解方法和系统，本发明通过虚单元填充方法以邻居矩阵的形式存储接触信息，实现了颗粒离散元程序的纯矩阵运算，便于程序编写与GPU并行加速，实现了高效的离散元邻居搜索及求解计算。

技术方案：为实现上述目的，本发明采用的技术方案为：

一种基于完全矩阵计算的离散元邻居搜索方法，包括以下步骤：

步骤1，单元编号及信息存储：将单元顺序编号，记从1到m，同时添加虚单元，其编号m+1，通过编号唯一标识单元；通过虚单元填充得到接触搜索的邻居矩阵表示，邻居矩阵为m行，每一行记录单元可能接触的单元编号，长短不一部分用虚单元编号填充，单元信息存储为单元属性数组，记为aX、aY、aZ和aR，其中，编号m+1单元为虚单元，aX、aY、aZ为三维坐标，aR为半径；

步骤2，网格划分及编号：

步骤21，输入不包括虚单元的单元三维坐标及半径数组，记三维坐标数组为X,Y,Z，半径数组为R，其均为m行的一维数组；

步骤22，根据三维坐标数组确定网格边界及边长：

p方向网格两侧边界：

Pmin＝min(P)-minR/2

Pmax＝max(P)+minR/2

其中，[Pmin，Pmax]为p方向两侧边界，[min(P),max(P)]为p方向坐标数组最小最大值，minR为最小单元半径，P表示x,y,z方向；

记[Xmin,Xmax],[Ymin,Ymax],[Zmin,Zmax]为x,y,z方向网格两侧边界；

网格元胞边长：

gSide＝2*maxR+dSide

其中，gSide为网格元胞边长，maxR为单元最大半径，dSide为可调值；

以最小值边界为起点，gSide为网格元胞边长，划分空间网格；则p方向元胞数PNum：

其中，ceil表示向上取整运算，元胞沿p方向顺序编号，记从1到PNum；