許 斌，李忠科，宋大虎

（第二炮兵工程大學(xué) 計(jì)算機(jī)教研室，陜西 西安710025）

0 引 言

使用掃描儀獲取模型表面的三維數(shù)據(jù)時(shí)，由于物體表面反射性、掃描設(shè)備的缺陷、模型表面的自遮擋等原因，點(diǎn)云數(shù)據(jù)通常含有孔洞，進(jìn)而導(dǎo)致重構(gòu)后的網(wǎng)格曲面也存在著孔洞?？锥吹拇嬖跇O大地影響了模型的完整性和顯示效果，在工業(yè)應(yīng)用場(chǎng)合則會(huì)導(dǎo)致后續(xù)操作無(wú)法進(jìn)行，因此，對(duì)這些孔洞進(jìn)行修復(fù)是曲面網(wǎng)格模型數(shù)據(jù)處理的一個(gè)關(guān)鍵步驟。

從實(shí)際應(yīng)用角度考慮，一個(gè)好的曲面修復(fù)算法應(yīng)盡可能地滿(mǎn)足下面兩個(gè)條件：①修復(fù)曲面盡量接近原始真實(shí)形狀，與周?chē)W(wǎng)格有很好的光順連接；②算法有較強(qiáng)的魯棒性，即能夠處理各種類(lèi)型孔洞?？傮w來(lái)說(shuō)，孔洞的修復(fù)問(wèn)題可以歸結(jié)為空間多邊形的三角剖分問(wèn)題，很多學(xué)者都對(duì)這個(gè)問(wèn)題進(jìn)行了深入研究，下面對(duì)其中一些具有代表性的算法進(jìn)行簡(jiǎn)要介紹。文獻(xiàn) ［1］提出了一種采用體數(shù)據(jù)場(chǎng)融合進(jìn)行孔洞修補(bǔ)的算法，該算法首先迭代計(jì)算擴(kuò)展體數(shù)據(jù)場(chǎng)的描述范圍，然后在此基礎(chǔ)上建立一個(gè)能夠描述整個(gè)孔洞區(qū)域及其周?chē)娴膱?chǎng)函數(shù)來(lái)完成孔洞缺失數(shù)據(jù)的修復(fù)。文獻(xiàn) ［2］使用RBF隱式曲面來(lái)對(duì)孔洞曲面進(jìn)行修復(fù)，算法中首先計(jì)算孔洞邊界頂點(diǎn)的鄰域頂點(diǎn)，然后根據(jù)其幾何屬性插值生成孔洞區(qū)域的RBF隱式曲面，最后對(duì)該曲面進(jìn)行網(wǎng)格化處理并與原始網(wǎng)格曲面進(jìn)行融合得到孔洞修復(fù)后的網(wǎng)格曲面。上述幾種算法的優(yōu)點(diǎn)是能夠修復(fù)較大的孔洞，但是對(duì)于缺失部分包含豐富幾何細(xì)的模型來(lái)說(shuō)上述幾種修復(fù)方法就存在局限性，因?yàn)樗鼈儾荒芎芎玫鼗謴?fù)模型缺失部分的幾何細(xì)節(jié)。

文獻(xiàn) ［3］提出了一種基于PDE （partial differential equations）方法的孔洞修復(fù)方法，該方法首先分析模型的幾何特征，然后采用正交梯度的方法來(lái)修復(fù)孔洞模型。Remacle等［4］提出了一種在孔洞邊界采樣點(diǎn)和缺失部分曲面定點(diǎn)之間建立調(diào)和映射的孔洞修補(bǔ)算法，最后通過(guò)最小化一個(gè)帶約束的全局能量函數(shù)來(lái)實(shí)現(xiàn)對(duì)缺失部分集合特征的修復(fù)。Park等［5］提出一種基于局部參數(shù)化的孔洞模型幾何和紋理的修復(fù)算法，該算法首先通過(guò)孔洞檢測(cè)算法檢測(cè)出孔洞，然后對(duì)待修復(fù)區(qū)域進(jìn)行局部參數(shù)化，再通過(guò)自選擇或者人機(jī)交互在模型上選擇一塊相似的區(qū)域來(lái)填補(bǔ)孔洞，最后通過(guò)解Poisson方程將選擇區(qū)域變形移植到孔洞區(qū)域，從而完成孔洞模型幾何和紋理的修復(fù)。

