约束数据域Delaunay四面体网格生成算法

1、 收稿日期:2004207213.作者简介:关文革(19662 , 女, 博士研究生; 北京, 中国矿业大学资源与安全学院(100083 .E 2m ail :gwg -007基金项目:教育部跨世纪优秀人才基金资助项目(200023 ; 教育部青年骨干教师基金资助项目(2000265 ; 河北省科技厅攻关资助项目(032135125 .约束数据域Delaunay 四面体网格生成算法关文革1,2武强1贾丽萍2刘明海2(1中国矿业大学资源与安全学院, 北京100083;2石家庄经济学院, 河北石家庄050031摘要:提出了一种快速Delaunay 四面体网格生成的分治算法三角剖分, 然后从边界三角

2、形开始递归生成四面体网格. 四面体, 边界三角形都将成为内部四面体的面, , , 而且.关键词:; ; 网格生成; 边界一致中图分类号:文献标识码:A 文章编号:167124512(2005 0520067203Algorithm of mesh generation of Delaunaytetrahedral in constrained domainGuan Wen ge W u Qiang Jia L i pi ng L i u M i nghaiAbstract :A fast dividing and conquering algorithm of mesh generation

3、of Delaunay tetrahedra in con 2strained domain was presented. After the boundary of the given domain was divided into Delaunay triangles from a given sample point set , tetrahedral meshes was made by selecting a suit point from given point cloud to have Delaunay property. The surface triangles becom

4、e a direct consequence of interior tetrahedron. The algorithm does not require any surface conforming checks to avoid penetrated surface boundaries and over 2lapped tetrahedrons.K ey w ords :constrained domain ; Delaunay tetrahedron ; mesh generation ; boundary conformanceG uan Wen ge Doctoral candi

5、date ; School of Resource &Security , China University of Mining and Tech 2nology , Beijing 100083, China.目前, 四面体网格生成技术广泛用于地质领域、有限元分析及科学计算可视化等许多领域. 该类算法研究已取得许多重要进展, 主要有局部变换法、Delaunay 三角化增量算法、四叉树/八叉树算法、单元生长法、网格前沿法等. 这些算法1,2对二维任意域Delaunay 三角化和三维凸域De 2launay 四面体网格生成的实现比较成熟, 但在三维限定区域四面体网格生成算法理论上还没有定

6、论. 本文提出一种高效的Delaunay 四面体网格生成算法, 生成的四面体符合空球原则, 边界三角形是内部四面体的面, 不需边界一致性检查, 可避免四面体穿过边界和狭长四面体的产生, 且算法方便编程.1基本概念四面体网格生成:将空间中任意分布的离散点用直线连在一起, 形成在空间既不重叠又无间隙的紧邻四面体集. Delaunay 规则:在三维情况下网格符合空球原则, 即每个四面体外接球面内不包含其他点. 在二维情况下网格符合空圆原则, 即每个三角形的外接圆弧内不包含其他点3. 约束Delaunay 剖分:在实际应用中, 部分散乱点之间常常存在某种约束关系, 在进行Delaunay 剖分时, 剖

7、分结果既要符合Delaunay 规则, 又要保证离散点之间的约束关系.第33卷第5期华中科技大学学报(自然科学版 Vol. 33No. 52005年5月J. Huazhong Univ. of Sci. &Tech. (Nature Science Edition May 2005 2约束数据域四面体生成算法描述2. 1数据结构考虑到算法需要, 本文仅对三角形和四面体数据结构简单描述. 三角形三个顶点的PID 按逆时针排布; 四面体的四个顶点A B CD 则按照如下原则排列:A B C 是四面体的一个面, 从其相对的顶点D 看按逆时针排列, D 是按右手法则大拇指指向的该面对应顶点.2

8、. 2数据域剖分主算法tri -Extend已知三维约束数据域S p s v |v R 3组成, . 首先将数据域S Delaunay 三角剖分, 构成三角形集合a fl =f |f =v i |i =1, 2, 3,v i R 3.本文算法采用文献4, 5中提供的算法实现.算法tri -Extend 是一个递归的分治算法, 输入数据是简单数据域S i 的散乱点坐标集合p s 和边界Delaunay 三角形集合, 输出结果是剖分所得四面体集合. 为提高四面体剖分效率, 算法的核心思想是在递归过程中依次利用垂直于X 轴、Y 轴或者Z 轴的平面a fl -face 将给定数据域空间分为两个半空间p

9、 s 1和p s 2, 并将初始边界三角形和四面体生成过程中产生的新三角形, 按照与平面a fl -face 交叉、在平面a fl -face 左侧、右侧的位置不同依次存放在三角形链表a fl0, a fl1和a fl2中. 函数每次递归仅对与平面a fl -face 相交的三角形构造Delaunay 四面体, 然后递归调用函数, 依次在半空间p 1和p 2中构造四面体.算法中a fl 是主程序传来的简单数据域边界三角形邻接面表, a fl -face 可以是垂直于X 轴、Y 轴或者Z 轴的平面, 初始值为垂直X 轴的平面. 边界三角形顶点顺序从给定数据域内部观看按照逆时针方向排列, 由四面体

10、产生的三角形从四面体外部观看按逆时针排布. 算法步骤如下:步骤1置三角形链表a fl0, a fl1, a fl2初值为空.步骤2使用一个平面a fl -face 将点集p s 分割为三部分p s 0, p s 1和p s 2.步骤3将a fl 分割为三个三角形邻接面表a fl0, a fl1和a fl2.步骤4对a fl0中每个三角形f 1调用簇状四面体集构造算法MakeSimplex (f , p s , a fl0 生成四面体簇t -link ,ten -lin k =ten -lin k t -link ; 对生成的新四面体除f 1之外的其他面f 调用函数Update (f , a f

11、l0 ; 求a fl -face 的新值.步骤5如果平面a fl1非空, 则递归调用函数tri -Extend (int s -id , p -bound , p s 1, a fl1 . 步骤6如果平面a fl2非空, 则递归调用函数tri -Extend (int s -id , p -bound , p s 2, a fl2 . 2. 3簇状四面体集算法MakeSimplex判定点P 与给定A B 所在平面的位置关系公式. D (x , y , , 已A 1, y , (, y 2, z 2 , C (x 3,3A B C 所在平面的位置:假定A B C 顶点顺序按逆时针方向排列, 按照

12、右手法则确定, 点D 在右手大拇指所指向方向, 表示在平面的上方; 在右手大拇指所指向方向逆向, 表示在平面下方.设H =xy z 1x 1y 1z 11x 2y 2z 21x 3y 3z 31, 则当H =0时, 表示D 点在平面上; 当H >0时, 表示D 点在平面上方; 当H <0时, 表示D 点在S 面的下方.MakeSimplex 函数算法思想是:依据给定的三角形f 构造一个与当前已经构造的四面体集合和已知的当前数据域边界三角形没有重叠和交叉的Delaunay 四面体. 构造过程实际上是寻找与三角形f 对应的另一个顶点V 的过程, 记T (f ,V 为三角形f 与顶点V

13、构成的四面体. 寻找的顶点V 符合以下原则:a . 顶点V 在三角形f 所在平面上方; b . 四面体T (f , V 是Delaunay 四面体, 即外接球内不包含其他顶点;c . 在寻找顶点V 的过程中同时记录下与V顶点在同一个外接球球面上的其他顶点集合, 存放在链表v -link 中. 四面体T (f , V 构造成功后, 会生成除f 之外的新三角形, 继续为这些新产生的三角形在链表v -link 中寻找合适的点构造新四面体.3算法正确性分析3. 1几个引理设点P 和点Q 位于三角形f 所在平面异侧, 以三角形f 所在平面为界将空间分为两个半空间, 用H (f , P 表示P 一侧的半空

14、间. T (f , P 表示三角形f 与P 构成的四面体.86华中科技大学学报(自然科学版 第33卷 引理1设点P 和点Q 在三角形f 的异侧, 则四面体T (f , Q 的外接球包含Q 点的充要条件是四面体T (f , P 的外接球包含P 点.引理2设点P 和点Q 在三角形f 的同侧, 如果点P 落在四面体T (f , Q 内部或者T (f ,Q 的外接球内部, 则四面体T (f , P 的外接球在半空间H (f , P 一侧的半球一定落在T (f ,Q 的外接球在半空间H (f , P 一侧的半球内.引理3若四面体T (f , P 是三角形f 与点P 构成的Delaunay 四面体, 则在

15、半空间H (f ,P 一侧T (f , P 以与三角形f 构成Delaunay H (f , P 与f 然在T (f , P .引理4若四面体T (f , P 是三角形f 与点P 构成的Delaunay 四面体, 设点P 与点Q 位于f 异侧, 若点Q 在T (f , P 外接球的球面上, 则T (f , Q 必然是Delaunay 四面体. 3. 2算法思想正确性分析a . 点V 在三角形f 所在平面上方是选取顶点的必要条件. 证明如下:三角形f 顶点按逆时针排列, f 或者是当前区域外边界三角形, 或者是一个已经构成的四面体的某个侧面, 若顶点V 在平面的下方, 四面体T (f , V 必

16、定在给定数据域边界以外或者与已经构成的四面体相交; 若在平面的上面, 则不能构成四面体. 因此顶点V 在三角形f 所在平面上方.b . 在有限时间内可生成四面体网格. 设p s是空间中的有限点集, 三角形f 已经与点D 构成Delaunay 四面体T (f , D , 则依据f 构造的另一个四面体的新顶点V 肯定在以f 为界的与点D 异侧的半空间内. 由于四面体T (f , D 是Delau 2nay 的, 因此其外接球内必然不包含V 点, 由引理1可知T (f , V 的外接球也必不包含D 点. 由引理2可知半空间H (f , D 一侧的点必然不包含在T (f , V 的外接球内, 因此在寻

17、找新点V 构造第二个四面体时, 只需测试半空间H (f , V 一侧的点是否包含在T (f , P 的外接球内即可.设点V 和V 1在四面体T (f , D 的同侧, 若四面体T (f , V 外接球内部包含点V 1, 由引理2可知, 为T (f , V 1 的外, 选择V 1作为一, , 因此最V , 使得四面体T (f , V 是Delaunay 四面体.由引理3和引理4可知, 应该将位于当前正在测试四面体外接球球面上的点用链表保存, 然后选择形状最佳的四面体.参考文献1Sloa S W. A fast algorithm generating constrained De 2launay

18、 triangulation of Boundary SurfacesJ.Comput 2ers &Structures , 1991, 39(5 :4935002Xu Y A , Y ang Q. The algorithm of 3D constrained De 2launay triangulation J .Journal of S oftware , 2001, 12(1 :1031103Victor J D T. Delaunay triangulations in creation :Anoverview and a linear 2time algorithm J.International Journal of G eographical Information Systems , 1993, 7(6 :5015244Wu Q , Xu H. An ap