用于保特征四邊形網格生成的改進Morse算法

籍 冉 冉， 鄭 曉 朋， 雷 娜， 羅 鐘 鉉

( 大連理工大學 軟件學院, 遼寧 大連 116620 )

0 引 言

數字幾何用于描述空間中的幾何物體，在工業設計、影視、游戲等領域都有廣泛的應用.數字幾何首先需要在計算機中對空間幾何模型進行描述，從而在此基礎上進一步進行建模、分析、仿真等.三角形網格是表述中最常用的一種網格形式，生成技術也比較成熟，但在許多應用領域，四邊形網格卻具有無可替代的優勢.在數值模擬仿真中，四邊形網格的穩定性、精度、收斂速度等都優于三角形網格；四邊形網格也可以用于生成樣條曲面；四邊形網格在模型參數化、紋理貼圖等領域也有著廣泛的應用.在四邊形網格生成中保持模型的特征至關重要.

曲面四邊形網格生成的方法很多，主要有基于方向場的方法[1-2]、參數化方法[3-4]、“分而治之”的方法[5-6]和基于Morse理論的方法.Dong等[7]基于Morse理論，分析拉普拉斯矩陣的特征向量得到四邊形網格，算法首先計算輸入網格的拉普拉斯矩陣，求解拉普拉斯矩陣的特征向量，選擇合適的特征向量作為特征函數；然后利用該函數生成Morse-Smale復形完成對三角網格的四邊單元剖分；最后在此剖分的基礎上進行全局分片參數化，并在參數域上進行規則采樣獲得最終的四邊形網格.在之后的工作中，Huang等[8]、Zhang等[9]和Ling等[10]加入了方向、特征對齊和各向異性的約束，通過優化求解帶有約束的廣義特征問題來求解特征函數.Huang等[8]將譜分析與離散余弦變換的關系應用在特征函數的優化上，建立了特征函數與其對應Morse復形的內在聯系，通過改變拉普拉斯算子中的面積項來控制單元密度；通過在廣義拉普拉斯特征問題中引入額外的能量項來控制Morse復形邊的方向與位置，方向由用戶輸入的方向場控制，特征對齊對應于網格邊的具體位置.Zhang等[9]提出一種基于駐波構造的各向異性四邊形網格生成技術，通過挖掘四邊形網格和駐波的內在關系，將所有的用戶需求同時轉化為駐波的各種物理屬性約束，通過輸入密度場實現了各向異性和密度控制，同樣也通過輸入方向場來實現方向和特征對齊約束.Ling等[10]通過引入邊界條件來實現特征對齊控制，通過求解Helmh?ltz方程來得到特征函數，在Helmh?ltz方程中加入了密度場的控制.這些算法均需將所有用戶控制和質量需求集中在一個全局優化問題中統一求解，算法復雜，隨著模型規模的增大，計算代價將急速上升，噪聲增多，也需要后期復雜的優化.2018年，Fang等[11]將參數化和Morse理論結合起來，算法通過用戶輸入標架場來進行參數化，在參數化的過程中也輸入特征線來進行特征約束，在參數化發生退化的區域用Morse理論來生成四邊形網格，通過將兩種方法優勢互補，解決了參數化的退化問題和Morse方法在規模較大時運算量大的缺陷.

本文將曲率信息加入到拉普拉斯矩陣中，通過迭代的方法計算特征函數，在迭代的過程中將特征線的信息加入，使得臨界點能夠準確定位到特征線上，并通過臨界點來構造保特征的Morse-Smale復形，進而參數化生成四邊形網格.

1 技術路線

基于Morse理論生成滿足特征的四邊形網格，算法步驟概括為5步，圖1是原始算法與本文算法的對比，本文接下來將對每一步進行更加詳細的介紹.

(a) 原始算法步驟

(b) 本文算法步驟

2 Morse理論與Morse-Smale復形

2.1 Morse理論

(a) 極小值點

(b) 鞍點

(c) 極大值點

Banchoff[13]用分段線性函數替代光滑函數f將Morse理論推廣到三角網格上，f的值定義在三角網格頂點上，通過判斷頂點的一環領域內的點來判斷臨界點.三角網格上光滑的Morse函數通常是求解一個拉普拉斯方程得到，即

(1)

式中：Ni為頂點i的一環領域的所有頂點的集合，wij表示與頂點i相連接的邊(i，j)的權重，在三維流形曲面網格上，通過用拉普拉斯矩陣的形式來表示拉普拉斯算子，即

(2)

式中：αij、βij為邊(i，j)的兩個三角形的兩個對角.-Lf=λf，特征函數即為計算特征值所對應的特征向量，特征值0=λ1≤λ2≤…≤λn所對應的特征向量e1，e2，…，en，隨著特征值的增大，頻率越來越高，臨界點的個數越來越多，在模型上的可視化如圖3所示.拉普拉斯方程的解稱為調和函數，因此拉普拉斯矩陣的特征向量也是離散情形下球面調和函數的基.

2.2 Morse-Smale 復形