近年來(lái)，神經(jīng)網(wǎng)絡(luò)及支持向量機(jī)等智能化的機(jī)器學(xué)習(xí)方法被大量應(yīng)用到網(wǎng)格曲面的孔洞修補(bǔ)中，用以推斷缺失部分曲面的幾何特征。文獻(xiàn) ［6］使用混合徑向基函數(shù)模擬網(wǎng)格曲面缺失部分頂點(diǎn)坐標(biāo)與原始網(wǎng)格曲面頂點(diǎn)坐標(biāo)之間的函數(shù)關(guān)系，該算法能夠較好地恢復(fù)修補(bǔ)曲面的幾何特征；文獻(xiàn) ［7］在此基礎(chǔ)上將遺傳算法、模擬退火法和BP網(wǎng)絡(luò)方法相結(jié)合，在一定程度上提高了算法的效率，但是生成網(wǎng)格的質(zhì)量不高；文獻(xiàn) ［8］中的方法也存在類(lèi)似問(wèn)題。

綜合上述各種算法可以看出，由于孔洞部分的幾何信息已經(jīng)完全缺失，怎樣根據(jù)模型已有部分的幾何特點(diǎn)推斷出缺失部分的幾何信息，是網(wǎng)格曲面孔洞修復(fù)的根本難點(diǎn)所在。以神經(jīng)網(wǎng)絡(luò)和支持向量機(jī)等機(jī)器學(xué)習(xí)系統(tǒng)為基礎(chǔ)的修補(bǔ)算法正是針對(duì)這一核心難點(diǎn)提出的。本文采用在統(tǒng)計(jì)學(xué)習(xí)理論的基礎(chǔ)上發(fā)展起來(lái)的新一代機(jī)器學(xué)習(xí)方法［9］支持向量機(jī) （support vector machines，SVM）作為幾何信息回歸手段，結(jié)合當(dāng)前非常流行的網(wǎng)格曲面Laplacian坐標(biāo)，提出了一種全新的，針對(duì)三角網(wǎng)格曲面模型的孔洞修復(fù)算法，該算法突出的優(yōu)點(diǎn)就是能最大限度地使修補(bǔ)后的曲面在幾何風(fēng)格上與模型曲面原有部分保持一致。

1 主要數(shù)學(xué)工具

本文中使用的兩種核心數(shù)學(xué)工具是網(wǎng)格曲面的Laplacian坐標(biāo)和最小二乘支持向量機(jī)，其作用分別是描述網(wǎng)格曲面的幾何細(xì)節(jié)，和依據(jù)孔洞周?chē)糠智娴膸缀涡畔?duì)孔洞缺失部分曲面的幾何信息進(jìn)行回歸推斷。

1.1 網(wǎng)格曲面上的Laplacian坐標(biāo)

網(wǎng)格上的Laplacian坐標(biāo)是一種描述網(wǎng)格模型幾何細(xì)節(jié)的數(shù)學(xué)工具，頂點(diǎn)Vi處的Laplacian坐標(biāo)的計(jì)算公式如式（1），Vj是頂點(diǎn)Vi的局部片中的任意一點(diǎn)

Laplacian坐標(biāo)根據(jù)權(quán)重的不同分別是同一Laplacian坐標(biāo)和余切Laplacian坐標(biāo)，本文采用余切Laplacian坐標(biāo)

式 （2）、式 （3）是權(quán)重系數(shù)計(jì)算公式，采用式 （2）和式 （3）計(jì)算出來(lái)的就是頂點(diǎn)Vi處的余切Laplacian坐標(biāo)，記做δi。如式 （4）所示余切Laplacian在數(shù)值上線(xiàn)性逼近網(wǎng)格頂點(diǎn)處的平局曲率，方向上逼近該點(diǎn)的主曲率法向，可以非常準(zhǔn)確地描述網(wǎng)格頂點(diǎn)的局部片的幾何特性

