許華強(qiáng)， 徐 崗， 胡維華， 汪國(guó)昭

（1. 杭州電子科技大學(xué)計(jì)算機(jī)學(xué)院圖形圖像研究所，浙江 杭州 310018；2. 浙江大學(xué)數(shù)學(xué)系圖像圖形研究所，浙江 杭州 310027）

隨著CAD/CAE 領(lǐng)域的快速發(fā)展和降低設(shè)計(jì)成本的需要，設(shè)計(jì)與分析的一體化已成為大勢(shì)所趨[1]。目前，有限元分析已經(jīng)成為CAE 領(lǐng)域?qū)AD 產(chǎn)品進(jìn)行仿真模擬的主流方法。然而，當(dāng)今設(shè)計(jì)與模擬分析階段幾何數(shù)據(jù)模型的不一致性，嚴(yán)重影響了求解精度和計(jì)算效率。據(jù)統(tǒng)計(jì)，數(shù)據(jù)交換可以占到CAE 軟件大約70%的運(yùn)行時(shí)間，而整個(gè)CAD/CAE 過(guò)程80%左右的工作量用于離散網(wǎng)格生成階段。此外，由于有限元分析中的計(jì)算域只是對(duì)CAD模型的逼近表示，這必然導(dǎo)致模擬分析及優(yōu)化結(jié)果存在一定的誤差。對(duì)于某些模擬仿真中的物理問(wèn)題而言，由此所導(dǎo)致的誤差可能會(huì)非常顯著，因此，如何實(shí)現(xiàn)真正意義上的CAD/CAE 系統(tǒng)的無(wú)縫融合，便成為CAD/CAE 向前發(fā)展過(guò)程中亟需解決的問(wèn)題。

等幾何分析方法的提出為實(shí)現(xiàn)CAD/CAE系統(tǒng)的無(wú)縫融合提供了新途徑[2]。相對(duì)于經(jīng)典的有限元分析方法，該方法的主要優(yōu)勢(shì)在于：

1）具有與CAD階段統(tǒng)一的幾何數(shù)學(xué)模型，或者兩個(gè)階段的幾何數(shù)學(xué)模型具有互用性，無(wú)需進(jìn)行容錯(cuò)轉(zhuǎn)換；

2）采用高階基函數(shù)作為形函數(shù)，具有良好的收斂性；

3）求解所得到的物理場(chǎng)具有高階連續(xù)性，對(duì)某些薄殼問(wèn)題的求解可達(dá)到比較高的精度；

4）實(shí)現(xiàn)了CAD 模型與計(jì)算域的精確轉(zhuǎn)換，提高了模擬分析及優(yōu)化設(shè)計(jì)的準(zhǔn)確度；

5）在細(xì)化過(guò)程中（節(jié)點(diǎn)插入或升階操作）可以使計(jì)算域的幾何保持不變[3]，這一特性適用于流體和固體的耦合問(wèn)題。等幾何分析方法的提出為已趨成熟的幾何造型領(lǐng)域注入了新的活力，目前已成為幾何設(shè)計(jì)與計(jì)算領(lǐng)域新的關(guān)注熱點(diǎn)，由于等幾何樣條分析的研究目前處于起步階段，該方法在應(yīng)用的廣度與模擬的自適應(yīng)性等方面仍存在一定的局限性。