Morse函數f，如果它的穩定流形和不穩定流形之間只有橫向相交，則f是Morse-Smale函數[12]，穩定流形與不穩定流形相交形成的胞腔稱為Morse-Smale復形.Morse函數f定義在三角

網格上，f的值定義在三角網格頂點上，通過比較頂點的值與一環領域內點的值來判斷臨界點的類型，如圖4所示.

Morse-Smale復形是在一個流形上標量函數的胞腔剖分.如圖5所示，從鞍點出發沿梯度上升或者下降最快的方向，直到達到極大值點或者極小值點，一般來說會從鞍點引出兩條上升梯度線和兩條下降梯度線，這些路徑把流形分成若干四邊形區域，所有的四邊形區域都有兩個鞍點、一個極大值點、一個極小值點.

Morse-Smale復形的生成直接關乎到最終生成的四邊形網格的質量，而生成Morse-Smale復形的關鍵在于求解Morse函數，在計算Morse函數時通常會充分考慮特征對齊約束的要求，但由于離散化中的精度損失以及噪聲點的影響，臨界點不能完全準確地定位在特征線上，導致Morse-Smale復形中的邊沒能與特征對齊，因此需要對臨界點進行局部擾動，使其準確定位到特征線上.本文針對這一點提出了相應的改進算法，下面將詳細說明.

3 加入特征約束的Morse-Smale復形

生成滿足模型特征的網格在曲面四邊形化中是非常重要的一個約束要求.而生成滿足模型特征的Morse-Smale復形，關鍵在于求解拉普拉斯矩陣的特征函數，即Morse函數.本文給出如何求解特征函數，使得臨界點準確地定位在特征線上.

3.1 加入曲率信息的拉普拉斯矩陣

在計算三角網格上標量函數時，通過拉普拉斯矩陣的某個特征值的特征向量來定義，即

-Lf=λf

(3)

通過借鑒之前的研究[8，14]，對拉普拉斯矩陣的對角線進行擾動.在拉普拉斯方程中加入了模型曲率的信息，使用三角網格曲面模型頂點上的高斯曲率，并歸一化：

(4)

其中θi為點v的一環領域的三角形中所對應的頂點為v的角的大小.曲率反映了模型表面的彎曲程度，刻畫了模型的幾何特征信息，求解特征函數的方程定義如下：

-Lf=λ(I+K)f?-(I+K)-1Lf=λf

(5)

其中矩陣K為對角矩陣，對角線上的元素ki等于第i個頂點vi處的高斯曲率.

拉普拉斯矩陣在加入曲率信息后，臨界點會因為加入的曲率值而進行相應的位置擾動，臨界點能夠更加精確地定位到模型特征的地方，生成的Morse-Smale復形更加貼合模型的特征，如圖6所示.

3.2 特征函數的迭代算法

通過計算拉普拉斯矩陣L的某個特征值λk的特征向量來作為特征函數f，在這里一般使用的是ARPACK稀疏矩陣特征問題求解器來求解，求解出來的只能是特定的數值λk所對應的特征函數，特征函數的選取不具有靈活性，具體體現在以下兩個方面：

(a) 加入曲率信息之前的結果

(b) 加入曲率信息之后的結果

(1)2.1節分析的臨界點個數隨著特征值的數值增大而增多，如果在求解中λk所對應的特征函數fk的臨界點個數偏少，而λk+1所對應的特征函數fk+1的臨界點個數偏多，根據原始的方法不能求得fk和fk+1之間的特征函數；

(2)想要對特征值λk所對應的特征函數fk做少許的擾動，或者是對局部頂點上的特征函數做擾動，使得臨界點移動到模型特征的地方，根據原始的方法也沒有辦法實現.

從以上兩點可以看出，特征函數的選取受到局限，因此本文提出了一種迭代算法來計算特征函數.

(6)

對于任意的頂點vi，式(6)的第i行可以寫為

(7)

這里λ不僅限于特征值，可以是任意數值，則相對應求解出來的f就更具有靈活性，則迭代公式即為式(7).根據拉普拉斯算子的離散形式(1)，定義全局M上的調和能量[15]

(8)

當迭代求解趨于穩定時，調和能量逐漸減小并趨于0.想要對特征值λk所對應的特征函數fk做相應擾動，就將fk作為迭代初值.算法流程如下：

輸入： 拉普拉斯矩陣L.

輸出： 特征函數f.

Step 1計算矩陣L的特征值λk所對應的特征函數fk，將fk作為迭代的初值f0.

Step 2計算調和能量E0.

Step 3對所有網格上的點進行更新：

Step 4計算調和能量Em+1，當|Em+1-Em|<ε時，停止迭代，否則返回Step3.

3.3 加入特征信息

在模型表面生成四邊形網格時，希望生成的網格的邊能與模型的特征線對齊，在基于Morse理論生成四邊形網格的過程中，Morse-Smale復形的邊一定是四邊形網格的邊，所以生成的Morse-Smale復形的邊與模型的特征線對齊.在本文算法中，關鍵是讓Morse-Smale復形的臨界點準確定位到模型的特征線上，步驟如下：