圖1 頂點(diǎn)Vi的局部片

在計(jì)算網(wǎng)格頂點(diǎn)的Laplacian坐標(biāo)時(shí)，本文使用n×n階的Laplacian矩陣L，n為需要計(jì)算的網(wǎng)格區(qū)域內(nèi)的頂點(diǎn)數(shù)

式 （5）中的Lij的形式如式 （6）所示

式 （8）中的Vd是頂點(diǎn)坐標(biāo)在某一坐標(biāo)軸方向的分量，式 （7）中的Δd是Laplacian坐標(biāo)在某一坐標(biāo)方向的分量。

1.2 最小二乘支持向量機(jī)

最小二乘支持向量機(jī) （least－squares support vector machines）是由Suykens提出的，這種方法采用最小二乘線(xiàn)性系統(tǒng)作為損失函數(shù)，突出優(yōu)點(diǎn)是求解速度快，魯棒性強(qiáng)。最小二乘支持向量機(jī)在本文中的應(yīng)用是一個(gè)典型的回歸問(wèn)題，回歸的目標(biāo)是最優(yōu)化函數(shù)f（x）。本文給定訓(xùn)練樣本集（如式 （9）所示）為孔洞周?chē)鷧^(qū)域網(wǎng)格頂點(diǎn)的某一空間坐標(biāo)和Laplacian坐標(biāo)在同一坐標(biāo)軸方向的分量

式中：d∈ ｛x ，y ，z｝，l——訓(xùn)練樣本集中的樣本數(shù)量。平衡考慮回歸函數(shù)的復(fù)雜度和擬合誤差，依照結(jié)構(gòu)風(fēng)險(xiǎn)最小化原理，將本文中的回歸問(wèn)題描述為下面的約束優(yōu)化問(wèn)題

式 （11）中的w是權(quán)矢量，ei是松弛變量，γ為正則化參數(shù)，b為偏差量，為求解上述優(yōu)化問(wèn)題，需要建立Lagrange函數(shù) （如式 （12）所示）

式中：a＝ （a1，a2，…，al）為 Lagrange乘子，在 KKT 最優(yōu)化條件 （Karush－kuhn－Tuchker）下優(yōu)化函數(shù)各變量的偏導(dǎo)數(shù)，并在鞍點(diǎn)處得到零值［10］（如式 （13）所示）

將式 （13）寫(xiě)成矩形形式，并消去w和e可以得到式（14）

式中：Ωij＝K（xi，xj），求解得到式 （15）中的回歸函數(shù)

2 算法具體流程

本文算法分為以下4個(gè)部分：①孔洞邊界處理：包括孔洞邊界檢測(cè)與邊界邊預(yù)處理；②孔洞平面網(wǎng)格填充：包括建立孔洞特征平面、孔洞平面缺失采樣點(diǎn)估計(jì)與孔洞平面的網(wǎng)格生成；③三維缺失數(shù)據(jù)獲取：包括訓(xùn)練樣本集的采集和最小二乘支持向量機(jī)回歸函數(shù)求解；④空間網(wǎng)格生成：包括三維空間中孔洞區(qū)域修補(bǔ)網(wǎng)格的生成。

2.1 孔洞邊緣檢測(cè)

在三角網(wǎng)格曲面上孔洞通常是由那些只包含于一個(gè)三角形的邊組成的，這些邊都屬于模型的邊界邊。為了能夠加快邊界檢測(cè)速度，本文采用層次空間剖分技術(shù)來(lái)對(duì)點(diǎn)云數(shù)據(jù)進(jìn)行剖分。目前層次空間剖分技術(shù)主要包括八叉樹(shù)、四叉樹(shù)、BSP樹(shù)和K－D樹(shù)等，其中八叉樹(shù)表現(xiàn)出了優(yōu)良的特性。因此，采用八叉樹(shù)來(lái)加快邊界檢測(cè)的速度。邊界檢測(cè)算法的具體步驟如下：

