有限元分析中的二维Delaunay三角网格剖分代码实现
有限元分析中的二维Delaunay三角网格剖分
有限元分析中的二维Delaunay三角网格剖分摘要本文从有限元分析出发，引出三角网格剖分的概念。
随后着重介绍了二维平面点集的Delaunay三角剖分。
给出了一些重要的Delaunay三角形的定理和性质，也体现出了Delaunay三角剖分的优点。
接着重点分析了构造二维Delaunay三角形的空洞算法，并用程序完成了它。
最后又分析了算法中的不足，并给出论文改进的方法。
关键词：Delaunay三角形，V oronoi图，网格剖分III1 第一章绪论1.1网格剖分的背景有限元分析是数学的一个分支。
其思想是将复杂的问题简单化，然后进行处理。
处理办法是将整个研究对象分成一些有限的单元，然后对每个小单元做相应的处理，最后整合起来去逼近原来的整个对象。
所以我们可以看到，有限元分析中将单元剖分的越小，得到的近似值就会越逼近真实值。
但是往往我们需要处理的对象很复杂，需要的计算量也很大，人工很难完成。
在早起年代，这个问题也阻止了有限元分析的发展。
近年来，随着计算机的发展，带动了一些需要大量计算的科学领域的发展。
有限元分析就是其中一种，因为当计算机取代人力之后，其快速的计算能力作用愈发凸显，人们只需要控制相应的算法即可。
作为最常用的处理手段，被大大的发展了之后，有限元分析也被应用于诸多方面。
早期的有限元分析主要应用与航空航天和地质、地球物理方面，现在越来越多的在工程分析计算和计算流体力学中看见。
图 1.1图 1.2常见的有限元分析可以分为六大步骤：问题及求解域的定义、求解域的网格剖分、确定状态变量及控制方法、单元推导、总装求解和结果解释。
上述步骤又可被分为三大阶段：前置处理、计算求解和后置处理。
而在前置处理中网格剖分作为最重要又最复杂的一个步骤，其处理结果制约着有限元的最后逼近结果。
网格剖分有很多形式：二维的主要剖分形状有三角形、四边形，三维的有四面体、六面体。
在有限元分析中网格剖分有如下要求：1、节点合法性。
指每个单元的节点最多只能是其他单元的节点或边界点，而不能是内点。
Delaunay三角剖分及matlab实例
Delaunay三⾓剖分及matlab实例鉴于Delaunay三⾓剖分在点云拟合⽅⾯有不错的应⽤，现对该算法的原理进⾏简单的汇总~----------------------------原理部分------------------------1、三⾓剖分与Delaunay剖分的定义如何把⼀个离散⼏何剖分成不均匀的三⾓形⽹格，这就是离散点的三⾓剖分问题，散点集的三⾓剖分，对数值分析以及图形学来说，都是极为重要的⼀项处理技术。
该问题图⽰如下：1.1 三⾓剖分定义【定义】三⾓剖分：假设V是⼆维实数域上的有限点集，边e是由点集中的点作为端点构成的封闭线段，E为e的集合。
那么该点集V的⼀个三⾓剖分T=(V,E)是⼀个平⾯图G,该平⾯图满⾜条件：1、除了端点，平⾯图中的边不包含点集中的任何点。
2、没有相交边。
//边和边没有交叉点3、平⾯图中所有的⾯都是三⾓⾯，且所有三⾓⾯的合集是散点集V的凸包。
//：⽤不严谨的话来讲，给定⼆维平⾯上的点集，凸包就是将最外层的点连接起来构成的凸多边型，它能包含点集中所有的点。
1.2 Delaunay三⾓剖分的定义在实际中运⽤的最多的三⾓剖分是Delaunay三⾓剖分，它是⼀种特殊的三⾓剖分。
先从Delaunay边说起：【定义】Delaunay边：假设E中的⼀条边e(两个端点为a,b)，e若满⾜下列条件，则称之为Delaunay边：存在⼀个圆经过a,b亮点，圆内（注意是园内，圆上最多三点共圆）不含点集V中任何其他的点，这⼀特性⼜称空圆特性。
【定义】Delaunay三⾓剖分：如果点集V的⼀个三⾓剖分T只包含Delaunay边，那么该三⾓剖分称为Delaunay三⾓剖分。
1.3 Delaunay三⾓剖分的准则要满⾜Delaunay三⾓剖分的定义，必须符合两个重要的准则：1、空圆特性：Delaunay三⾓⽹是唯⼀的（任意四点不能共圆），在Delaunay三⾓形⽹中任⼀三⾓形的外接圆范围内不会有其它点存在。
Delaunay三角剖分
public Tin(Point_T p1, Point_T p2, Point_T p3) { pthree[0] = p1; pthree[1] = p2; pthree[2] = p3; lthree[0] = new Line_T(p1, p2); lthree[1] = new Line_T(p2, p3); lthree[2] = new Line_T(p3, p1); } }
envpts.Add(pts[0]); envpts.Add(pts[pts.Count - 1]); othpts.Remove(pts[0]); othpts.Remove(pts[pts.Count-1]); pts.Sort(comxandy); if(!envpts.Contains(pts[0])) { envpts.Insert(1, pts[0]); } if (!envpts.Contains(pts[pts.Count - 1])) { envpts.Add(pts[pts.Count - 1]); } othpts.Remove(pts[0]); othpts.Remove(pts[pts.Count-1]); //构建以x-y，x+y最大最小值组成的初始矩形框 int i = 0; int tag = 0; bool over = true; while(i<envpts.Count) { Line_T cline; if (i==envpts.Count-1) { cline = new Line_T(envpts[i], envpts[0]); } else { cline = new Line_T(envpts[i], envpts[i + 1]); } double dismax=0; for (int j = 0; j < othpts.Count ;j++ ) { if (IsLeftPoint(othpts[j], cline)) { double distance = PointToLine(othpts[j], cline); if (distance > dismax) { dismax = distance; tag = j;
三角剖分
for(j=O；j‘pNT；j++)
{dc．MoveTo(．)；
de．LineTo()；
dc．LineTo()
dc．LineTo()； }
6.Delaunav剖分是二维Delaunav三角形化和三维De—launay四面体化的总称.
在Ma￡lab中．剖分的命令是tri=delaunay(x，y)，给出显示的数据结构是每个三角形编号和每个三角形顶点的编号，画剖分三角形的命令是triplot(tri，x，y)，用tsearch(x，y，tri，xi，yi)来找点(xi，yi)所在的三角形，返回的点所在三角形的编号，要找到某个三角形顶点的坐标用的则是x(tri(l，2))，指的是上边剖分的数据结构tri中第一个数据结构中三角形第二个顶点的x的坐标。
3.两种以上语言混合编程是不好的。如果预先把Mathb／Pdetool工具箱库函数对求解域Delaunay三角化过程独立完成，把三角化的结果文件的形式用于单一的高级程序语言开发的有限元前处理程序中，既避免了有限元前处理程序中两种语言同时运行。又充分利用了Matlab／பைடு நூலகம்detool工具箱中优良的Delaunay剖分算法，缩短了开发周期，提高了前处理程序的可移植性和适用性
4.pdetool工具箱已为我们提供了几何形状描述、网格剖分、网格加密、三角形
单元局部调整和图形显示等一系列相关的库函数标准格式，且相关源代码都是公开的，这就为基于Matlab／pdetool的求解域三角化提供了必要的条件。
5.获得了特定求解域三角化剖分结果的数据文件后。就要把它用到有限元前处理程序中，此时，前处理程序的编制就非常简单，只需要先把3个文件复制到有限元程序的相应目录中，然后在程序实体中首先读这些数据文件。以Visual C++程序为例，运用VC中库函数ifstream可以把数据文件读进来，但要注意，它默认的是按行读的，而数据文件中，每一列为一个操作单元，所以在获得数据文件时．就先把数据进行转置再存储。程序中读完数据后通过简单的连线语句就可以实现剖
基于Delaunay三角剖分处理二维欧式空间MTSP的近似算法
基于Delaunay三角剖分处理二维欧式空间MTSP的近似算法寿涛;刘朝晖【摘要】This paper discussed Multi Travelling Salesman Problem (MTSP)on 2D Euclidean space.This problem could be simplified to solve several TSP by Delaunay Triangulation.It could be proven that the approximate ratio of Tree Decomposed Algorithm was 2 and the core proof was based on empty circle property of Delaunay edge.The paper testified the performance and efficiency of the algorithm by some numerical examples.%考虑了在二维欧式平面内的多旅行商问题,通过Delaunay三角剖分的方法,将问题转化为求解多个旅行商问题.树分解算法的核心是Delaunay边的空圆性质并且可以证明该算法的近似比为2.最后,通过数值模拟验证了算法的有效性.【期刊名称】《华东理工大学学报（自然科学版）》【年(卷),期】2017(043)006【总页数】4页(P895-898)【关键词】MTSP;Delaunay三角剖分;近似算法【作者】寿涛;刘朝晖【作者单位】华东理工大学数学系,上海200237;华东理工大学数学系,上海200237【正文语种】中文【中图分类】O224多旅行商问题(MTSP)是TSP问题的推广[1]。
通常可以把MTSP问题拆分成2个子问题,即:先确定每个旅行商访问客户点集;再对每个旅行商访问的点求解TSP问题。
对于MTSP问题的研究,最早可以追溯到1990年Li等[2]的5近似算法。
Delaunay三角剖分的最优化网格节点生成算法研究
Delaunay三角剖分的最优化网格节点生成算法研究张晶飞;李射;崔向阳【摘要】针对任意域Delaunay三角剖分存在的局部网格质量不佳问题,提出了一种改进的Delaunay算法.利用边界三角形单元节点和重心的关系是否满足右手定则来判断初始三角形单元是否位于剖分域内的三角形重心法,并保留剖分域内的三角形单元;对待插入节点进行最优化处理以获得高质量网格,避免产生畸形单元;算例结果表明,所提方法可以适应复杂几何边界区域的划分,并可获得质量较高的三角形网格.【期刊名称】《电子设计工程》【年(卷),期】2019(027)009【总页数】7页(P10-16)【关键词】Delaunay三角剖分;任意域;有限元网格;节点【作者】张晶飞;李射;崔向阳【作者单位】湖南大学汽车车身先进设计制造国家重点实验室,湖南长沙410082;湖南大学汽车车身先进设计制造国家重点实验室,湖南长沙410082;湖南大学汽车车身先进设计制造国家重点实验室,湖南长沙410082【正文语种】中文【中图分类】TP391对任意平面而言，Delaunay三角化法[1-3]（DT）、栅格法和推进波前法是目前最为流行的3种平面网格划分方法[4]。
其中Delaunay三角化是计算几何的重要研究内容之一[5]，其拥有优异的特性及较完善的理论体系，在很多领域得到了广泛应用。
该方法也是最为流行的全自动网格生成方法之一，它的“最大-最小角”特性可以自动避免生成小内角长薄单元[6]，因此特别适用于有限元网格剖分。
但是在实际应用中仍然存在一些问题亟待解决，比如产生域外三角形、确定新节点插入位置等问题。
许多学者对产生域外三角形的问题进行了研究并提出了解决方案，如L。
A.Piegl[7]的边界外法向法、凸分解法[8]、边界指示法[9]、矢量三角形法[10]和解决域法[11]。
凸分解法需要将凹域分解成多个凸域，但难以处理多联通域问题。
矢量三角形法可以解决简单的凹多边形，然而对于自相交多边形效果不佳。
基于Delaunay三角网格剖分算法在三维造型中的研究
基于Delaunay三角网格剖分算法在三维造型中的研究作者：王牌来源：《科学与财富》2014年第06期摘要：在对三维图像进行有限元数值模拟解析时，为了对连续的计算区域进行数值计算，达到模拟仿真的效果，必须先对三维图像进行网格剖分。
Delaunay三角网格剖分算法是生成网格的一种有效方法。
本文介绍了Delaunay三角网格剖分算法，以及在约束条件下的网格细分，最后给出了该算法在三维实体造型中的应用。
关键词：三角剖分；网格生成；网格细分Abstract： In the simulation analysis of the 3D finite element numerical， in order to carry out the numerical calculation for the calculation of continuous area， achieve the simulation results， we must first on the 3D mesh. Delaunay triangulation algorithm is an effective method to generate mesh. This paper introduces the Delaunay triangulation algorithm， and in the condition of mesh subdivision， finally the application of the algorithm in 3D solid modeling are given in this paper.Keywords： triangulation，mesh generation，mesh subdivision1、引言网格生成是有限元模拟计算的先决条件，有限元计算的效率和精确度在很大程度上受生成的网格质量的影响。
Delaunay三角剖分的快速实现
第25卷第3期2005年5月海 洋 测 绘HY DROGRAPH I C S URVEYI N G AND CHARTI N GVol 125,No 13M ay,2005收稿日期:2005203215作者简介:汪连贺(19732),男,天津人,工程师,主要从事GI S 软件开发研究。
D e l aunay 三角剖分的快速实现汪连贺,董 江(天津海事局海测大队,天津 300220) 摘要:主要探讨了以Delaunay 三角剖分的逐点插入法为基础构建不规则三角网的方法,并在程序设计中对该算法进行了改进,极大地提高了Delaunay 三角网的构建效率。
关键词:不规则三角网;Delaunay 三角网;优化策略中图分类号:P208 文献标识码:B 文章编号:167123044(2005)03200512031 引 言不规则三角网的自动联结(又称三角剖分),无论是在地理信息系统中,还是在计算机辅助设计系统中,都占有十分重要的位置。
目前,三角剖分以Delaunay 法的应用最为广泛。
Delaunay 三角剖分的方法主要有如下三种:三角网生长法、逐点插入法和分割—归并法。
本文在成图系统中采用了Delaunay 三角剖分的逐点插入法,并对该算法进行了改进,极大地提高了三角网的构建效率。
2 D e launa y 三角网的性质(1)Delaunay 三角网是唯一的;(2)三角网的外边界构成了点集P 的凸多边形“外壳”;(3)没有任何点在三角形的外接圆内部,反之,如果一个三角网满足此条件,那么它就是Delaunay 三角网;(4)如果将三角网中的每个三角形的最小角进行升序排列,则Delaunay 三角网的排列得到的数值最大,从这个意义上讲,Delaunay 三角网是“最接近于规则化的”的三角网。
3 逐点插入法的基本步骤(1)定义一个包含所有数据点的初始多边形;(2)在初始多边形中建立初始三角网,然后迭代以下步骤,直至所有数据点被处理;(3)插入一个数据点P,在三角网中找出包含P 点的三角形,把P 点与三角形的三个顶点相连,生成新的三角形;(4)用LOP 算法优化三角网。
Python实现DelaunayTriangulation
Python实现DelaunayTriangulationDelaunay三角剖分是一种广泛用于计算机图形学、地理信息系统和计算几何等领域的技术。
它将一组离散点集划分成一组互不重叠的三角形，同时满足一定的性质。
在Python中实现Delaunay三角剖分可以使用几何库scipy和matplotlib。
Scipy库中的Delaunay函数可以实现Delaunay三角剖分，而matplotlib库可以用于可视化结果。
下面是一个简单的例子，展示了如何使用Python实现Delaunay三角剖分：```pythonimport numpy as npimport matplotlib.pyplot as pltfrom scipy.spatial import Delaunay#生成随机点points = np.random.rand(30, 2)# 执行Delaunay三角剖分tri = Delaunay(points)# 绘制Delaunay三角剖分结果plt.triplot(points[:,0], points[:,1], tri.simplices)plt.plot(points[:,0], points[:,1], 'o')plt.show```在这个例子中，我们首先生成了30个随机点，然后使用Delaunay函数执行Delaunay三角剖分，并传入这些点作为参数。
最后，使用matplotlib库的triplot函数绘制Delaunay三角剖分结果，然后再绘制原始点。
这段代码可以生成一个简单的Delaunay三角剖分示意图，并且可以通过调整点的数量或位置来实现不同的效果。
实际应用中，Delaunay三角剖分常用于三维地理信息系统中的高程插值、人脸识别中的特征点定位、计算几何中的最近邻等领域。
此外，还可以使用其他几何库，如CGAL和OpenCV等库，来实现更复杂的Delaunay三角剖分算法。
delaunay 三角剖分 步骤
delaunay 三角剖分步骤1. Delaunay三角剖分是用于将点集分割成不规则三角形的方法。
The Delaunay triangulation is a method for dividing a set of points into irregular triangles.2.首先选择一个点作为起始点。
First, select a point as the starting point.3.然后选择另外两个点与起始点构成一个三角形。
Then select two other points to form a triangle with the starting point.4.接着选择一个未被包含在任何三角形内的点。
Then select a point that is not included in any triangle.5.在所有的三角形中寻找能将这个新点包含进去的三角形。
Find a triangle among all the triangles that can include this new point.6.如果找到了这样的三角形，将这个三角形和新点围成的区域删除。
If such a triangle is found, remove the area enclosed by this triangle and the new point.7.在新的边缘上寻找新的三角形。
Find new triangles on the new edges.8.重复以上步骤，直到所有的点都被包含在三角形内。
Repeat the above steps until all points are included in triangles.9. Delaunay三角剖分具有无重叠、最小化夹角和最大化最小角的性质。
Delaunay triangulation has the properties of non-overlapping, minimizing angles, and maximizing minimum angles.10.可以使用Delaunay三角剖分来进行网格生成和空间分析。
delaunay-三角剖分算法
一、概述Delaunay 三角剖分算法是计算机图形学领域中常用的一种算法，它可以将给定的点集进行高效的三角剖分，用于构建网格、进行地理信息系统分析、建立三维模型等应用。
本文将对该算法的原理、实现和应用进行介绍。
二、算法原理1. 待剖分点集在进行Delaunay三角剖分之前，需要准备一个点集，这个点集是待剖分的对象。
点集的数量取决于具体的应用，可以是二维平面上的点，也可以是三维空间中的点。
2. Delaunay 三角形在进行三角剖分时，Delaunay 三角形是一种特殊的三角形，满足以下性质：a. 任意一个点要么位于Delaunay 三角形的外接圆内部，要么位于外接圆的边上；b. 任意两个Delaunay 三角形之间的外接圆不相交。
3. Delaunay 三角剖分Delaunay 三角剖分是将给定点集进行三角剖分的过程，它的目标是构建满足Delaunay 三角形性质的三角形集合。
三、算法实现1. 基于增量法的实现增量法是Delaunay 三角剖分的一种经典算法，它的基本思想是逐步增加点，并根据Delaunay 三角形的性质进行调整。
具体步骤如下： a. 初始化：选择一个超级三角形包含所有点集，作为初始三角剖分；b. 顺序插入点：逐个将待剖分点插入到当前三角剖分中，并进行调整；c. 边界检测：检测新增的边界是否需要进行修正；d. 优化处理：对新增点周围的三角形进行优化调整。
2. 时间复杂度分析增量法的时间复杂度主要取决于点集的数量和点的分布情况，一般情况下，其时间复杂度可以达到O(nlogn)。
四、算法应用1. 图形渲染在计算机图形学中，Delaunay三角剖分常用于构建网格、进行三维渲染等。
它可以有效地分割空间，使得渲染效果更加真实。
2. 地理信息系统地理信息系统中常常需要对地理数据进行空间分析，Delaunay三角剖分可以帮助构建地理网格，进行地形分析、资源评估等。
3. 三维建模在三维建模领域，Delaunay三角剖分可以用于构建复杂的三维模型，并支持模型的分析、编辑等功能。
用c++实现Delaunay三角剖分
Delaunay三角剖分void main(PointLink *ptUpOutline, PointLink *ptDownOutline){//首先分别对上下轮廓线的点集进行三角化，ptUpOutline, ptDownOutline为对应点集TrianglizationInPlan(ptUpOutline) ;TrianglizationInPlan(ptDownOutline) ;//根据两层三角化的结果进行Delaunay剖分，得到四面体链表GetTetrahedron( UpTriangle , DownTriangle ) ;//后处理用于删除不正确的四面体PostProcess() ;}其中：1．函数TrianglizationInPlan()用于平面三角化void TrianglizationInPlan( PointLink * Head ){//对不符合Delaunay剖分的边进行细分bool Over = 1;while( Over ){First = Head ;Second = First->next ;while( Second ){if( !NotIncludeOtherPoint( First->Point , Second->Point , Head ) ){PointLink* Point = new PointLink ;First->next = Point ;Point->next = Second ;Over = false ;}First = Second ;Second = Second->next ;}}//对平面进行三角化，保存结果到TriangleLink的链表Trianglization( Head ) ;}2．函数GetTetrahedron()进行Delaunay剖分得到四面体链表void GetTetrahedron( TriangleLink * First, TriangleLink * Second )Cur = First ;while( Cur ){Temp1 = Cur->Triangle ;Center = GetCenter( Temp1 ) ;Cur = ContourUp ;while( Cur ){//初始化minif( Cur == ContourUp ) {min = Center.distance( Cur->Point ) ;Value = min ;Result = Cur->Point ;}//找距离最小的点else {Value = Center.distance( Cur->Point ) ;if( min > Value ){min = Value ;Result = Cur->Point ;}}Cur = Cur->next ;}TetrahedronLink *T1 = new TetrahedronLink ;T1->volume.Point[0] = Cur->Triangle.Point[0] ;T1->volume.Point[1] = Cur->Triangle.Point[1] ;T1->volume.Point[2] = Cur->Triangle.Point[2] ;T1->volume.Point[3] = Result ;T1->type = 1 ;ResultCur->next = T1 ;T1->next = 0 ;ResultCur = T1 ;Cur = Cur->next ;}Cur = Second ;while( Cur ){Temp1 = Cur->Triangle ;Center = GetCenter( Temp1 ) ;//逐个比较找到离圆心最近的点Cur = ContourDown ;while( Cur ){//初始化minif( Cur == ContourDown ){min = Center.distance( Cur->Point ) ;Value = min ;Result = Cur->Point ;}else //找距离最小的点{Value = Center.distance( Cur->Point ) ;if( min > Value ){min = Value ;Result = Cur->Point ;}}Cur = Cur->next ;}TetrahedronLink *T2 = new TetrahedronLink ;T2->volume.Point[0] = Cur->Triangle.Point[0] ;T2->volume.Point[1] = Cur->Triangle.Point[1] ;T2->volume.Point[2] = Cur->Triangle.Point[2] ;T2->volume.Point[3] = Result ;T2->type = 2 ;ResultCur->next = T2 ;T2->next = 0 ;ResultCur = T2 ;Cur = Cur->next ;}EdgeCur1 = EdgeUp ;EdgeCur2 = EdgeDown ;Edge tempEdge1 , tempEdge2 ;while( EdgeCur1 ){tempEdge1 = EdgeCur1->line ;while( EdgeCur2 ){tempEdge2 = EdgeCur2->line ;//判断两条边投影是否相交if( EdgesIntersect( tempEdge1 , tempEdge2 ) ){TetrahedronLink *T12 = new TetrahedronLink ;T12->volume.Point[0] = tempEdge1.Start ;T12->volume.Point[1] = tempEdge1.End ;T12->volume.Point[2] = tempEdge2.Start ;T12->volume.Point[3] = tempEdge2.End ;T12->type = 12 ;ResultCur->next = T12 ;T12->next = 0 ;ResultCur = T12 ;}EdgeCur2 = EdgeCur2->next ;}EdgeCur1 = EdgeCur1->next ;}}3．函数EdgesIntersect()两条边投影相交的判断bool EdgesIntersect(Edge edge1, Edge edge2){C3DPoint Vector , Vector1 , Vector2 ;C3DPoint Middle = edge1.End + edge1.Start ;Middle.x /= 2 ;Middle.y /= 2 ;Middle.z = edge2.Start.z ;Vector = edge1.End - edge1.Start ;Vector1 = edge2.Start - Middle ;Vector2 = edge2.End - Middle ;if( ( Vector*Vector1 >= 0 && Vector*Vector2 >= 0 ) || ( Vector*Vector1 <= 0 && Vector*Vector2 <= 0 ))return true ;Middle = edge2.End + edge2.Start ;Middle.x /= 2 ;Middle.y /= 2 ;Middle.z = edge1.Start.z ;Vector = edge2.End - edge2.Start ;Vector1 = edge1.Start - Middle ;Vector2 = edge1.End - Middle ;if( ( Vector*Vector1 >= 0 && Vector*Vector2 >= 0 ) || ( Vector*Vector1 <= 0 && Vector*Vector2 <= 0 ))return true ;return false ;}4．函数PostProcess()用来删除外部四面体void PostProcess(){TetrahedronLink* Cur = Head ;TetrahedronLink* Pre , * p ;Triangle Temp ;while( Cur ){switch( Cur->type ) {case 1:Temp.Point[0] = Cur->volume.Point[0] ;Temp.Point[1] = Cur->volume.Point[1] ;Temp.Point[2] = Cur->volume.Point[2] ;//根据向量间的叉乘，判断是否在轮廓线外if( IsOutSide( ContourUp , Temp ) ){if( Cur == Head ){p = Head ;Head = Cur->next ;delete p ;}else{p = Cur ;Pre->next = Cur->next ;delete p ;}}break;case 2:Temp.Point[0] = Cur->volume.Point[0] ;Temp.Point[1] = Cur->volume.Point[1] ;Temp.Point[2] = Cur->volume.Point[2] ;if( IsOutSide( ContourDown , Temp ) ) {if( Cur == Head ){p = Head ;Head = Cur->next ;delete p ;}else{p = Cur ;Pre->next = Cur->next ;delete p ;}}break;default:break;}Pre = Cur ;Cur = Cur->next ;}}。
用VC语言实现任意多边形的Delaunay完全三角剖分算法
用VC语言实现任意多边形的Delaunay完全三角剖分算法①涂治红 桑农(华中科技大学图像识别与人工智能研究所,图像识别与智能控制国家教委重点实验室 武汉 430074)摘 要多边形三角剖分是计算几何的一个几何基元,它可以简化问题规模,在计算机图形学、模式识别等方面有重要的应用。
本文针对已有的Delaunay三角剖分算法的不足,提出新算法,并采用Visual C语言MFC类进行链表的管理,使得编程容易实现。
整个算法简洁通用。
最后给出了在实际中的应用。
关键词:任意多边形 Delaunay三角剖分 链表 MFC类中图分类号:TP301.6Delaunay T riangulation Algorithm of Arbitrary Polygons with Visual C LanguageTu Zhihong S ang N ong(Institute for Pattern Recognition and Artificial Intelligence,HUST,State Education Committee K ey Lab for Image Processing and Intelligent Control,Wuhan430074)Abstract:Triangulation of arbitary polygons is geometric primitives of computational geometry.It can predigest,dimensions. There are so many applications in graphics,pattern recognition and so on.This paper proposes an improved algorithm of Delaunay triangulation of the arbitrary polygon.The programming with Visual C language is relatively simple by using the MFC function to manage the lists.This algorithm is concise and general.The application of this algorithm is presented.K ey w ords:arbitrary polygon,Delaunay triangulation,list,MFCClass number:TP301.61 引言在计算机三维曲面造型,有限元计算和模式识别等领域里,经常要解决平面多边形的三角剖分问题。
2D-Delaunay 三角网格的数据结构与遍历
2D-Delaunay三角网格的数据结构与遍历高晓沨1，黄懿2（1．清华大学数学系，北京 100080；南开大学数学科学学院，天津 300071） 摘要：本文总结了二维Delaunay三角网格的Bowyer-Watson自动生成算法及其实现步骤，提出了一种类的结构、函数范例（采用Visual C++ 6.0编写程序），并讨论了遍历三角网格各种方法的优劣性，给出实验数据对比；最后得出结论，用广度优先的遍历方法创建网格是生成三角网格一种相对便利有效率的方法；另外，讨论了初始点加入顺序对程序运行时间的影响。
关键词：Delaunay三角网格 类结构 自动生成 广度优先遍历Data Structure and Traverse of 2D-Delaunay TriangulationGao Xiaofeng1, Huang Yi2(Department of Mathematics ,Tsinghua University; Beijing,100080College of Mathematics, Nankai University, Tianjin, 300071) Abstruct: The article summarized the realization of 2D-Delaunay triangulation, thesteps of creating Bowyer-Watson automatic mesh generator, and then constructed a kindof Class structure as well as functions of this algorithm (using Visual C++ 6.0), andfinally discussed the advantage and disadvantage of different methods to traverse thetriangle mesh, using data examination as contrast. Lastly, the author got theconclusion that Width First Traversal method is more effective and convenient. Besides,we discussed the effect between the order of original point set and running time ofthe program.Key Word: Delaunay Triangulation, Class Structure,automatic generation, Width FirstTraversal1．引言近年来，平面任意点集的三角网格化(triangulation)问题一直是人们密切关注的问题。
基于Delaunay三角剖分的二维交互建模研究
第43卷第2期物探化探计算技术Vol.43No.2 2021年3月COMPUTING TECHNIQUES FOR GEOPHYSIC A L AND GEOCHEMIC A L EXPLORATION Mar.2021文章编号：1001-1749(2021)02-025605基于Delaunay三角剖分的二维交互建模研究杨强（中国石化石油物探技术研究院，南京211103）摘要：准确合理的地质模型是地震正演模拟的基础。
简单的层状结构建模方法无法描述复杂的地质模型的拓扑结构，基于块体的地质建模方法可以描述断层、尖灭、透镜体等多种复杂地质结构，但建模方法实现比较难。
这里提出了一种解决方案，该方法基于Delaunay的地质建模方法，采用带约束的Delaunay三角剖分对复杂模型三角化，在此基础上，以原始线段为索引对剖分结果进行递归查找，形成多边形拓扑结构图和封闭多边形组，再对封闭多边形进行网格化。
通过理论和实际模型，验证了本方法的可行性。
关键词：地震正演；建模；Delaunay三角剖分中图分类号：TP391.41文献标志码：A 0引言地震勘探中，正演模拟不但在地震数据采集中得到应用，在地震资料处理和地震资料解释中也是重要的验证技术手段，是进行地震反演的基础。
而正演模拟的基础就是需要准确且合理的地质模型，因此建模是地震正演的重要基础。
传统上地质模型都沿用Cenveny提出的模型结构，即所谓层状结构模型。
它要求每一个分界面都必须从模型体的左边界贯穿到模型体的右边界，分界面按顺序由上到下依序排列，不得交叉。
对于逆断层、尖灭、透镜体等复杂地层情况，只能人为地简化模型,从而满足层状模型。
在实际的地球物理勘探中，层状结构模型有两个严重不足：①在采集、处理、解释各阶段的模型大都有逆断层、尖灭、透镜体等复杂地质元素，层状结DOI：10.3969/j.issn.1001-1749.2021.02.16构模型无法准确描述复杂的地质拓扑结构，从而无法得到网格化模型;②处理中的叠前深度偏移速度建模以及地震解释构造建模等，都要求能够交互修改模型，交互编辑中的层位修改及移动必须遵循层状规则，这对交互式复杂建模而言很困难，也不便利。
基于MATLAB实现二维delaunay三角剖分
基于MATLAB实现二维delaunay三角剖分34基于MATLAB 实现二维delaunay 三角剖分刘锋涛凡友华（哈尔滨工业大学深圳研究生院深圳 518055）【摘要】在已知凸多边形的顶点坐标的前提情况下，利用MATLAB 中的meshgrid 函数产生多边形附近矩形区域内的网格点的坐标，然后再利用inpolygon 函数判断哪些点位于多边形内和哪些点位于多边形的边界上。
在此基础上利用delaunay 函数来完成delaunay 三角剖分。
【关键词】delaunay 三角剖分；MATLAB网格划分是有限元分析前处理中的关键步骤，网格划分的密度以及质量对有限元计算的精度、效率以及收敛性有着重要的影响作用。
自20世纪70年代开始，关于有限元网格生成方法的研究已经取得了许多重要成果，提出许多有效的算法。
Ho-Le 对网格生成方法进行了系统的分类[1]。
许多学者也对网格生成的方法进行了综述，如我国的学者胡恩球等[2]、关振群等[3]。
Delaunay 三角剖分(简称DT)是目前最流行的通用的全自动网格生成方法之一。
DT 有两个重要特性：最大-最小角特性和空外接圆特性。
DT 的最大-最小角特性使它在二维情况下自动地避免了生成小内角的长薄单元。
因此特别适用于有限元网格生成。
大体上可将DT 算法分为三大类：分治算法，逐点插入法和三角网生长法。
经典DT 技术已经相当成熟，近年来的研究重点是约束DT 的边界恢复算法，以及如何克服算法退化现象所产生的薄元(sliver element)问题[3]。
然而，实现DT 有限元网格生成，对于非计算机图形学专业的工程师来说还是很复杂的。
在处理一些对有限元网格划分质量不过的问题时，如极限分析的有限元方法，可以采用一些更为简单的方法来实现。
在Matlab 计算软件中，已有一个成熟的函数delaunay 可以实现对一系列点的DT 划分。
因此，本文基于Matlab 的delaunay 等一些函数来完成一个凸多边形的DT 网格划分。
二维Delaunay三角剖分在地质建模中的稳定化研究
二维Delaunay三角剖分在地质建模中的稳定化研究张琳;孔祥铭;张尧尧【摘要】二维三角剖分是要将平面内的一些散点以及线段组成的集合剖分成不均匀的三角形网格,使得三角网格中只含有Denaulay边和Denaulay三角形.在地质建模中,由于数据规模和数据大小通常非常大,因此计算的浮点误差也会有很大的影响,在使用Triangle生成CDT时,往往因为这种误差,导致程序崩溃.文章将讨论这种误差的产生原因,进而提出改进算法,使Triangle更加稳定,满足实际需求.【期刊名称】《中国高新技术企业》【年(卷),期】2010(000)024【总页数】3页(P53-55)【关键词】三角剖分;Denaulay;vironoi图;Triangle;地质建模【作者】张琳;孔祥铭;张尧尧【作者单位】华北电力大学,北京102206;华北电力大学,北京102206;华北电力大学,北京102206【正文语种】中文【中图分类】TP391三角剖分问题是计算几何和计算机辅助设计与分析中的一个重要问题。
Delaunay 三角剖分是二维平面内的最优三角剖，它在有限元分析、信息可视化、计算机图形学等应用领域有着重要应用。
当前，对平面三角剖分及其约束条件下的三角剖分主要有分治、逐点插入三角网格插入法，大多约束三角剖分算法也是基于此基础之上。
在分析研究区域(二维) 的离散数据时，都可以尝试一下采用Delaunay 三角网或Voronoi 图的分析途径。
如GIS 中的网络分析。
描述地表形态的一种最佳方法，是地表(地貌和地物) 数字化表现的手段和分析工具。
当然，它的应用不仅适用于地学，而且活跃于所有与二维分析有关的领域，而且还将活跃于地图信息识别领域。
Triangle是用C语言编写的程序，可以生成二维的Delaunay三角剖分、Constraint Delaunay三角剖分和Voronoi图。
Triangle的计算速度快，存储效率高，并且稳定。
python中scipy库delaunay的simplices用法
python中scipy库delaunay的simplices用法全文共四篇示例，供读者参考第一篇示例：Scipy库是一个强大的Python科学计算工具箱，其中包含了许多实用的工具和算法。
其中的delaunay模块是用来进行Delaunay三角剖分的工具，它可以将给定的点集进行三角形网格化，从而用于空间插值、地理信息系统、有限元分析等领域。
在scipy库中，delaunay模块中的simplices属性是一个用于存储Delaunay三角形网格中每个三角形顶点索引的属性。
利用这个属性，我们可以很方便地获取Delaunay三角形网格中每个三角形的顶点索引，从而进行进一步的计算和分析。
我们需要导入scipy库和相关的模块：```pythonimport numpy as npfrom scipy.spatial import Delaunay```接着，我们可以生成一个随机的二维点集，并利用Delaunay函数进行Delaunay三角剖分：现在，我们可以通过simplices属性获取Delaunay三角形网格中每个三角形的顶点索引：输出结果类似如下：```[[0 1 3][1 2 3][1 4 5][1 3 4][0 3 6][2 7 3][1 5 8][4 3 7][8 6 0][7 3 5][6 3 8][3 7 5]]```每一行表示一个三角形，其中的三个数字分别代表该三角形的三个顶点在原始点集中的索引。
通过这些顶点索引，我们可以进一步计算三角形的边长、面积、周长等属性，或者对Delaunay三角形网格进行可视化展示。
除了simplices属性之外，Delaunay对象还提供了许多其他有用的属性和方法，比如`vertex_neighbor_vertices`方法可以用来获取每个点的相邻点的索引列表，`find_simplex`方法可以用来查找给定点所在的三角形索引等等。
这些属性和方法可以帮助我们更好地理解和利用Delaunay三角剖分的结果。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
有限元分析中的二维Delaunay三角网格剖分代码实现//二维平面点集的Delaunay三角剖分#include "stdafx.h"#include <iostream>#include <GL/glut.h>#include <ctime>#include <cmath>using namespace std;#define point_size 600#define pi 3.1415926struct point{float x,y;};struct triangle{point* Pv[3];float r_of_sqrt;point o_of_tr;};struct d_t_node{triangle Tnode;d_t_node*Pt_l[3];int position_in_dt;int zhuangtai_when_new;};point p1,p2,p3,p4;int n;point p[point_size];int dt_last=0;point p_in_dtriangle1[point_size+1];d_t_node Dtriangle[point_size];point p_in_dtriangle2[point_size+1];d_t_node *queue_t[point_size];point p_in_dtriangle3[point_size+1];int ps_last=0;int queue_t_last=0;point get_spoint_cin(point*p,int n);point get_spoint_rank(point*p,int n);void get_spoint_squre(point*p,int );void get_spoint_cir(point*,int);void read_spoint(point *p,int n);void Display();void init_D_triangle(d_t_node* D_tr);void Pjoin_ps(point*ps,point* p_new,int ps_last);void pjion_D_triangle(d_t_node*Dtriangle,int &dt_last ,point p_new);int findt_clude_p(d_t_node *Dtriangle,int dt_last,point p_new);void sort_to_line(point line_of_tr[point_size][2],int line_of_lin[point_size],int this_fh_when_new[point_size][2], int line_of_tr_last);float hls(float a,float b,float c,float d);point get_o(point p1,point p2,point p3 );float long_of_lines(point point1,point point2);bool p_in_q(d_t_node**queue,point p_new);void put_bug(d_t_node );void put_this(int psize[point_size][2],int last);void read_triangle(d_t_node Dtriangle[point_size] );void put_out_line(point line_of_tr[point_size][2],int line_of_tr_last);int main(int argc, char* argv[]){for(int ii=0;ii<point_size;ii++){for(int jj=0;jj<3;jj++){switch(jj){case 0:Dtriangle[ii].Tnode.Pv[jj]=&p_in_dtriangle1[ii+1];break;case 1:Dtriangle[ii].Tnode.Pv[jj]=&p_in_dtriangle2[ii+1];break;case 2:Dtriangle[ii].Tnode.Pv[jj]=&p_in_dtriangle3[ii+1];break;}}};//点集为圆上的点，输入圆的个数/* int m;cin>>m;n=m*36+1;get_spoint_cir(p,m); *///点集为正方形网格的交点，输入网格的行数/* int m;cin>>m;n=m*m;get_spoint_squre(p,m); *///点集为随机生成的点，输入点的数量cout<<" 请输入点数大小n=";cin>>n;get_spoint_rank(p,n);glutInit(&argc, argv);glutInitDisplayMode(GLUT_RGB | GLUT_SINGLE);glutInitWindowPosition(100, 100);glutInitWindowSize(300, 300);glutCreateWindow("第一个OpenGL程序");glColor3f(0.0,0.0,1.0);init_D_triangle(Dtriangle);for(int i=0;i<n;i++){pjion_D_triangle(Dtriangle,dt_last,p[i]);}glClearColor(1.0,1.0,1.0,1.0);glutDisplayFunc(Display);glutMainLoop();return 0;};void get_spoint_cir(point*p,int nn){float r;p[0].x=0.5;p[0].y=0.5;//int xi,xii;for(int iii=0;iii<nn;iii++ ){//xi=nn-1;xii=xi-iii;r=(float(iii)+1.0)*0.5/float(nn+1);for(int jjj=0;jjj<36;jjj++){p[iii*36+jjj+1].x=float(r*cos(float(jjj)*pi/18.0))+0.5;p[iii*36+jjj+1].y=float(r*sin(float(jjj)*pi/18.0))+0.5;}}};void get_spoint_squre(point*p,int nn){for(int i=0;i<nn;i++)for(int j=0;j<nn;j++)p[i*nn+j].x=float(j+1.0)/float(nn+1),p[i*nn+j].y=float(i+1.0)/float(nn+1);};void put_out_line(point line_of_tr[point_size][2],int line_of_tr_last){for(int i=0;i<line_of_tr_last;i++){cout<<"("<<line_of_tr[i][0].x<<" , "<<line_of_tr[i][0].y<<")-> ";cout<<"("<<line_of_tr[i][1].x<<" , "<<line_of_tr[i][1].y<<") ";}cout<<endl;};void read_triangle(d_t_node Dtriangle[point_size]){for(int i=0;i<dt_last;i++){glBegin(GL_LINE_LOOP);glVertex2f((Dtriangle[i].Tnode.Pv[0]->x-0.5)*2.0,(Dtriangle[i].Tnode.Pv[0]->y-0.5)*2);glVertex2f((Dtriangle[i].Tnode.Pv[1]->x-0.5)*2.0,(Dtriangle[i].Tnode.Pv[1]->y-0.5)*2);glVertex2f((Dtriangle[i].Tnode.Pv[2]->x-0.5)*2.0,(Dtriangle[i].Tnode.Pv[2]->y-0.5)*2);glEnd();};};void read_spoint(point *p,int n){for(int i=0;i<n;i++)glVertex2f(p[i].x,p[i].y);};void Display(){glClear(GL_COLOR_BUFFER_BIT);glPolygonMode(GL_FRONT_AND_BACK,GL_FILL);glPointSize(4.0);read_triangle(Dtriangle);glFlush();}；void put_this(int psize[point_size][2],int last){for(int i=0;i<last;i++)cout<<psize[i][0]<<" ";cout<<endl;for( i=0;i<last;i++)cout<<psize[i][1]<<" ";};void put_bug(d_t_node Dtriangle1){for (int i=0;i<3;i++){cout<<"point : "<<Dtriangle1.Tnode.Pv[i]->x<<" , "<<Dtriangle1.Tnode.Pv[i]->y<<endl;}cout<<"the_o: "<<Dtriangle1.Tnode.o_of_tr.x<<" , "<<Dtriangle1.Tnode.o_of_tr.y<<endl;cout<<"the_r "<<Dtriangle1.Tnode.r_of_sqrt <<endl;cout<<"poistion_in_Dt "<<Dtriangle1.position_in_dt<<endl;cout<<"zhuangtai: "<<Dtriangle1.zhuangtai_when_new<<endl;for (i=0;i<3;i++){if(Dtriangle1.Pt_l[i]!=NULL){cout<<" _ "<<Dtriangle1.Pt_l[i]->position_in_dt;}else cout<<" _NULL ";}cout<<endl;};bool p_in_q(d_t_node*queue_t,point p_new){float p_o,r_t;p_o=long_of_lines(p_new,queue_t->Tnode.o_of_tr);r_t=queue_t->Tnode.r_of_sqrt;if(p_o<=r_t) return 1;else return 0;};float long_of_lines(point point1,point point2){float long_r;long_r=sqrt((point1.x-point2.x)*(point1.x-point2.x)+(point2.y-point1.y)*(point2.y-point1.y)) ;return long_r;};float hls(float a,float b,float c,float d){float e;e=a*d-b*c;return e;};point get_o(point p1,point p2,point p3 ){point o;float a,b,c,d,e;a=(p2.x*p2.x-p3.x*p3.x+p2.y*p2.y-p3.y*p3.y)/2.0;b=p2.y-p3.y;c=(p3.x*p3.x-p1.x*p1.x+p3.y*p3.y-p1.y*p1.y)/2.0;d=p3.y-p1.y;e=hls(p2.x-p3.x,p2.y-p3.y,p3.x-p1.x,p3.y-p1.y);o.x=hls(a,b,c,d)/e;b=(p2.x*p2.x-p3.x*p3.x+p2.y*p2.y-p3.y*p3.y)/2.0;a=p2.x-p3.x;d=(p3.x*p3.x-p1.x*p1.x+p3.y*p3.y-p1.y*p1.y)/2.0;c=p3.x-p1.x;o.y=hls(a,b,c,d)/e;return o;};void sort_to_line(point line_of_tr[point_size][2],int line_of_lin[point_size],int this_fh_when_new[point_size][2], int line_of_tr_last){point change[2];int change1;int change_n1;for(int i=1;i<line_of_tr_last;i++){for(int j=i+1;j<line_of_tr_last;j++){if(line_of_tr[i-1][1].x==line_of_tr[j][0].x&&line_of_tr[i-1][1].y==line_of_tr[j][0].y ) {change[0].x=line_of_tr[j][0].x;change[0].y=line_of_tr[j][0].y;change[1].x=line_of_tr[j][1].x;change[1].y=line_of_tr[j][1].y;line_of_tr[j][0].x=line_of_tr[i][0].x;line_of_tr[j][1].x=line_of_tr[i][1].x;line_of_tr[j][0].y=line_of_tr[i][0].y;line_of_tr[j][1].y=line_of_tr[i][1].y;line_of_tr[i][0].x=change[0].x;line_of_tr[i][0].y=change[0].y;line_of_tr[i][1].x=change[1].x;line_of_tr[i][1].y=change[1].y;change1=this_fh_when_new[i][0];this_fh_when_new[i][0]=this_fh_when_new[j][0];this_fh_ when_new[j][0]=change1;change1=this_fh_when_new[i][1];this_fh_when_new[i][1]=this_fh_when_new[j][1];this_fh_ when_new[j][1]=change1;change_n1=line_of_lin[i];line_of_lin[i]=line_of_lin[j];line_of_lin[j]=change_n1;break;}}}};int findt_clude_p(d_t_node* Dtriangle,int dt_last,point p_new){for(int i=0;i<dt_last;i++){point o_of_t;float r_of_t,r_pando;float a1=Dtriangle[i].Tnode.Pv[0]->x,a2=Dtriangle[i].Tnode.Pv[0]->y;float b1=Dtriangle[i].Tnode.Pv[1]->x,b2=Dtriangle[i].Tnode.Pv[1]->y;float c1=Dtriangle[i].Tnode.Pv[2]->x,c2=Dtriangle[i].Tnode.Pv[2]->y;o_of_t=get_o(*Dtriangle[i].Tnode.Pv[0],*Dtriangle[i].Tnode.Pv[1],*Dtriangle[i].Tnode.Pv[2 ]);r_of_t=sqrt((o_of_t.x-a1)*(o_of_t.x-a1)+(o_of_t.y-a2)*(o_of_t.y-a2));r_pando=sqrt((o_of_t.x-p_new.x)*(o_of_t.x-p_new.x)+(o_of_t.y-p_new.y)*(o_of_t.y-p_new. y));if(r_pando<=r_of_t) return i;}return 0;};void pjion_D_triangle(d_t_node*Dtriangle,int &dt_last,point p_new){queue_t_last=0;d_t_node* t_clude_p=NULL;int here=findt_clude_p(Dtriangle,dt_last, p_new);t_clude_p=&Dtriangle[here];queue_t[queue_t_last]=t_clude_p,queue_t[queue_t_last]->zhuangtai_when_new=1;queue_t_last++;int queue_t_first=0;point line_of_tr[point_size][2];int line_of_tr_last=0;int line_of_lin1[point_size];int k;int this_fh_when_new[point_size][2];while(queue_t_first!=queue_t_last){for(int j=0;j<3;j++){if(queue_t[queue_t_first]->Pt_l[j]!=NULL&&p_in_q(queue_t[queue_t_first]->Pt_l[j],p_new)){if(queue_t[queue_t_first]->Pt_l[j]->zhuangtai_when_new!=1){queue_t[queue_t_last]=queue_t[queue_t_first]->Pt_l[j],queue_t[queue_t_last]->zhuangtai_when_new=1;queue_t_last++;}}else{if(j<2){line_of_tr[line_of_tr_last][0].x=queue_t[queue_t_first]->Tnode.Pv[j]->x;line_of_tr[line_of_tr_last][0].y=queue_t[queue_t_first]->Tnode.Pv[j]->y;line_of_tr[line_of_tr_last][1].x=queue_t[queue_t_first]->Tnode.Pv[j+1]->x;line_of_tr[line_of_tr_last][1].y=queue_t[queue_t_first]->Tnode.Pv[j+1]->y;if(queue_t[queue_t_first]->Pt_l[j]!=NULL)cout<<queue_t[queue_t_first]->Pt_l[j]->position_in_dt<<endl;if(queue_t[queue_t_first]->Pt_l[j]!=NULL){line_of_lin1[line_of_tr_last]=queue_t[queue_t_first]->Pt_l[j]->position_in_dt;for(int l=0;l<3;l++){if(queue_t[queue_t_first]->Pt_l[j]->Pt_l[l]!=NULL&&queue_t[queue_t_first]->Pt_l[j]->Pt_l[ l]->position_in_dt==queue_t[queue_t_first]->position_in_dt)k=l;}this_fh_when_new[line_of_tr_last][0]=queue_t[queue_t_first]->Pt_l[j]->position_in_dt;this_fh_when_new[line_of_tr_last][1]=k;}else{line_of_lin1[line_of_tr_last]=-1;this_fh_when_new[line_of_tr_last][0]=-1;this_fh_when_new[line_of_tr_last][1]=-1;}line_of_tr_last++;}else{line_of_tr[line_of_tr_last][0].x=queue_t[queue_t_first]->Tnode.Pv[j]->x;line_of_tr[line_of_tr_last][0].y=queue_t[queue_t_first]->Tnode.Pv[j]->y;line_of_tr[line_of_tr_last][1].x=queue_t[queue_t_first]->Tnode.Pv[0]->x;line_of_tr[line_of_tr_last][1].y=queue_t[queue_t_first]->Tnode.Pv[0]->y;if(queue_t[queue_t_first]->Pt_l[j]!=NULL){line_of_lin1[line_of_tr_last]=queue_t[queue_t_first]->Pt_l[j]->position_in_dt;for(int l=0;l<3;l++){if(queue_t[queue_t_first]->Pt_l[j]->Pt_l[l]!=NULL&&queue_t[queue_t_first]->Pt_l[j]->Pt_l[l]->p osition_in_dt==queue_t[queue_t_first]->position_in_dt)k=l;}this_fh_when_new[line_of_tr_last][0]=queue_t[queue_t_first]->Pt_l[j]->position_in_dt;this_fh_when_new[line_of_tr_last][1]=k;}else{line_of_lin1[line_of_tr_last]=-1;this_fh_when_new[line_of_tr_last][0]=-1;this_fh_when_new[line_of_tr_last][1]=-1;}line_of_tr_last++;}}}queue_t_first++;}sort_to_line(line_of_tr,line_of_lin1,this_fh_when_new,line_of_tr_last);int last_k,first_k,i;for( i=0;i<queue_t_last;i++){int k=queue_t[i]->position_in_dt;Dtriangle[k].Tnode.Pv[0]->x=p_new.x;Dtriangle[k].Tnode.Pv[0]->y=p_new.y;Dtriangle[k].Tnode.Pv[1]->x=line_of_tr[i][0].x,Dtriangle[k].Tnode.Pv[1]->y=line_of_tr[i][0] .y;Dtriangle[k].Tnode.Pv[2]->x=line_of_tr[i][1].x,Dtriangle[k].Tnode.Pv[2]->y=line_of_tr[i][1] .y;Dtriangle[k].Tnode.o_of_tr=get_o(*Dtriangle[k].Tnode.Pv[0],*Dtriangle[k].Tnode.Pv[1],*Dt riangle[k].Tnode.Pv[2]);Dtriangle[k].Tnode.r_of_sqrt=long_of_lines(Dtriangle[k].Tnode.o_of_tr,*Dtriangle[k].Tnode .Pv[0]);Dtriangle[k].position_in_dt=k;if(i>0){Dtriangle[last_k].Pt_l[2]=&Dtriangle[k];Dtriangle[k].Pt_l[0]=&Dtriangle[last_k];} if(line_of_lin1[i]!=-1){Dtriangle[k].Pt_l[1]=&Dtriangle[line_of_lin1[i]];Dtriangle[this_fh_when_new[i][0]].Pt_l[this_fh_when_new[i][1]]=&Dtriangle[k];}else Dtriangle[k].Pt_l[1]=NULL;Dtriangle[k].zhuangtai_when_new=0;if(i==0)first_k=k;last_k=k;}for(i=0;i<line_of_tr_last-queue_t_last;i++){Dtriangle[i+dt_last].Tnode.Pv[0]->x=p_new.x,Dtriangle[i+dt_last].Tnode.Pv[0]->y=p_new.y;Dtriangle[i+dt_last].Tnode.Pv[1]->x=line_of_tr[i+queue_t_last][0].x,Dtriangle[i+dt_last].Tn ode.Pv[1]->y=line_of_tr[i+queue_t_last][0].y;Dtriangle[i+dt_last].Tnode.Pv[2]->x=line_of_tr[i+queue_t_last][1].x,Dtriangle[i+dt_last].Tn ode.Pv[2]->y=line_of_tr[i+queue_t_last][1].y;Dtriangle[i+dt_last].Tnode.o_of_tr=get_o(*Dtriangle[i+dt_last].Tnode.Pv[0],*Dtriangle[i+dt _last].Tnode.Pv[1],*Dtriangle[i+dt_last].Tnode.Pv[2]);Dtriangle[i+dt_last].Tnode.r_of_sqrt=long_of_lines(Dtriangle[i+dt_last].Tnode.o_of_tr,*Dtri angle[i+dt_last].Tnode.Pv[0]);Dtriangle[i+dt_last].position_in_dt=i+dt_last;if(line_of_lin1[i+queue_t_last]!=-1){Dtriangle[i+dt_last].Pt_l[1]=&Dtriangle[line_of_lin1[i+queue_t_last]];Dtriangle[this_fh_when_new[i+queue_t_last][0]].Pt_l[this_fh_when_new[i+queue_t_last][1] ]=&Dtriangle[i+dt_last];}else Dtriangle[i+dt_last].Pt_l[1]=NULL;if(i>=0){Dtriangle[i+dt_last].Pt_l[0]=&Dtriangle[last_k];Dtriangle[last_k].Pt_l[2]=&Dtriangle[i+dt_last]; }Dtriangle[i+dt_last].zhuangtai_when_new=0;last_k=i+dt_last;}Dtriangle[last_k].Pt_l[2]=&Dtriangle[first_k];Dtriangle[first_k].Pt_l[0]=&Dtriangle[last_k];dt_last=dt_last+line_of_tr_last-queue_t_last;};void Pjoin_ps(point*ps,point* p_new,int ps_last){ps[ps_last]=*p_new;ps_last++;};void init_D_triangle(d_t_node* D_tr){d_t_node *q1=new d_t_node;d_t_node *q2=new d_t_node;D_tr[0].Tnode.Pv[0]->x=0.0; D_tr[0].Tnode.Pv[0]->y=0.0;D_tr[0].Tnode.Pv[1]->x=1.0; D_tr[0].Tnode.Pv[1]->y=0.0;D_tr[0].Tnode.Pv[2]->x=0.0; D_tr[0].Tnode.Pv[2]->y=1.0;D_tr[1].Tnode.Pv[0]->x=1.0; D_tr[1].Tnode.Pv[0]->y=0.0;D_tr[1].Tnode.Pv[1]->x=1.0; D_tr[1].Tnode.Pv[1]->y=1.0;D_tr[1].Tnode.Pv[2]->x=0.0; D_tr[1].Tnode.Pv[2]->y=1.0;D_tr[0].Tnode.o_of_tr=get_o(*D_tr[0].Tnode.Pv[0],*D_tr[0].Tnode.Pv[1],*D_tr[0].Tnode.P v[2]);D_tr[1].Tnode.o_of_tr=get_o(*D_tr[1].Tnode.Pv[0],*D_tr[1].Tnode.Pv[1],*D_tr[1].Tnode.P v[2]);D_tr[0].Tnode.r_of_sqrt=long_of_lines(D_tr[0].Tnode.o_of_tr,*D_tr[0].Tnode.Pv[0]);D_tr[1].Tnode.r_of_sqrt=long_of_lines(D_tr[1].Tnode.o_of_tr,*D_tr[1].Tnode.Pv[0]);D_tr[0].Pt_l[0]=NULL;D_tr[0].Pt_l[1]=&D_tr[1];D_tr[0].Pt_l[2]=NULL;D_tr[1].Pt_l[0]=NULL;D_tr[1].Pt_l[1]=NULL;D_tr[1].Pt_l[2]=&D_tr[0];dt_last=2;D_tr[0].position_in_dt=0;D_tr[1].position_in_dt=1;D_tr[0].zhuangtai_when_new=0, D_tr[1].zhuangtai_when_new=0;};point get_spoint_rank(point*p,int n){point back;back.x=0.0;back.y=0.0;srand(time(0));for(int i=0;i<n;i++){p[i].x=rand()/float(RAND_MAX);if(back.x*back.x<p[i].x*p[i].x) back.x=p[i].x;p[i].y=rand()/float(RAND_MAX);if(back.y*back.y<p[i].x*p[i].y) back.y=p[i].y;cout<<p[i].x<<" "<<p[i].y<<endl;}return back;};point get_spoint_cin(point*p,int n){float M=0,N=0;point back;for(int i=0;i<n;i++){cin>>p[i].x;if(M*M<p[i].x*p[i].x) M=p[i].x;cin>>p[i].y;if(N*N<p[i].y*p[i].y) N=p[i].y;}back.x=M;back.y=N;return back;};。