复杂几何模型接触关系并行识别方法及应用

2022-06-14刘田田刘伟杰

计算机仿真 2022年5期

刘田田，刘伟杰，杨洋，郑澎

(1. 中物院高性能数值模拟软件中心，北京 100088；2. 北京应用物理与计算数学研究所，北京 100094；3. 中国工程物理研究院计算机应用研究所，四川绵阳 621900)

1 引言

高性能计算机能力的不断提升使得结构分析计算能针对数千部件的复杂系统进行模拟[1]。复杂系统的多部件之间通过粘接、铰接、铆接等方式装配在一起，存在着接触、摩擦等相互作用，获得零部件之间的相互作用是模拟计算过程中的重要一环[2]，也是决定结构分析结果准确可靠的关键之一[3，4]。

在进行接触分析计算时，需要先确定哪些部件间会存在接触关系。在工业生产部门中，数值模拟分析人员得到的来自于设计部门的几何模型往往缺少装配信息，对于多部件、装配关系复杂的系统级模型，不得不花费大量的时间人工寻找接触关系。因此，发展自动识别接触关系的方法显得尤为重要。近年来，基于几何信息的接触关系识别方法得到各界学者的广泛研究[5-14]。但是，目前已有的方法大都只适用于简单模型，局限于效率或精度的原因，对复杂模型和高阶曲面不适用。

面向复杂模型，接触关系识别算法达到实用化标准需要解决的两个关键问题是：提高识别精度和提高识别效率。模型中存在几何体自身的误差和几何体之间的误差等多种误差来源，选取合适的计算容差是提高算法识别精度的关键。已有的方法大都采用全局容差进行接触识别[10]，但全局容差的识别精度不高。已有的方法大都通过空间筛选技术提高算法效率[15]。

本文面向复杂B-Rep几何模型，提出了一种适用于任意曲面类型的接触关系识别算法，并将该算法应用到真实工程模型的结构分析实例中。算法基于几何模型的离散数据进行接触识别，根据模型的几何特征自动计算适用于相邻几何体和几何面之间的局部容差，并采用基于局部容差的碰撞检测技术提高识别精度。采用多线程并行技术，组合使用空间划分技术和包围盒筛选技术提高识别效率。

2 接触关系识别算法

接触表示两个几何体之间有某些面、线或点重合，但是没有体积上的重合[16]。几何体间的接触关系有六类：面-面接触、面-线接触、面-点接触、线-线接触、线-点接触和点-点接触。不同类型的接触识别方法都是类似的，由于计算时大部分接触约束是通过面联接的，因此，本文只介绍面-面接触的识别方法。

2.1 算法概述

B＿Rep几何模型通常由大量的三维几何体组成，包含多种类型的曲面。为了适用于包含任意曲面类型的模型，算法基于几何面的离散数据 (三角面片) 进行接触探测。目前，商业和开源几何引擎都提供B＿Rep模型的离散化接口，获取离散数据非常简单。

对模型中的任意两个几何面，判断是否存在接触关系的方法是对两个几何面的离散数据进行碰撞检测。复杂工程模型的尺寸特征跨度大，为了提升算法精度，采用基于局部容差的碰撞检测。由于复杂装配模型中往往存在数千上万个几何面，离散数据可能包含数百万个三角面片，直接进行三角形接触检测计算量太大。本文采用两种策略提高接触关系的识别效率：

1)采用空间划分和包围盒筛选相结合的技术筛选几何面对；

2)采用多线程并行技术加速筛选和接触检测过程。

为方便算法描述，记{Bi,i=1,…，n}为模型中所有三维几何体的集合，{Fij，j=1,…,m}为第i个几何体包含的所有几何面的集合。Box(Bi)表示几何体Bi的包围盒，Box(Fij)表示几何体Bi中第j个面Fij的包围盒。Box(Bi)±ε表示Bi的包围盒向各方向外扩ε。

2.2 局部容差计算

局部容差是根据模型的局部特征长度计算得到的。几何体的特征长度定义为包围盒的最短边长

L(B)=min{lx，ly，lz}，

(1)

其中lx，ly，lz分别表示包围盒在x，y，z方向的边长。几何面的特征长度定义为组成该面的最短边长度或包围盒的对角线长度

(2)

其中k为该面包含的边个数。

定义几何体的特征容差εb如下

εb=max{L(Bi)/R，i=1，…，n}，

(3)

其中n表示模型包含的几何体个数，R是容差计算因子，取值大于0。R越大，局部容差越小，根据经验，R默认取值20。

定义几何面的特征容差εf如下

εf=max{L(Fj)/R，j=1，…，m}

(4)

其中m表示模型包含的几何面个数。

定义两个相邻几何面之间的局部容差如下

tol(Fj1，Fj2)=min{L(Fj1)/R，L(Fj2)/R}

(5)

2.3 面-面接触检测

记Fj1和Fj2是两个待检测的几何面，Sj1和Sj2分别为两几何面的离散数据，两面之间的局部容差根据式(5)计算得到。判断两个几何面之间是否存在接触关系，只需判断两个集合面的离散数据中是否存在三角形共面相交，即存在s1∈Sj1，s2∈Sj2满足