步驟1 計(jì)算每個(gè)三角網(wǎng)格的中心以及它和三角形頂點(diǎn)間的聯(lián)系；

步驟2 通過(guò)三角網(wǎng)格的中心建立八叉樹(shù)；

步驟3 搜索八叉樹(shù)，檢測(cè)與它相關(guān)的邊的拓?fù)溥B接性，如果某條邊只屬于一個(gè)三角形，則是邊界邊。

步驟4 通過(guò)頂點(diǎn)的不同幾何連通性分辨出不同的邊界回路，表示出不同的孔洞。

2.2 孔洞邊界預(yù)處理

由于在實(shí)際測(cè)量中孔洞區(qū)域的數(shù)據(jù)信息非常少，因此，重構(gòu)后孔洞邊界處的三角片形狀可能會(huì)不理想，會(huì)出現(xiàn)一些狹長(zhǎng)的三角片。這些三角片上的邊同其它邊界邊不協(xié)調(diào)，在孔洞填充時(shí)會(huì)造成不好的影響。為了降低這種影響，本文先對(duì)孔洞邊界進(jìn)行預(yù)處理，將一些比較長(zhǎng)的邊進(jìn)行剖分，然后再對(duì)孔洞進(jìn)行修補(bǔ)。

設(shè)vi（i＝1，2，…，n）為孔洞邊界上的n個(gè)頂點(diǎn)，ε為孔洞邊界上的邊，則孔洞邊界預(yù)處理可以分為三步。

步驟1 計(jì)算孔洞邊界上所有邊的平均長(zhǎng)度Lave。

步驟2 遍歷所有的邊，如果存在某一條邊εij＝ （vi，vj）的長(zhǎng)度Lεij＞k·Lave，k為比例系數(shù) （本文中取k＝2），則計(jì)算邊εij的中點(diǎn)vnew。

步驟3 將vnew增加為孔洞的邊界點(diǎn)，同時(shí)刪除原來(lái)的邊界邊εij與三角片vvivj，增加邊界邊vivnew、vivnew和三角片vivvnew、vnewvvj，如圖2所示。最后迭代執(zhí)行上一步操作，直到邊界邊沒(méi)有變化為止。

圖2 孔洞邊界預(yù)處理

2.3 建立孔洞特征平面

孔洞特征多邊形的平面化是通過(guò)將特征多邊形投影到特征平面上實(shí)現(xiàn)的，所以特征平面的選取是一個(gè)關(guān)鍵步驟，設(shè)Ve＝ ｛vi，i＝1，2，…，m｝為孔洞的邊界點(diǎn)集合，本文選取與Ve中的所有點(diǎn)的距離平方和最小的平面作為特征平面，這樣的特征平面的法向量如式 （16）所示

式中：N——特征平面法向量，nsum——孔洞周邊所有三角面片法矢量的向量和，這樣的法平面選取可以保證N和最后修復(fù)得到的孔洞部分填充曲面的主法向量 （所有新生三角面片的矢量和）偏差最小，進(jìn)而保證對(duì)平面三角化后得到的修補(bǔ)網(wǎng)格頂點(diǎn)進(jìn)行重新定位時(shí)Laplacian坐標(biāo)的旋轉(zhuǎn)失真最小。

2.4 特征面平面上的孔洞填充

為了提高最終生成的三角網(wǎng)格質(zhì)量，本文提出一種基于步長(zhǎng)均值采樣的孔洞平面缺失數(shù)據(jù)估計(jì)算法來(lái)估計(jì)孔洞平面內(nèi)缺失的數(shù)據(jù)，其基本思想是使用兩組步長(zhǎng)相等的平行線(xiàn)對(duì)孔洞平面均值采樣，孔洞平面內(nèi)的交點(diǎn)即為估計(jì)的數(shù)據(jù)。算法的具體步驟如下：