等幾何分析方法雖然無(wú)須由CAD 模型的邊界曲面表示生成離散的四面體網(wǎng)格作為計(jì)算域，但仍然需要由其邊界表示生成連續(xù)的參數(shù)樣條體表示以完成計(jì)算域的構(gòu)造。這一體參數(shù)化過(guò)程與有限元方法中的網(wǎng)格生成過(guò)程的不同之處在于：我們需要在精確保持邊界曲面的前提下，將CAD 模型的內(nèi)部用若干個(gè)三變量參數(shù)體填充起來(lái)，而有限元分析中的網(wǎng)格生成過(guò)程得到的則是該CAD模型的四面體網(wǎng)格逼近表示。CAD 模型的體參數(shù)化結(jié)果對(duì)模擬分析的精度有著重要影響[4]，這一問(wèn)題在Cottrell 和 Hughes 等所著的《Isogeometric Analysis: Toward Integration of CAD and FEA》一書中被列為等幾何分析面臨的最富有挑戰(zhàn)性的公開問(wèn)題之一。目前，國(guó)際上對(duì)這一問(wèn)題的研究屈指可數(shù)。Martin 等利用離散調(diào)和函數(shù)方法來(lái)完成對(duì)四面體網(wǎng)格的參數(shù)體逼近[5]；Escobar 等利用參數(shù)域剖分和插值逼近技術(shù)來(lái)構(gòu)造三角網(wǎng)格模型的三變量T 樣條體參數(shù)化[6]；Zhang 等利用參數(shù)映射和八叉樹技術(shù)來(lái)構(gòu)造三角網(wǎng)格模型的三變量T樣條體參數(shù)化[7]。上述3種方法雖然均可有效生成體參數(shù)化，但它們的缺點(diǎn)是顯而易見的，即其輸入是三角網(wǎng)格模型或者四面體網(wǎng)格模型，而參數(shù)樣條表示是CAD系統(tǒng)中的標(biāo)準(zhǔn)形式，顯然違背了等幾何分析方法的精髓——精確幾何表示。理想的方法應(yīng)該是以CAD 模型的邊界樣條曲面表示作為輸入，構(gòu)造其內(nèi)部空間的參數(shù)體表示。本文提出了基于模板的由邊界曲面生成參數(shù)體的新方法。首先給出了插值邊界曲面的Coons體的參數(shù)表達(dá)式。由該表達(dá)式的離散表示，內(nèi)部控制頂點(diǎn)可表為邊界控制頂點(diǎn)的線性組合。然后由其離散表示可得到Coons模板，并將Coons模板推廣到統(tǒng)一形式，這為內(nèi)部控制頂點(diǎn)的生成提供了更多選擇。本文通過(guò)幾個(gè)熱傳導(dǎo)問(wèn)題的例子對(duì)由不同的模版所得到的不同的體參數(shù)化結(jié)果及其對(duì)等幾何分析結(jié)果的影響進(jìn)行了比較分析。

1 Coons 參數(shù)體

給定6張邊界曲面S(0,v,w),S(1,v,w) ,S(u,0,w),S(u, 1,w),S(u,v, 0),S(u,v, 1), 插值上述6張曲面的Coons參數(shù)體S(u,v,w)可按如下方式構(gòu)造：

圖1 Coons參數(shù)體

其中，α(u)β(v)γ(w)為滿足如下條件的參數(shù)函數(shù)：

通過(guò)直接計(jì)算可證明參數(shù)體S(u,v,w)插值6張邊界參數(shù)曲面，它可以看作插值4條邊界曲線的Coons曲面的推廣。圖1給出了一個(gè)Coons參數(shù)體S(u,v,w)的例子。

若給定的邊界曲面為B樣條曲面，并且相對(duì)的B樣條曲面具有相同的次數(shù)、控制頂點(diǎn)數(shù)目和節(jié)點(diǎn)向量，α(u)β(v)γ(w)為參數(shù)多項(xiàng)式函數(shù)，則按照Coons 參數(shù)體的定義所生成的參數(shù)體仍為B樣條形式。若α(u) =u,β(v)=v，γ(w) =w,并且邊界控制頂點(diǎn)按圖2給出，則每一個(gè)控制頂點(diǎn)Pi,j,k可賦予一組參數(shù)值(u,v,w) =(i/l,j/m,k/n)， 內(nèi)部的控制頂點(diǎn)Pi,j,k可按照如下的離散形式構(gòu)造：

則相應(yīng)的B樣條參數(shù)體可記為

其中，Ni(u),Nj(v),Nk(w)為B樣條基函數(shù)，其節(jié)點(diǎn)向量由邊界B樣條曲面給出。

圖2 排列方式

2 B樣條參數(shù)體生成的模板方法

本文所提出的模板實(shí)際上表達(dá)了每一個(gè)內(nèi)部控制頂點(diǎn)與其相鄰控制頂點(diǎn)的線性關(guān)系。對(duì)B樣條Coons體的每一個(gè)內(nèi)部控制頂點(diǎn)Pi,j,k來(lái)說(shuō)，其與相鄰控制頂點(diǎn)之間的關(guān)系可表示為