①s1和s2的法向相反；

②s1和s2在容差范围内共面相交 (或包含)。

三角形之间的相交检测是图形学中的经典问题，有多种求解方法，如分离轴算法[17]、Moller算法[18]、Devillers算法[19]、Shen算法[20]、Tropp算法[21]等。采用的是分离轴算法。为了减少不必要的碰撞操作，对每个几何面的三角面片构建k-d树，基于k-d树进行邻域内的三角形搜索。

2.4 空间划分+包围盒分层筛选技术

由于几何模型的零件之间是紧密排列的，且进行的碰撞检测是静态的，为了提高效率，算法选用空间划分和包围盒分层筛选相结合的方法来减少面-面之间的接触检测。根据B-Rep模型的层级关系，如果两个几何体的包围盒不相交 (含相切)，那么这两个几何体之间一定不存在接触面对。如果两个几何面的包围盒不相交，那么这两个几何面一定不接触。几何面对筛选的算法如下：

Step 1：计算所有几何体、几何面的包围盒。

Step 2：对所有几何体的包围盒构建k-d树，记为BTree。对第每个几何体包含的几何面的包围盒构建k-d树，记为FTree(i)， i取值1～n。

Step 3：对任一几何体Bi，在BTree中查找与Box(Bi)±εb相交的包围盒列表{Box(Bik),k=1,…,K}。遍历Bik，为避免重复计算，若Bik的实体编号小于等于Bi，则跳过；否则，执行Step 4。

Step 4：对OBi的任一几何面Fij，在FTree(ik)中查找与Box(Fij)±εf相交的包围盒列表{Box(Fikl),l=1,…,L}。若Fikl的实体编号小于等于Fij，则跳过；否则，Fij与Fikl时可能存在接触关系的面对，需进行面-面接触检测。

图1 接触算法并行设计流程

2.5 多线程并行

串行算法难以满足复杂模型在计算效率方面的应用需求。目前的处理器都是多核系统，采用多线程分配任务的策略是提高算法性能的有效手段。本文采用OpenMP实现多线程并行算法。

算法中可以被并行分解的任务内容包括：包围盒计算、几何特征长度计算、多个k-d树的构建、包围盒筛选和面对接触检测。因此，并行设计流程如图1所示，其中I、II、III、IV表示循环。由于几何体包围盒筛选、几何面包围盒筛选和接触检测之间存在层级依赖，为了充分利用并行资源，将并行粒度放在几何体包围盒筛选的循环上。

3 算法分析

利用图2中的装配模型对算法效率和算法精度进行了测试，该模型包含76个体， 14467个面，含多种类型的曲面，模型内包含505对面-面接触关系。如果对所有的几何面进行两两接触检测，则需要检测209，294，089次。通过包围盒筛选后，只需要检测557，367次，检测次数降到了0.27%。几何面对筛选和面-面接触检测过程中，算法采用k-d树进行邻近搜索， k-d树的构建复杂度是O(NlogN)，搜索的复杂度是O(logN)，大大提高了计算效率。

图2 公开发动机模型

利用该模型对算法的并行效率进行分析，分别采用1/2/4/8线程进行测试，结果如表1所示。随着线程数增加，计算时间逐渐下降。8核并行的运行时间是串行时间的1/3左右，表明多线程并行的有效性。

表1 不同线程数目算法运行时间

由于识别方法主要基于几何信息，并未考虑物理模型，因此对于识别结果，还需要分析人员进行手工剔除，把需要被剔除的接触关系称为“冗余接触关系”。分别测试了本文的接触识别方法和商业软件的自动识别方法，结果如表2所示。测试环境CPU型号为Intel(R) Xeon(R) CPU E5-2680 v4 @ 2.4GHz， 8核。从表中看出，采用本文的方法识别结果更精确，只识别出少量的冗余接触，效果好于商业软件ABAQUS和ANSYS。分析原因，主要是因为本文的识别方法是基于局部容差进行的，而商业软件采用一致化的全局容差，为识别最大特征，不得不设置最大容差值，从而识别出大量的冗余接触。对于模型尺寸跨度比较大的复杂模型，选取合适的全局容差非常困难，使用基于几何特征的局部容差计算效果明显优于全局容差。从表中的运行时间来看，本文中的算法效率略低于ANSYS，但明显优于ABAQUS。

表2 不同方法识别的接触关系结果对比

4 结构分析应用实例

本节以某土工离心机模型的静强度分析为例验证接触关系算法的有效性。离心机是一种利用旋转来模拟特定高重力场加速度环境的试验设备，广泛地应用于科研和国民经济建设的多个领域。土工离心机为岩土离心模型试验提供了重要的研究手段，广泛用于堤坝、边坡地基、挡土结构、路基、隧道、重要建筑物基础等方面的工程应用研究，是评估岩土工程设计结果的合理性、安全性和经济性等重要设备。

图3中的离心机模型，包含126个几何实体， 821个几何面。采用本文的接触关系识别算法，共检测到393对面-面接触关系，如图5所示，图中红色表示面-面接触。