步驟1 在特征平面內(nèi)求取孔洞特征多邊形的最小包圍盒；

步驟2 計(jì)算投影后邊界邊的平均邊長(zhǎng)珚Lb；

步驟3 以珚Lb為步長(zhǎng)，用兩組互成60°的平行線(xiàn)對(duì)特征多邊形的包圍盒均勻采樣，如圖3所示。

圖3 孔洞特征多邊形均值采樣

步驟4 采用多邊形內(nèi)外點(diǎn)判斷方法對(duì)步驟3中的交點(diǎn)進(jìn)行判斷，保留特征多邊形的內(nèi)點(diǎn)；

步驟5 分別計(jì)算步驟4中所得估計(jì)點(diǎn)與邊界邊和邊界點(diǎn)間的距離le和lp，如果le＜0.2珚Lb或是lp＜0.2珚Lb，則刪除對(duì)應(yīng)的估計(jì)點(diǎn)。

步驟6 最終所剩下的點(diǎn)則為所要計(jì)算的缺失數(shù)據(jù)點(diǎn)，如圖4（a）所示。

圖4 孔洞平面網(wǎng)格填充

在特征平面內(nèi)通過(guò)采樣得到缺失數(shù)據(jù)點(diǎn)的坐標(biāo)值之后，需要生成三角網(wǎng)格以方便后續(xù)的處理。直接在局部坐標(biāo)系下將上一步求得的數(shù)據(jù)點(diǎn)與孔洞的邊界點(diǎn)一起進(jìn)行Delaunay三角剖分，得到三角網(wǎng)格，完成孔洞平面的網(wǎng)格填充，填充后的網(wǎng)格如圖4（b）所示。

2.5 為最小二乘支持向量機(jī)選取學(xué)習(xí)樣本

迄今為止所有的孔洞修復(fù)算法都是依據(jù)已有網(wǎng)格曲面的幾何特征去推測(cè)缺失部分曲面的幾何特征，本文算法就是選取孔洞周?chē)狞c(diǎn)及改點(diǎn)上Laplacian坐標(biāo)作為最小二乘支持向量機(jī)的訓(xùn)練樣本集，在實(shí)現(xiàn)過(guò)程中以孔洞邊界為中心不斷向外擴(kuò)展，逐層采集相鄰幾層的三角形頂點(diǎn)，直到采集到滿(mǎn)足需要的層數(shù)為止。

假設(shè)Ve為孔洞邊界B 上的頂點(diǎn)集合，SVj（j＝0，1，…k）為孔洞外面第j層的頂點(diǎn)集合，樣本采集的具體算法如下：

步驟1 在邊界上取出一點(diǎn)vi，然后尋找與它相關(guān)的三角面片頂點(diǎn)，如果該頂點(diǎn)不在Ve中，則將它加入到孔洞的下一層頂點(diǎn)集SV1中；

步驟2 遍歷數(shù)組Ve中的每一個(gè)頂點(diǎn)，重復(fù)步驟1，這樣就采集到了由孔洞多邊形向外的第一層三角片頂點(diǎn)，并存放在頂點(diǎn)集SV1中；

步驟3 對(duì)頂點(diǎn)集SVj（j＝1，2…k－1）中的頂點(diǎn)，尋找與它相關(guān)的三角形的頂點(diǎn)，并將不在SVj和SVj－1中的點(diǎn)加入到頂點(diǎn)集SVj＋1中；

步驟4 迭代執(zhí)行步驟3，直至采集達(dá)到所設(shè)定的層數(shù)為止，則SVj（j＝0，1，…k）中所有的點(diǎn)及其Laplacian坐標(biāo)即構(gòu)成了訓(xùn)練樣本集，采樣點(diǎn)分布如圖5所示。

圖5 為最小二乘為支持向量機(jī)選取的樣本點(diǎn)集

2.6 填充部分頂點(diǎn)的重新定位