上述關(guān)系可以表示為圖3(a)中簡(jiǎn)潔的模板表示。對(duì)一個(gè)具有(l-1)× (m-1)× (n-1)個(gè)內(nèi)部控制頂點(diǎn)的B樣條參數(shù)體來(lái)說(shuō)，我們將得到(l- 1)× (m-1)× (n -1 )個(gè)線性方程。其所組成的線性系統(tǒng)的系數(shù)矩陣為一稀疏矩陣。在B樣條體的內(nèi)部，其表現(xiàn)為內(nèi)部待求控制頂點(diǎn)之間的關(guān)系，而在其邊界，則表現(xiàn)為邊界已知控制頂點(diǎn)與內(nèi)部待求控制頂點(diǎn)之間的關(guān)系。

等式(2)所表示的離散Coons模板可以擴(kuò)展到如下形式：

上述關(guān)系可以表示為圖3(b)中簡(jiǎn)潔的模板表示。該模板表示方法提供了3個(gè)參數(shù)α,,βγ，增加了B樣條參數(shù)體生成的自由度。為使內(nèi)部待求的控制頂點(diǎn)位于邊界已知的控制頂點(diǎn)所組成的凸包內(nèi)，參數(shù)α,,βγ應(yīng)滿足如下約束條件

圖3 模板表示方法

該模板表示形式可以記為Mα,β(Pi,j,k)。對(duì)于給定的邊界B樣條曲面，用戶可選擇不同的參數(shù)α,β來(lái)構(gòu)造B樣條參數(shù)體，γ可根據(jù)α,β的值由式（4）求出。例如，離散Coons模板可記為其中離散Laplacian 模板可記為，其中離散調(diào)和模板可記為其中

圖4給出了在邊界B樣條曲面相同的情形下，依據(jù)上述3種不同的模板方法所生成的B樣條參數(shù)體的例子。為將3種模板進(jìn)行比較，本文給出了所生成的B樣條參數(shù)體的等參數(shù)曲面分布，可以發(fā)現(xiàn)離散調(diào)和模板所生成的B樣條參數(shù)體的等參數(shù)曲面分布最為均勻。另外，本文將生成的體參數(shù)化結(jié)果在三維等幾何分析中進(jìn)行了測(cè)試。所求解的物理仿真問(wèn)題為如下的三維熱傳導(dǎo)問(wèn)題：

圖4 體參數(shù)化實(shí)例

3 結(jié) 論

三維等幾何分析中的計(jì)算域通常表示為參數(shù)樣條體。本文提出了一種由邊界曲面生成參數(shù)體的模板方法。首先給出了插值邊界參數(shù)曲面的Coons體的參數(shù)表達(dá)式。由其離散表示可得到Coons模板，并將Coons模板推廣到統(tǒng)一形式，這為內(nèi)部控制頂點(diǎn)的生成提供了更多自由度。通過(guò)模板方法，內(nèi)部控制頂點(diǎn)可通過(guò)求解一個(gè)稀疏線性系統(tǒng)得到。本文通過(guò)熱傳導(dǎo)問(wèn)題的實(shí)例對(duì)由不同的模板所得到的體參數(shù)化結(jié)果及其對(duì)等幾何分析的影響進(jìn)行了比較分析。

[1]高曙明, 何發(fā)智. 異構(gòu)CAD系統(tǒng)集成技術(shù)綜述[J]. 計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào), 2009, 21(5): 561-568.

[2]Hughes T J R, Cottrell J A, Bazilevs Y. Isogeometric analysis: CAD, finite elements, NURBS, exact geometry, and mesh refinement [J]. Computer Methods in Applied Mechanics and Engineering, 2005,194(39-41): 4135-4195.

[3]徐 崗, 王毅剛, 胡維華. 等幾何分析中的r-p型細(xì)化方法[J]. 計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào), 2011,(12): 2019-2024.

[4]Xu Gang, Mourrain B, Duvigneau R, et al.Parameterization of computational domain in isogeometric analysis: methods and comparison [J].Computer Methods in Applied Mechanics and Engineering, 2011, 200(23-24): 2021-2031.

[5]Martin T, Cohen E, Kirby R M. Volumetric parameterization and trivariate B-spline fitting using harmonic functions [J]. Computer Aided Geometric Design, 2009, 26(6): 648-664.

[6]Escobar J M, Cascon J M, Rodr?guez E, et al. A new approach to solid modeling with trivariate T-splines based on mesh optimization [J]. Computer Methods in Applied Mechanics and Engineering, 2000, (45-46);2011, 3210-3222.

[7]Zhang Y, Wang W, Hughes T J R. Solid T-spline Construction from boundary representations for genus-zero geometry [J]. Computer Methods in Applied Mechanics and Engineering, 2012,(249-252):185-197.