Delaunay三角剖分算法介绍
一.什么是Delaunay三角剖分
從事數值計算相關領域的讀者,相信或多或少都聽說過“三角剖分”這個概念。在諸如有限元仿真,光線追蹤渲染等計算當中,都需要把幾何模型轉化為三角網格數據,即“三角網格生成”。在這個過程中,三角剖分是十分重要的一步。下面我們就具體了解一下三角剖分的概念。
三角剖分:
顧名思義,三角剖分(triangulation)就是對給定的平面點集,生成三角形集合的過程。考慮平面點集P={?P1...Pn?},我們希望得到三角形集合T=在{?t1...tm?},滿足:
a)所有三角形的端點恰好構成集合P。
b)任意兩個三角形的邊不相交(要么重合,要么沒有交點)。
c)所有三角形的合集構成P的凸包(convex hull)。
注:有時我們也稱三角剖分為“三角化”,二者表示同一概念。而后者有時也指從多邊形或者多面體生成三角網格的過程,即“三角網格生成”。為避免混淆,本文中不使用“三角化”這一詞語,而專用“三角剖分”和“三角網格生成”來指代兩個不同的過程。
一般來說給定一個點集,往往存在不止一個三角剖分。以下圖為例,我們只要把第一個三角剖分中的兩個相鄰三角形的公共邊做一次翻轉,就可以得到一個新的三角剖分。
圖一:對于給定點集的兩種三角剖分
由于給定點集的三角剖分不唯一,我們希望從中挑選出一個“最優”的三角剖分。那么何為最優呢?這涉及到對三角形“質量”(quality)的評定。一般在數值計算中,我們不希望三角形過于狹長,也就是說越接近等邊三角形越好。以下是幾種常見的質量評定標準:
a) 最小角(minimum angle):即所有三角形的內角當中最小的角。
b) 縱橫比(aspect ratio):三角形最短邊與最長邊的比例。
c) 半徑比(radius ratio):三角形內接圓半徑的兩倍與外接圓半徑的比例。
Delaunay三角剖分
所有三角形的外接圓均滿足空圓性質的三角剖分,稱為一個Delaunay三角剖分。
空圓性質
即一個三角形(或邊)的外接圓范圍內(邊界除外),不包含點集P中的任何頂點。
對于上面的例子來說,左圖的三角剖分是Delaunay的,任何一個三角形的外接圓內部都不包含點集中的頂點;而右圖的三角剖分不是Delaunay的,因為下面的兩個三角形的外接圓內部都包含了頂點。
圖二:左圖的剖分滿足空圓性質,而右圖的剖分不滿足
二.Delaunay三角剖分的相關理論
這部分簡單介紹一下Delaunay三角剖分的相關理論,是為了更好的理解Delaunay三角剖分的最優性質和構造算法。對理論不感興趣的小伙伴可以直接跳到結論部分。先介紹一下Delaunay三角形和(局部)Delaunay邊的概念。
a) 對三角剖分中的一個三角形,如果其外接圓滿足空圓性質,稱為一個Delaunay三角形。
b) 對三角剖分中的一條邊,若存在一個外接圓滿足空圓性質,稱為一條Delaunay邊。(一條邊的外接圓可以是經過這條邊兩個頂點的任意圓,不唯一。)
c)對三角剖分中的一條屬于兩個三角形?t1,t2?的邊e,若存在一個外接圓不包含三角形?t1,t2中的任何頂點稱為局部Delaunay邊。若一條邊只屬于一個三角形,也屬于局部Delaunay邊。
由定義可知,Delaunay邊一定是局部Delaunay的,但反過來未必。以下圖為例,左側三角剖分中的藍色邊是局部Delaunay邊,因為可以找到一個外接圓不包含A,B,C,D中任何一個點。但它并非一條Delaunay邊,因為它的外接圓總是包含其他兩個頂點。而右側三角剖分中的藍色邊不是局部Delaunay邊,因為它的任意一個外接圓都至少包含A,B,C,D中的一個點。
3.左圖的邊BD是局部Delaunay邊,而右圖的邊AC不是局部Delaunay邊。
讀者們可以再去回看一下圖2,其中左側三角剖分中的藍色邊既是Delaunay邊也是局部Delaunay邊,而右側三角剖分中的藍色邊既非Delaunay也非局部Delaunay。可見判斷“局部”二字意味著它只與相鄰的兩個三角形有關。
以下結論被稱為Delaunay引理
對一個三角剖分T,以下三個命題相互等價
a) T中所有三角形均為Delaunay三角形。(Delaunay三角剖分的定義)
b) T中所有三角形的邊均為Delaunay邊。
c) T中所有三角形的邊均為局部Delaunay邊。
此處略去證明,感興趣的讀者可以閱讀Jonathan Schewchuk關于Delaunay triangulation的講義。
基于局部Delaunay邊的概念,可以定義下面的翻轉邊算法(flip algorithm)。
定理:對三角剖分T中的一條邊e,若它不是局部Delaunay的,則可以被翻轉成為一條局部Delaunay邊e’。
邊的翻轉:假設一條邊e屬于兩個三角形?t1,t2?,這兩個三角形的合集所構成的四邊形有兩條對角線,e為其中一條。現在用另一條對角線e’替換e,得到兩個新的三角形?t1′,t2′?,并用這兩個新的三角形替換原三角剖分中的?t1,t2?,得到新的三角剖分T’。這個操作稱為邊的翻轉。
圖2和圖3從右到左都是一次翻轉邊操作。可以看到,每進行一次翻轉邊操作時,都把三角剖分“局部”地改善了。更具體的講,有以下結論:
a) 最小角增大:翻轉邊e后的兩個相鄰三角形的6個內角中的最小角大于等于原先的兩個三角形6個內角中的最小角。
b) 外接圓半徑縮小:翻轉邊e后的兩個相鄰三角形的外接圓半徑最大值小于等于原先的兩個三角形外接圓半徑的最大值。
c) 四點共圓的情況:若一條邊e的兩個相鄰三角形的頂點A,B,C,D四點共圓,則無論翻轉與否,這條邊都是局部Delaunay邊。
圖4.經過一次翻轉邊操作后,兩個相鄰三角形的最小角增大了,而外接圓半徑縮小了,這說明這次翻轉邊操作改善了三角剖分的質量。
根據前面介紹的Delaunay引理,對任意一個三角剖分T,只要持續進行上述的翻轉邊操作,最終總可以轉化為一個Delaunay三角剖分。這實際上提供了一種構造Delaunay三角剖分的思路(Lawson算法,詳見第三部分)。而根據上面提到的關于翻轉邊的結論,可以得到下面三條關于Delaunay三角剖分的最優性質。
Delaunay三角剖分的最優性質
a) 最大化最小角性質:
在給定點集P的所有三角剖分當中,Delaunay三角剖分得到的最小角(所有三角形的內角中的最小值)是最大的。
b) 最小化外接圓性質:
在給定點集P的所有三角剖分當中,Delaunay三角剖分得到的最大外接圓半徑(所有三角形的外接圓半徑中的最大值)是最小的。
c)若點集P中任意四點不共圓,則存在唯一的Delaunay三角剖分T。若點集P中四點A,B,C,D共圓,且△ABC,△BCD屬于Delaunay三角剖分T,那么將邊BC翻轉后得到的三角剖分T’(包含△ABD,△ACD)同樣是一個Delaunay三角剖分。
其中最大化最小角性質正是我們在第一部分中提到的,三角形質量的度量標準之一。由于具備這些最優性質,在許多三角剖分的應用場景中,Delaunay三角剖分總是會被優先考慮。
三.Delaunay三角剖分的構造算法
大部分Delaunay三角剖分的構造算法都是通過逐點插入來實現的,這里我們介紹兩種: Lawson算法和Bowyer-Watson算法。
前面提到的翻轉邊算法已經可以幫助我們消除“壞邊”(指非局部Delaunay的邊)。而每當插入一個新的頂點時,將其與附近的舊頂點連接會加入一些新的邊,我們只要用翻轉邊算法來修正其中的“壞邊”,就可以將其重新“修復”為Delaunay三角剖分。順著這個思路我們就會得到下面的Lawson算法。
Lawson算法:
a) 先計算點集P的包圍盒(bounding box),將包圍盒的四個頂點加入P中得到P’。根據包圍盒生成兩個超三角形(super triangles),構成初始三角剖分?T0?。由于只包含兩個直角三角形,?T0?是(包圍盒四個頂點的)一個Delaunay三角剖分。
b) 將點集P中的頂點逐一插入現有的三角剖分?Ti?中,并進行如下調整:
?設插入的頂點v位于三角形t中,將v與三角形的三個頂點連接,使t分裂為3個三角形?t1,t2,t3。
?分別檢查?t1,t2,t3?是否滿足空圓性質,若不滿足則進行翻轉邊操作,直到沒有壞邊為止。此時得到一個包含頂點v的新Delaunay三角剖分。
c) 當最后一個頂點插入到三角剖分中,并且完成所有翻轉邊操作后,我們得到了點集P’的一個Delaunay三角剖分。現在刪除第一步中加入的包圍盒的四個頂點,并且去除所有與它們連接的三角形,則剩下的三角形就構成點集P的Delaunay三角剖分T。
在Lawson算法的基礎上,可以再做一點改進:在插入新的頂點之后,不要立刻與舊頂點連接生成新的三角形,再去逐一翻轉壞邊。而是先去尋找哪些邊是需要翻轉的,直接將這些邊刪掉形成一個“空穴”(cavity),然后在空穴當中連接新的邊。這就是Bowyer-Watson算法的思路。下面來詳細描述一下調整的過程。
Bowyer-Watson算法
a) 同Lawson算法
b) 將點集P中的頂點逐一插入現有的三角剖分中,并進行如下調整:
?在現有三角剖分中,所有外接圓包含頂點v的三角形的合集構成一個“星形多邊形”(star shaped polygon)。星形多邊形的含義是多邊形的任何一個頂點到v的連線都在多邊形內部。
?對于上述星形多邊形,將其內部的三角形全部刪除,形成一個“空穴”。將空穴邊界的頂點與新添加的頂點v連接得到新的三角形,替代剖分中被刪除的三角形,此時得到一個包含頂點v的新Delaunay三角剖分。
c) 同Lawson算法
當然,上述算法的正確性是有一些理論來證明的,這里就不詳細探討了,感興趣的讀者可以查看相關講義。下圖展示了Bowyer-Watson算法的調整過程。
圖5.左圖:原Delaunay三角剖分及新插入的頂點v,中間:所有外接圓包含頂點v的三角形構成一個星形多邊形,右圖:將星形多邊形內部的三角形刪除,并將v與空穴頂點重連后得到的新三角剖分,這是一個包含頂點v的Delaunay三角剖分。
由于Delaunay三角剖分的唯一性,兩種算法在效果上沒什么區別。但Bowyer-Watson算法因為縮短了搜索壞邊的過程,效率更高一些,所以在實際應用中我們往往會使用Bowyer-Watson算法。
四.基于Delaunay三角剖分的網格生成算法
三角網格生成(mesh generation)算法有很多種,這里我們介紹基于Delaunay三角剖分的算法(其他還有波前法,八叉樹方法等等)。
圖6.左圖的幾何模型由點、線、面等幾何元素組成,右圖為基于Delaunay三角剖分生成三角網格,全部由三角形組成。
在實際應用中,三角網格生成的對象大多是三維空間中的幾何曲面。常見的幾何曲面,如平面,圓柱面,球面,B樣條曲面,都是在二維空間參數化。因此實際上我們是在二維的參數空間上生成三角網格,最后再根據參數表達式投射到三維空間中。這里我們就簡述一下對于給定的二維參數區域,如何基于Delaunay三角剖分去生成一個三角網格。
生成三角網格與構造三角剖分的主要區別在于,點集P并非給定,需要自己去生成。而一旦確定了點集P,剩下的工作就是構造三角剖分了。實際上,生成點集P與構造三角剖分的過程是同時進行的。其過程如下:
Delaunay平面區域網格生成算法
a) 計算區域的包圍盒,生成兩個超三角形作為初始三角網格。根據邊界曲線生成邊界點。
b) 將邊界點逐一插入到三角網格中。每一步插入后,用Bowyer-Watson算法得到新的Delaunay三角網格。
c) 將第一步插入的輔助頂點刪除(同時刪除與其連接的三角形),得到關于全部邊界點的Delaunay三角網格(亦即邊界點的Delaunay三角剖分)。
d) 在區域內部的三角形,根據一定的條件,將重心(centroid)插入三角網格,并根據Bowyer-Watson算法調整,得到新的Delaunay三角剖分。直到不再有需要插入的重心點,此時得到區域的完整Delaunay三角網格。
下面用一個簡單的例子來直觀地展示這個過程:
圖7.從上至下,從左至右依次為:1)參數區域 2)由包圍盒生成的初始三角網格及邊界曲線上生成的邊界點3)將邊界點逐一插入并調整后得到的Delaunay三角網格 4)去除輔助頂點及相連的三角形得到的關于邊界點的Delaunay三角網格5)加入內部點之后得到的完整Delaunay三角網格。
每個曲面都可以用上述算法在參數區域上生成其Delaunay三角網格,然后投射到三維空間中。然而實際的模型往往包含多個曲面,因此要保證邊界的一致性。下面我們介紹完整算法:
三維模型的Delaunay三角網格生成算法
a) 遍歷模型的所有邊(edge),生成邊界點。
b) 遍歷模型的所有面(face),收集其邊界點在參數區域的坐標。然后在參數區域上按照上面的算法生成Delaunay三角網格。
c) 將每個面在參數區域上生成的內部點和三角形投射到三維空間,添加到三角網格中。
這里要注意的是,每個面上的點需要兩份編號,一個局部(local)編號用于參數區域上的三角網格,及一個全局(global)編號用于三維空間的整體三角網格。在網格生成過程中需要時刻保持這兩種編號的對應關系。
五.結果展示
以下結果為齒輪和斯坦福兔子的幾何模型以及shonMesh中生成的Delaunay三角網格。目前相關功能還在開發中,歡迎對網格生成算法感興趣的讀者與作者聯系。感謝閱讀!
總結
以上是生活随笔為你收集整理的Delaunay三角剖分算法介绍的全部內容,希望文章能夠幫你解決所遇到的問題。
 
                            
                        - 上一篇: 高性能Mysql运维应用实战-高俊峰-专
- 下一篇: MSP430F2111IPWR 超低功耗