將第2.5小節(jié)中得到的訓(xùn)練樣本集代入最小二乘支持向量機(jī)，按照2.2小節(jié)中的方法可以得到3個(gè)回歸函數(shù)：fx（x），fy（x），fz（x）。

設(shè)定Vnew為2.4小節(jié)中填充的新頂點(diǎn)集合。式 （17）中的V′nd是對(duì)Vnew中的頂點(diǎn)重新定位后的頂點(diǎn)坐標(biāo)在某一坐標(biāo)軸方向上的分量，重新定位的過(guò)程是對(duì)網(wǎng)格頂點(diǎn)坐標(biāo)的3個(gè)分量分別進(jìn)行求解，求解的是式 （18）中所示的n×n的線(xiàn)性系統(tǒng)，n為Vnew中頂點(diǎn)數(shù)

式 （18）中的LVn是在2.4小節(jié)中剖分得到的修補(bǔ)三角網(wǎng)格上，按照2.1小節(jié)中的方法建立的余切Laplacian矩陣，Δ′dn是重新定位后的修補(bǔ)區(qū)域網(wǎng)格頂點(diǎn)余切Laplacian坐標(biāo)的某一分量構(gòu)成的向量，是按照式 （19）構(gòu)造的。式（19）中的fd（v′nd）是由回歸函數(shù)推測(cè)出的缺失區(qū)域網(wǎng)格頂點(diǎn)的Laplacian坐標(biāo)的某一分量。

3 算法評(píng)價(jià)及實(shí)驗(yàn)結(jié)果

為了檢驗(yàn)本文算法的實(shí)際性能，以Visual C＋＋2010為開(kāi)發(fā)平臺(tái)，結(jié)合OpenGL圖形庫(kù)實(shí)現(xiàn)了文中算法，并在head模型和happyvrip模型上進(jìn)行了對(duì)比試驗(yàn)，試驗(yàn)結(jié)果如下：

圖6 中的模型就是head模型，在頭發(fā)部分人為制造了一個(gè)孔洞，作為待修補(bǔ)孔洞。試驗(yàn)中選取的對(duì)比算法為傳統(tǒng)的直接修復(fù)算法和文獻(xiàn) ［11］中的算法，之所以選擇文獻(xiàn) ［11］算法是因?yàn)樵撍惴ê捅疚乃惴ǘ际腔谧钚《酥С窒蛄繖C(jī)的，區(qū)別主要是文獻(xiàn) ［11］算法直接將網(wǎng)格頂點(diǎn)的一個(gè)歐氏坐標(biāo)分量做回歸函數(shù)的輸出，而本文引入了網(wǎng)格頂點(diǎn)的余切Laplacian坐標(biāo)作為回歸函數(shù)的輸出，再以由回歸函數(shù)推斷得到的缺失部分網(wǎng)格曲面的Laplacian坐標(biāo)為基礎(chǔ)建立線(xiàn)性系統(tǒng)求解得到缺失部分的網(wǎng)格頂點(diǎn)坐標(biāo)。

圖6 原始head模型和制造孔洞后的head模型

圖7 中是使用原始的直接修復(fù)算法得到的修付效果，明顯比較粗糙。圖8中是使用文獻(xiàn) ［11］中算法得到的修復(fù)效果，效果很好，修復(fù)區(qū)域的網(wǎng)格曲面質(zhì)量也不錯(cuò)。圖9是使用本文算法得到的修復(fù)效果，修復(fù)區(qū)域的網(wǎng)格質(zhì)量也很好，同時(shí)可以看到相對(duì)于文獻(xiàn) ［11］中的算法，使用本文算法的到的修補(bǔ)曲面與原始曲面的融合度明顯更高，效果與沒(méi)有制造孔洞前的原始模型曲面幾乎沒(méi)有差別。