(1)首先在模型的表面畫出特征線，在本文中，為了驗證算法的效果，只選取了部分特征線，特征線的集合為Se，定義這些邊的特征屬性作為輸入，如圖7所示

圖7 加入特征線信息

(2)在3.2節的迭代算法流程中，針對Step 3，如果[vi，vj]∈Se，μij>1，否則μij=1.針對不同的模型μij的取值不同，在本文的模型中μij∈[1，2].特征函數的迭代公式為

(9)

即增大特征線的權重.通過增大特征線的權重，就可以使得原本在特征線附近的臨界點，準確定位到特征線上，并且采用迭代的算法，對不是特征線上的點沒有產生太大的影響，如圖8給出了增加特征線權重前后的結果對比.

圖8 加入特征線的臨界點前后對比圖

從圖中可以看出，臨界點都準確地定位到了特征線上，并且其他位置的特征點沒有發生太大的變化.從之前的分析也可以知道，Morse-Smale復形的邊都是從鞍點出發，連向極大值點或極小值點，因此，只要臨界點在特征線上，Morse-Smale復形的邊大多就會沿著特征線，少數情況下在特征線上存在兩個相鄰的臨界點為極值點，兩個點之間沒有復形的邊，此時可以通過對偶操作或者補齊操作來優化，補齊操作的優化過程在下一節會詳細介紹，優化使得生成的復形的邊可以沿著模型的特征線，如圖9所示，藍線是補齊算法補齊的特征線.

圖9 沿特征線的復形結果圖

4 Morse-Smale復形的優化

4.1 噪聲點的消除

初步生成的Morse-Smale復形會產生很多的噪聲點，如圖10所示.

通過持久性簡化[16]來消除不穩定的點對.在生成的Morse-Smale復形中，對于兩個相鄰的臨界點vi、vj，它們對應于特征函數所對應的標量場的一組拓撲特征，持久性值定義為兩個點的函數值之差的絕對值，即|f(vi)-f(vj)|，其是描述移除這個拓撲特征時，特征函數的變化量，持久性值較小的點對屬于不穩定點對，要進行消除.一對相鄰的臨界點，必然一個是鞍點，另一個是極值點，定義在它們上的消除操作是指將鞍點和極值點消除，同時將所有連向該極值點的梯度線連向鞍點的另一個同屬性的極值點.如圖11所示.

(a) 存在大量噪聲的結果圖

(b) 消除噪聲的結果圖

圖11 持久性簡化操作

4.2 幾何優化

希望生成的Morse-Smale復形是近似于正方形或者矩形，但是由于根據特征函數計算梯度線過程中的不可控制性，初始復形的形狀可能會出現菱形、梯形等，如果不對其優化，會對后續參數化生成四邊形網格的質量產生影響，因此對全局的復形進行形狀的優化.

(1)文獻[14]中提到通過參數化和松弛迭代的方法來重新定位復形的邊緣點，文獻[8]中由于加入了特征線，為了不使復形的邊緣偏離特征線，對文獻[14]的方法進行了擴充，固定特征線上的參數化值，借鑒文獻[8]的算法，可以對復形的邊緣進行優化；

(2)本文迭代算法中，需要對特征線上相鄰的兩個極值點之間進行補齊操作，由此一個復形區域被分割成了兩個三角形，針對因為這種情況生成的三角形區域，可以對這些區域采用“分而治之”的思想，將1個三角形剖分成3個四邊形.

5 四邊形網格生成

優化完成后的Morse-Smale復形可以保證生成高質量的四邊形網格.將每一個復形區域參數化到平面上，這一步可以通過調和映射完成[15]，根據密度要求d，在二維平面上生成規則的d×d維的四邊形網格，再將其映回到模型表面，輸出四邊形網格.

本文給出了幾種模型迭代結果，如圖12(c)，臨界點可以很準確地與特征線對齊，通過補齊操作可以實現最終的保特征四邊形網格生成；圖12(d)給出了補齊操作的結果.幾何優化部分和最后的分區域參數化四邊形網格生成都是非常成熟的算法，在本文中沒有顯示，可以在以后的優化中實現.

圖12 4個模型的結果

6 結 語

本文在Morse理論的基礎上，提出了一種特征函數的新計算方法.通過迭代算法計算特征函數，使得特征函數的選取更加具有靈活性；在迭代算法中加入特征約束，可以將臨界點準確定位在特征線上，生成Morse-Smale復形；最后對初始形成的Morse-Smale復形進行拓撲優化和幾何優化，進而生成保特征的四邊形網格.本文所展示的結果與文獻[8-11]的結果相比，存在不足，但是本文提出了計算特征函數的一種特殊的新算法，并且這個算法實現了臨界點對齊這一好的結果，對后續的保特征的四邊形網格生成起到了關鍵作用，在優化完成后也可以生成質量較好的四邊形網格，并且本文提出的算法簡單，易于實現，輸入信息較少，計算效率高.

后續的研究工作包括以下內容：(1)實現對目前四邊形網格的優化，并將本文算法應用到實際問題中，驗證其有效性；(2)將本文的算法推廣到高維，生成三維的Morse-Smale復形，進而生成六面體網格.