圖10 是使用本文算法對(duì)happyvrip模型的修復(fù)效果，這個(gè)試驗(yàn)的目的是檢驗(yàn)本文算法在高曲率曲面上的修復(fù)效果，在具體實(shí)驗(yàn)中，為了應(yīng)付高曲率環(huán)境，加大了2.5小節(jié)中支持向量機(jī)的訓(xùn)練樣本的選取范圍 （head試驗(yàn)中的選取范圍是孔洞外的3層點(diǎn)，happyvrip試驗(yàn)中為7層點(diǎn)），從試驗(yàn)結(jié)果看本文算法對(duì)高曲率曲面的修復(fù)效果也是比較理想的。

4 結(jié)束語(yǔ)

本文算法的創(chuàng)新點(diǎn)在于引入了網(wǎng)格頂點(diǎn)的余切Laplacian坐標(biāo)作為回歸函數(shù)的輸出，并通過(guò)求解線(xiàn)性系統(tǒng)的方式對(duì)修復(fù)區(qū)域的網(wǎng)格頂點(diǎn)進(jìn)行二次定位，這樣做的優(yōu)勢(shì)是發(fā)揮了Laplacian坐標(biāo)對(duì)網(wǎng)格曲面幾何細(xì)節(jié)的刻畫(huà)能力，可以使修復(fù)后的網(wǎng)格曲面與原始曲面融合的更加自然。

本文算法的局限性在于對(duì)面積過(guò)大的孔洞的修復(fù)效果不理想，這也是目前所有的曲面修復(fù)算法面對(duì)的共同問(wèn)題，對(duì)付這個(gè)問(wèn)題的根本途徑是在構(gòu)建修復(fù)區(qū)域網(wǎng)格時(shí)不僅要以孔洞周?chē)娴膸缀翁匦宰鳛橐罁?jù)，還要考慮模型的整體幾何特征，這也是本文下一步研究和改進(jìn)的方向。

［1］Geuzaine C，Remacle J.Gmsh：A three－dimensional finite element mesh generator with built－in pre－and post－processing facilities ［J］.International Journal for Numerical Methods in Engineering，2009，79 （11）：1309－1331.

［2］Wu X，Wang M，Han B.An automatic hole－？lling algorithm for polygonal meshes ［J］.Computer－Aided Design and Applications，2008，5 （6）：889－899.

［3］Marchandise E，Piret C.CAD and mesh repair with radial basis functions［J］.Journal of Computational Physics，2012，231 （5）：2376－2387.

［4］Remacle J，Geuzaine C，Compère G，et al.High quality surface remeshing using harmonic maps ［J］.International Journal for Numerical Methods in Engineering，2010，83 （4）：403－425.

［5］Marchandise E，Carton de Wiart C，Vos W，et al.High quality surface remeshing using harmonic maps—Part II：Surfaces with high genus and of large aspect ratio ［J］.International Journal for Numerical Methods in Engineering，2011，86（11）：1303－1321.

［6］Wright G B，F(xiàn)lyer N，Yuen D.A hybrid radial basis function－pseudo－spectral method for thermal convection in a 3dspherical shell［J］.Chemistry Geophysics Geosystems，2010，11 （7）：1－18.

［7］CHEM Jie，GAO Chenghui，HE Bingwei.Automatic recognition of boundary features of non－closed triangulation model［J］.Machinery Design ＆ Manufacture，2011 （11）：89－94.

［8］Macdonald C B，Ruuth S J.The implicit closest point method for the numerical solution of partial differential equations on surfaces［J］.SIAM Journal on Scienti？c Computing，2009，31（6）：4330－4350.

［10］WU Dehui.Intelligent prediction model for surface roughness in milling ［J］.Computer Integrated Manufacturing Systems，2007，13 （6）：1137－1141 （in Chinese）.

［11］LIU Deping，YU Shuijing.Hole repairing in triangular meshes based on least－squares support vector machines ［J］.Computer Integrated Manufacturing Systems，2009，15 （9）：1868－187（in Chinese）.［劉德平，于水晶.基于最小二乘支持向量機(jī)的三角網(wǎng)格修補(bǔ)算法 ［J］.計(jì)算機(jī)集成制造系統(tǒng)，2009，15 （9）：1868－187.］