01  前奏

要理解等幾何分析，首先要了解有限元方法發(fā)展的來龍去脈。

達芬奇說過，數(shù)學(xué)是一切科學(xué)的基礎(chǔ)。我們所生活的這個世界，大到宇宙天體的宏觀尺度，小到肉眼不可見粒子的微觀尺度，但凡總結(jié)出的原理規(guī)律，都可以用一定的數(shù)學(xué)語言描述；同樣，工程中的許多問題都可以用微分方程（常微分方程、偏微分方程）、加以特定的邊界條件（力學(xué)條件、位移條件、溫度條件等）表達。例如，學(xué)過材料力學(xué)和結(jié)構(gòu)力學(xué)的人知道，簡支梁的自由彎曲變形可以用以下公式描述：



只要知曉施加在梁上的特定邊界條件（例如端點位移為0，或受力點載荷值），便可通過求解該常微分方程，精確得到某種工況下（載荷）和梁的彎曲程度（位移/撓度）的關(guān)系。


（簡支梁的自由彎曲。注：上文公式僅考慮歐拉梁，忽略鐵木辛柯梁的剪切效應(yīng)）

稍微復(fù)雜點的，例如求解物體受載后的內(nèi)部應(yīng)力/應(yīng)變，彈性力學(xué)告訴我們，可以通過構(gòu)造三大平衡方程（體現(xiàn)力和位移的關(guān)系特性）+物理方程（體現(xiàn)物體本身應(yīng)力、應(yīng)變的材料特性）+幾何方程（體現(xiàn)位移和應(yīng)變的關(guān)系特性）的閉環(huán)系統(tǒng)，外加邊界條件，即可聯(lián)立求解方程組，得到域內(nèi)任意點的力或位移等信息。


（物體受載內(nèi)部微觀單元體的三向應(yīng)力狀態(tài)）

但是，能用公式描述出萬事萬物來是一方面，進一步能否求解出方程是另一方面。實際上，只有極少數(shù)簡單邊界條件或簡化模型，才能通過嚴格的數(shù)學(xué)公式推導(dǎo)得到精確解（即解析解）。大多數(shù)情況下，我們根本沒有辦法解出方程，就好比我們描述某個人沿著某條路跑步，知道他用了多久跑了多遠，但是無法得到過程中的每一秒他在哪里、速度是多少等等。因此，人們只能曲線救國，通過一定的數(shù)值方法求取其近似解。

早在3世紀，魏晉數(shù)學(xué)家劉徽首創(chuàng)割圓術(shù)，是人類數(shù)學(xué)史上“化整為零”較早的思想萌芽。17世紀，牛頓和萊布尼茨獨立提出微積分方法，證明了整體對局部的可加性，從此人們對描述客觀世界復(fù)雜連續(xù)性找到了一個嶄新的思路，即化整為零——無窮逼近——化零為整的方法。之后，高斯提出的加權(quán)余量法、拉格朗日提出的泛函分析，再到后來瑞利、里茲首次提出對全定義域用位移函數(shù)來表述，都對后來的有限元提出建立了理論基礎(chǔ)。


（割圓術(shù)：割之彌細，失之彌少，割之又割，以至于不可割，則與圓合體，而無所失矣。——劉徽）

有限元方法（FEM，F(xiàn)inite Element Method）首先應(yīng)用在航空航天領(lǐng)域。1950年，艾格瑞斯（Argyris）提出力學(xué)分析的結(jié)構(gòu)位移方法，被譽為有限元的發(fā)明人和創(chuàng)始人；1956年，波音公司在紐約舉行的航空學(xué)會上介紹了將矩陣位移法推廣到平面應(yīng)力問題的求解上，提出把整個結(jié)構(gòu)劃分成一個個三角形和矩形單元進行分析，這實際上就是有限元方法的早期應(yīng)用；1960年，克拉夫（Clough）首次提出有限元這一術(shù)語（《The Finite Element in Plane Stress Analysis》）；1967年，辛克維奇（Zienkiewicz.O.C.，又有人叫杰卡維奇）出版了第一本有限元專著《The Finite Element Method》。與國外有限元發(fā)展各階段同時期甚至更早，胡海昌（廣義變分原理）、錢偉長（拉格朗日乘子法）、馮康（先于國外獨立提出有限元分析收斂性的理論）等國內(nèi)學(xué)者也對有限元進行了深入研究。隨著國內(nèi)外學(xué)者對單元構(gòu)造、微分方程的數(shù)值求解方法、收斂性和精度等有限元基礎(chǔ)理論研究的完善，以及計算機算力的提高和相關(guān)CAE軟件的大力發(fā)展，工程上真正開啟了有限元方法馳騁CAE領(lǐng)域的時代。


（有限元方法發(fā)展的先驅(qū)，從左到右：John Argyris——獲得16個國家授予的杰出科技與工程獎；Ray W.Clough——美國兩院院士，我國外籍院士；Zienkiewicz.O.C.——號稱5院院士，我國外籍院士，其撰寫的有限元開山制作至今仍是新手入門的不二之選；馮康：中科院院士，我國有限元理論的拓路者。很可惜，以上大佬均未獲諾貝爾獎，因為有限元屬于計算數(shù)學(xué)領(lǐng)域，而諾獎不設(shè)數(shù)學(xué)獎。）


（史上第一本有限元著作：《The Finite Element Method》——Zienkiewicz.O.C.）

計算固體力學(xué)告訴我們，有限元分析的一般流程為：結(jié)構(gòu)離散（網(wǎng)格劃分）——單元構(gòu)造（確定位移插值方法）——單元分析（局部剛度矩陣）——整體結(jié)構(gòu)分析（所有單元組裝成整剛矩陣）——施加邊界條件（邊值、初值）——構(gòu)造方程、解方程（數(shù)值方法求近似解）。以工程上最常用的平面四節(jié)點矩形單元為例，我們選取四個節(jié)點的位移作為基本量，其余位置（邊界、域內(nèi)）點的位移通過形函數(shù)和四個節(jié)點的位移進行雙線性插值近似得到：


（平面四節(jié)點單元，應(yīng)變?yōu)榫€性）


即


其中形函數(shù)為：



通過某點的位移求導(dǎo)即可得到該點的應(yīng)變，通過材料本構(gòu)特性即可由應(yīng)變得到應(yīng)力，最終得到以位移為基本量表示的所有力學(xué)特性（當然除了四個節(jié)點是精確的，其余點都是近似插值得到）。將其代入上文提到過的彈性力學(xué)三大方程，即可求解出該單元的剛度矩陣（當然單元剛度矩陣也可由虛功原理或變分原理求導(dǎo)得出，具體不再展開論述）。


（筆者師公所著，有限元和計算固體力學(xué)的經(jīng)典教材）

我們說過，有限元是通過將真實結(jié)構(gòu)離散成很多簡單的單元求解的，按照上文得到某個單元的單剛矩陣后，將所有單元的剛度矩陣組裝成真實結(jié)構(gòu)的整剛矩陣，再代入邊界條件，即可求解出該結(jié)構(gòu)所有單元節(jié)點的位移，進而求得所有的力學(xué)量。



當然，實際工程涉及的模型遠比本例復(fù)雜，所用的單元類型也不僅限于本例中的平面矩形。以汽車為例，整車離散后的單元個數(shù)、節(jié)點個數(shù)是千萬級的，裝配所有單元節(jié)點后所構(gòu)造的剛度矩陣，其維數(shù)規(guī)模遠超人力所能計算的范圍；不同的單元類型（選擇什么樣的形函數(shù)、什么樣的積分點進行插值），決定了求解的精度和效率；不同的矩陣運算和積分算法，決定了剛度矩陣的規(guī)模；不同的邊界條件和工況，決定了采用何種求解算法（隱式、顯示等）。

時至今日，在有限元理論發(fā)展和CAE軟件高度集成的雙buff加持下，今天的工程師已不需要從零開始構(gòu)造單元、組裝矩陣、求解方程（當然這也是不可能的，單靠人力怎么可能求解超大規(guī)模矩陣）。我們需要做的，只不過將設(shè)計師做好的三維數(shù)據(jù)，通過前處理工具進行網(wǎng)格劃分、選擇單元類型和積分方式，然后做點必要的連接、接觸等邊界條件設(shè)置，再丟給求解器進行計算即可，依靠工程經(jīng)驗，幾乎不需要過多理解有限元本身的豐富內(nèi)涵。無論如何，我們所做的所有努力，實質(zhì)都是在解牛頓第二定律所表示的靜力學(xué)/動力學(xué)方程，只不過公式里的向量是用非常復(fù)雜形式的偏微分表示，由于我們無法直接求解方程，只能采用數(shù)值方法（有限元），將其轉(zhuǎn)化成有限個單元組成的超大型矩陣運算。其中涉及的矩陣運算、數(shù)據(jù)存儲、方程求解、各種各樣單元構(gòu)造的方式和形函數(shù)的推算，其實都凝聚了幾代學(xué)者的努力，這里就不展開講了，詳情可查閱附錄相關(guān)資料。



（LS-DYNA中選擇某個單元類型，其實就是在選擇前輩們?yōu)槲覀兺茖?dǎo)出的形函數(shù)、積分點、積分方式等，不同的單元類型，直接決定了模擬真實物體的力學(xué)特性逼真度如何，也導(dǎo)致所得結(jié)果的精度和求解效率有很大不同）

為什么要花一些時間介紹有限元呢？因為等幾何分析理論的提出，正是伴隨著有限元的優(yōu)勢和缺點而誕生的。

02  萌芽

有限元的優(yōu)勢在于，無論多么復(fù)雜的形體，都可以離散成簡單的單元進行逼近，且邊界條件非常容易施加（例如施加本質(zhì)邊界條件，只需約束單元節(jié)點的位移即可）。然而這恰恰帶來了有限元與生俱來的缺點：Sandia國家實驗室曾做過一項統(tǒng)計，工程上80%以上的精力都要用來對模型前處理，即網(wǎng)格劃分、單元構(gòu)造、模型搭建的過程，工程師真正用來聚焦于問題本身和結(jié)構(gòu)優(yōu)化的時間少之又少；其次，有限元畢竟是一種近似的數(shù)值解法，原本連續(xù)光滑的模型被硬生生分割成有限個單元，不但精度大打折扣，而且單元之間的光滑性、連續(xù)性無法保證，需要后期進行特殊處理；想要提高結(jié)果的精度，要么用階次更高的高精度單元（p-收斂），要么將網(wǎng)格劃分的更加細密、更加逼近原始光滑模型（h-收斂），而這無疑增加了人的工作量和電腦的求解時間；不同的網(wǎng)格劃分方式、單元類型的選擇，直接影響求解精度和效率，工程分析經(jīng)驗性太大，即同一個工況的問題交給不同的工程師來做，可能得到相差甚遠的結(jié)論；對于大變形問題、裂紋擴展問題等容易產(chǎn)生網(wǎng)格畸變和網(wǎng)格重構(gòu)的障礙；最后，一旦上游的設(shè)計思路變更導(dǎo)致CAD模型更改，有限元模型必須跟著同步更新（例如重新劃分網(wǎng)格、重新搭建模型），即伴隨不同版本的CAD數(shù)據(jù)必須同步更新不同版本的CAE模型。

以上所說的每一條痛點，相信在座從事CAE仿真的各位攻城獅都深有體會——每次網(wǎng)格劃分質(zhì)量太差、模型怎么調(diào)試都算不通、項目節(jié)點逼近又一輪分析任務(wù)在等待，都想質(zhì)問自己一聲當初為何腦子進水要從事仿真分析工作。那么，有沒有什么方法，能夠保留有限元的優(yōu)勢、而又同時解決它的缺點甚至取代有限元呢？

學(xué)過數(shù)值分析和線性代數(shù)的人知道，有限元本質(zhì)上還是一種解方程的數(shù)值方法，除此之外，還有很多近似解法，例如：離散元法（DEM，Discrete Element Method，解決不連續(xù)介質(zhì)問題），邊界元法（BEM，Boundary Element Method，只需離散邊界，域內(nèi)無需離散，降維使問題更簡單），無網(wǎng)格法（Mesh-less Method，無需劃分網(wǎng)格，在邊界和域內(nèi)“撒點”，全局插值），還有有限差分法、有限體積法、譜方法等。但是，須知道任何一種數(shù)值求解方法都必須依附于計算機，其中最需要考慮的問題便是效率和精度。盡管當今CAE軟件也越來越集成其他的數(shù)值計算方法，但這些方法都由于缺乏精度或效率某一方面、或不具有普適性、或?qū)嶋H應(yīng)用代價太大，終究無法撼動有限元法的霸主地位，僅僅停留在某些特殊問題的研究范圍。大多方法，仍然無法擺脫前處理的“魔咒”，即在設(shè)計領(lǐng)域和分析領(lǐng)域有一道天然的鴻溝，設(shè)計師做出來的模型無法直接交給工程師仿真分析，必須通過前處理打通CAD-CAE之間的橋梁。


（有限元法（左）和無網(wǎng)格法（右）的模型）

于是，基于打通從CAD領(lǐng)域和CAE領(lǐng)域、使得模型“做出即可用”的想法——設(shè)計師畫好的CAD模型可以交給工程師直接仿真分析，無需太多的前處理步驟，且具有更高的分析精度、比有限元更接近原始模型的逼近度，2005年，美國德克薩斯大學(xué)計算與應(yīng)用數(shù)學(xué)系主任、航空航天工程與工程力學(xué)教授Thomas J.R.Hughes，提出了等幾何分析（IGA，Isogeometric Analysis，又有人叫同幾何分析）方法。


（Thomas J.R.Hughes老爺子（左）和他出版的等幾何分析經(jīng)典教材（中）。右圖是2018年8月大連理工大學(xué)舉辦的International Symposium on Isogeometric Analysis and Mesh Generation大會上，老爺子精神矍鑠的視頻連線，祝賀大會召開）

03  真容

IGA是打通CAD和CAE的橋梁，其提出是兩個領(lǐng)域理論成果的結(jié)合。自然的，除了CAE方面的知識（上文介紹的有限元），還需了解CAD領(lǐng)域的內(nèi)容，為此先介紹下文的NURBS相關(guān)知識。

非均勻有理B樣條（Non-uniform Rational B-spline，NURBS），是基于樣條曲線（B樣條）的有理形式，在大多數(shù)工業(yè)產(chǎn)品數(shù)據(jù)交換標準（例如IGES、STEP）中，被作為定義幾何形狀的主要甚至唯一數(shù)學(xué)方法。什么意思呢？看下圖的B樣條曲線：


這條曲線的數(shù)學(xué)形式可以寫成參數(shù)表示：




這里需要明確幾個概念：

定義域：即參數(shù)u所在的自變量區(qū)間U∈[0,1]，任意一條參數(shù)曲線，都可以定義在[0,1]上，對應(yīng)的都能計算得到曲線上一點坐標，你可以理解為當自變量u從0到1滑過每個位置，相應(yīng)的就畫出一條曲線。

節(jié)點：定義域[0,1]之間劃分成不等距的子區(qū)間，每個自變量參數(shù)ui稱為節(jié)點。

基函數(shù)：即公式中的N。

控制頂點：即公式中的P，也就是圖中的包絡(luò)曲線的P0到P6的控制點。

階次：即基函數(shù)N的下標p。曲線/曲面都是有階次的，例如：最常見的兩點線段是線性，也就是一階；拋物線以及這里的圖示是二階曲線。曲線階次越高，表現(xiàn)的越光順，可導(dǎo)性越高（一階導(dǎo)為曲線切線，二階導(dǎo)為曲線曲率）。但是由于人眼可見性，曲線/曲面超過三階次，人眼很難分辨出來，可是現(xiàn)實中看到階次更高的形體，會不自覺的感到很舒服，高階曲面在特定的光影環(huán)境中也表現(xiàn)的更為炫麗。飛機高速飛行需要絕對光順的表面、汽車的造型設(shè)計和曲面設(shè)計，其實都是圍繞這方面做文章。


（飛機和汽車的光順曲面，極致優(yōu)美）

具體的計算機圖形學(xué)知識不需要過多了解（想深入了解的參考下方兩本教材），你只需要知道，任意一條復(fù)雜曲線或一張復(fù)雜曲面，都可以用上述公式這種“基函數(shù)×控制頂點坐標”的求和形式描述出來。我們在各種造型軟件（CATIA/UG等）看到的所有模型，其底層都是這種數(shù)學(xué)描述形式，通過一定的渲染和求交技術(shù)，最終展現(xiàn)在我們眼前的是可視化的點、線、面、體。


（CATIA的創(chuàng)成式設(shè)計或自由造型模塊，其中的樣條線創(chuàng)建按鈕，即為創(chuàng)建樣條曲線。此外，做好的數(shù)據(jù)保存成stp或igs不同格式，其實遵循的就是不同的工業(yè)數(shù)據(jù)交換標準STEP或IGES）

（計算機圖形學(xué)經(jīng)典教材：皮格爾著《The NURBS Book》、施法中著《CAGD&NRBS》。CAD領(lǐng)域的發(fā)展史，更是一部集數(shù)學(xué)、數(shù)值方法、圖形學(xué)、幾何學(xué)、藝術(shù)學(xué)、心理學(xué)等泛領(lǐng)域的精彩史詩，期間涌現(xiàn)出大批牛人，諸如高斯、孔斯、德·卡斯特里奧、貝齊爾、德布爾、考克斯等學(xué)者，以及各大軟件研發(fā)廠商互相吞并、競爭，各種技術(shù)相互較量淘汰，最終才有了我們今天高度靈活的造型系統(tǒng)）



（筆者基于Visual Studio MFC開發(fā)的創(chuàng)建NURBS曲線程序，可以實現(xiàn)基于圖形區(qū)的人機界面自由交互，也能自由對曲線升階、降階、增減頂點等）

到這里，細心的讀者相信已經(jīng)發(fā)現(xiàn)了，描述曲線/曲面坐標的數(shù)學(xué)公式，和我們上文中介紹有限元位移場的插值公式，居然如此相似！


（左：有限元單元位移插值公式；右：B樣條參數(shù)表達式）

實際上，不只是公式表達形式相似，其所蘊含的物理意義也能一一對應(yīng)：



我們知道，有限元單元的特性（單元階次、積分點、積分形式），完全取決于其插值形函數(shù)，即不同的形函數(shù)對應(yīng)了不同的單元類型，你可以簡單的認為，一個有限元模型中任意點的位移，都是通過周邊節(jié)點的位移和形函數(shù)加權(quán)相乘湊出來的；同樣的，對NURBS曲線曲面來說，其特性（階次、可導(dǎo)性、局部性，仿射不變性）也完全取決于其插值基函數(shù)。

為什么有限元分析要先進行前處理劃分網(wǎng)格？其實就是把原本的CAD模型（NURBS參數(shù)公式）轉(zhuǎn)化成CAE模型（有限元位移插值公式）——比如一張光滑曲面片，在CATIA中使用NURBS公式表述，經(jīng)過網(wǎng)格劃分、單元構(gòu)造后，原本的光滑曲面被許多線性小單元來逼近擬合，其中每個小單元又可以用有限元位移公式表述，這就完成了從CAD語言到CAE語言的“翻譯”。






（天才的靈感火花是不可捉摸的，你永遠搞不懂多種多樣的基函數(shù)、形函數(shù)，如此簡單又復(fù)雜的公式，是被怎樣發(fā)現(xiàn)并推導(dǎo)出來的）

所以，如果能直接用描述CAD模型的數(shù)學(xué)公式作為CAE模型的公式，不就可以省去“翻譯”的過程，幫助工程師從前處理的泥潭中解放出來了嘛？IGA，便是直接將NURBS曲線曲面表達式中的基函數(shù)拿來，作為有限元位移場的插值形函數(shù)，使得CAD模型和CAE模型具有統(tǒng)一的數(shù)學(xué)表達形式。

04  發(fā)展

可能有人還是不明白，IGA直接利用NURBS基函數(shù)作為位移形函數(shù)，這有什么用呢？

首先，省去網(wǎng)格劃分：不同于FEM每次要將CAD模型離散成很小的網(wǎng)格單元，IGA將原本的復(fù)雜曲面分割成一塊塊的不那么復(fù)雜的子級曲面片（patch）。本質(zhì)上IGA還是基于有限元的“分割”思想，但是，相比較有限元動輒千萬個幾毫米的線性網(wǎng)格，IGA可以用很少很大的幾個patch來逼近原模型，前處理的工作量大大減少。更重要的是，前處理的減負不僅僅體現(xiàn)在某一次的模型分析中，須知實際工程中的模型迭代更新是很頻繁的，以汽車為例，往往從預(yù)研到SOP期間要經(jīng)過3到5輪分析，每次分析都意味著有限元網(wǎng)格的重新劃分和建模，而使用IGA方法，等CAD模型更新后，其patch也會同步更新（因為patch用的形函數(shù)就是NURBS基函數(shù)），無需或者簡單調(diào)整模型即可直接進行分析。


（NURBS子級曲面片，圖引自Piegl L A, Tiller W. Geometry-based triangulation of trimmed NURBS surfaces）

其次，幾何精確性：由于patch使用了和原模型相同的數(shù)學(xué)表達式（NURBS基函數(shù)），幾乎可以說是原模型的“分裂體”，逼近程度非常高；patch之間的連接也能保證一定的光順性，無需擔(dān)心有限元單元之間出現(xiàn)的應(yīng)力不連續(xù)、單元畸變，無需在后處理時特意處理單元之間的連接問題。


（CAD曲面（左）、FEM逼近曲面（中）、IGA逼近曲面（右））

最后，結(jié)果高精度：我們知道有限元要想提高結(jié)果的精度，勢必要加密網(wǎng)格（甚至重新劃分網(wǎng)格），導(dǎo)致計算時長增加。而采用IGA單元，由于NURBS基函數(shù)具有天然的局部性和遞推性，低階次的單元可以很方便的構(gòu)造高階單元（通過插入節(jié)點、基函數(shù)升階或者兩者混合使用的方式），這就意味著IGA提高精度無需重新分割成子級patch，無需改變單元本身的形狀，這也是“等幾何分析”名稱的由來。相比較有限元加密網(wǎng)格得到的的h-收斂，IGA通過采用高階單元得到的p-收斂具有更高的精度。

Hughes最初提出IGA是以NURBS基函數(shù)作為有限元形函數(shù)，而Bazilevs則證明了IGA具備收斂性和穩(wěn)定性的特征，這實際上論證了IGA應(yīng)用于工程實際的可能，從此打開了IGA蓬勃發(fā)展的大門。直到今天的短短十幾年間，針對IGA的研究成果層出不窮，包括對IGA魯棒性研究、構(gòu)造更高精度的單元、應(yīng)用于相域建模、將IGA和有限元或無網(wǎng)格法耦合分析等。此外，IGA誕生于NURBS的參數(shù)表達式，但不同于幾何形狀的顯式或隱式表示，參數(shù)曲線曲面并不適合進行復(fù)雜的拓撲運算，例如曲面造型中最常見也最難的“求交裁剪”問題：


（兩個參數(shù)曲面求交，很容易導(dǎo)致交點漏求、多求、求錯的情況）


（曲面求交表面上看令人滿意，然而放大局部會發(fā)現(xiàn)有縫隙（gap）或重疊（overlay）這種現(xiàn)實物體不可能發(fā)生的問題，這會導(dǎo)致難以構(gòu)造IGA單元）

外行人乍一看覺得，不就是兩張曲面交在一起，求出它們的交線不就行了嘛？問題在于，怎么求，或者說怎么又快又準的求出。當前IGA僅能對比較簡單的形狀進行模擬，對于復(fù)雜型面，尤其是通過曲面求交得到的裁剪曲面，要構(gòu)造IGA單元并不是那么容易。



（曲面求交涉及的知識和要解決的問題非常繁雜，圖示所列的是基于包圍盒和碰撞檢測算法求交，其中的碰撞檢測類似于LS-DYNA中的接觸算法，都是現(xiàn)在大火的人工智能、圖像識別技術(shù)常用的方法）


（為了解決傳統(tǒng)IGA的B樣條基函數(shù)出現(xiàn)縫隙問題，有學(xué)者提出利用T樣條代替B樣條作為基函數(shù)，圖引自Sederberg T W, Zheng J, Bakenov A, et al. T-splines and T-NURCCs）

由于NURBS基函數(shù)是有理多項式，不像FEM形函數(shù)那樣滿足正交條件和插值特性，因此其邊界條件施加遠不如FEM方便（FEM只需對網(wǎng)格節(jié)點施加位移邊界添加即可）。如果模型的邊界是通過求交裁剪得到的不規(guī)則曲線，則更為困難。除此之外，盡管IGA單元具有數(shù)量少的優(yōu)勢，但由于其高階特性，同等規(guī)模下IGA的計算效率并不比FEM高出多少，而工程上FEM的計算精度也完全滿足要求，沒必要刻意追求更高精度，因此IGA在實際工程應(yīng)用上還有一段路要走。


（裁剪邊界的IGA單元構(gòu)造）

無論如何，IGA自誕生以來，僅僅經(jīng)過十幾年時間的發(fā)展，便大有取代FEM之勢，盡管當前工程上的主流分析方法仍然是FEM。隨著IGA相關(guān)理論的完善和諸多問題（如裁剪曲面的單元構(gòu)造、邊界條件施加、計算效率問題等）的解決，相信會有越來越多的軟件樂意將這一新的分析方法集成進去（當前ANSYS、Altair、DYNA已經(jīng)將IGA集成到軟件中了），并在更多領(lǐng)域開展應(yīng)用。未來的CAE仿真領(lǐng)域，IGA勢必能和FEM分一杯羹，或兩者耦合分析發(fā)揮最大效能——畢竟，不用前處理就能直接分析的模型，有哪個工程師不愛呢？

05  題外話：國產(chǎn)CAE軟件開發(fā)之路

無論是IGA方法還是FEM方法，現(xiàn)如今都被集成到各種各樣工業(yè)軟件中廣泛應(yīng)用于航空航天、船舶工程、汽車研發(fā)制造、家具家電、醫(yī)藥器材等工程領(lǐng)域，然而從事相關(guān)行業(yè)的工程師或?qū)W者會發(fā)現(xiàn)，當前主流的工業(yè)軟件幾乎清一色被國外公司壟斷（DYNA/ANSYS/ABAQUS/NASTRAN/ADINA等），筆者曾經(jīng)做過一定的國產(chǎn)CAD/CAE軟件研發(fā)工作，結(jié)合自己的研發(fā)體會，在本文最后介紹下國產(chǎn)CAE軟件的坎坷歷程。

第一，很多人有個錯覺，認為工業(yè)軟件和我們生活中常用的app（微信、美團、王者榮耀等）并無二致，這好比是把一輛汽車的含量和一件襯衫相比。工業(yè)軟件的本質(zhì)和內(nèi)核，從來就不屬于IT產(chǎn)品，舉個例子，現(xiàn)在比較熱門的機器學(xué)習(xí)、圖像識別、語音識別、大數(shù)據(jù)，其內(nèi)核全部屬于工業(yè)軟件。工業(yè)軟件最終呈現(xiàn)在我們面前的只是一個經(jīng)過渲染和圖形處理的集成平臺的前端接口，它底層所含有的架構(gòu)、數(shù)據(jù)結(jié)構(gòu)、算法、邏輯，凝聚了各行各業(yè)最前沿的基礎(chǔ)理論研究成果，包括但不限于計算機圖形學(xué)、幾何學(xué)、固體力學(xué)、流體力學(xué)、電磁學(xué)、熱力學(xué)、聲學(xué)、振動學(xué)、光學(xué)……，這意味著一個工業(yè)軟件的開發(fā)，是背后各行各業(yè)團隊高度合作的成果，當各領(lǐng)域最前沿的知識固化到一行行代碼中，才能真正體現(xiàn)其應(yīng)用價值、真正實現(xiàn)知識的代際傳承。世界上最大的軟件公司不是微軟、不是谷歌、不是蘋果，而是世界上最大的軍火商洛克希德公司，多少年來其一直致力于工業(yè)軟件的研發(fā)和更新

第二，我國工業(yè)軟件研發(fā)之路充滿了坎坷和無奈。初始，我們也是從模仿國外軟件研發(fā)開始起步（最開始主要集中在二維CAD制圖領(lǐng)域）。七五規(guī)劃到十五規(guī)劃的20年間，在國家機械部、國家科委的扶持下，國產(chǎn)工業(yè)軟件（CAD/CAE）開發(fā)迎來了一波蜜月期，期間涌現(xiàn)出一批以唐榮錫、鐘萬勰院士等為代表的先鋒人物，并開發(fā)出以北航（熊貓CAD）、大連理工（JIFEX）、中科院（飛箭）、中航工業(yè)（HAJIF）等為代表的一批優(yōu)秀國產(chǎn)軟件；然而，十二五規(guī)劃之后，中國開始兩化融合，科技部不再分管這方面工作，而工信部認為工業(yè)軟件研發(fā)并不屬于基礎(chǔ)科研領(lǐng)域，不予資助，工業(yè)軟件研發(fā)從此變成了“爹不疼娘不愛”的野孩子。什么概念呢？以全球最大的仿真軟件開發(fā)商ANSYS為例，據(jù)統(tǒng)計，其每年在軟件研發(fā)投入上達20億人民幣；而我國在十五規(guī)劃到十二五規(guī)劃的15年期間，對工業(yè)軟件研發(fā)的總投入不足2億人民幣。這意味著國外一家廠商1年的研發(fā)費用投入吊打我們整個國家15年總投入的10倍，更遑論研發(fā)團隊和高精尖人才的投入。

第三，工業(yè)軟件的研發(fā)只是起步工作，更重要的是后續(xù)一系列更新迭代。工業(yè)軟件誕生的一般流程是“工程需求——理論研究——軟件研發(fā)——工程試用和反饋——商業(yè)化——市場推廣——更新迭代”。軟件開發(fā)難度自然不言而喻，但若沒有一個良好的生態(tài)循環(huán)保證，其最終的下場就是被淘汰或并購。近十幾年來，國外各大軟件供應(yīng)商針對基礎(chǔ)理論和技術(shù)的更新拓展競爭達到了水深火熱的地步，各大廠商瘋狂并購、重組，其間每一次的技術(shù)整合都帶來軟件的進化——無論是人機交互還是底層算法。國外軟件起步早，早期打入我國市場后并不注重盜版打擊，一方面是由于我國知識產(chǎn)權(quán)保護法不完善，另一方面有一定的刻意成分，因為借助正版或盜版軟件的廣泛傳播使用，廠商不但省去了宣傳推廣成本，還能及時收集軟件使用中的問題，最終的目的只有一個——使自家軟件不斷迭代更新、不斷提升產(chǎn)品競爭力。相較而言，我國普遍在整個生命周期的第一步“工程需求”就折戟沉沙，試想一下：作為企業(yè)人員，與其自己花幾十年搞一套程序出來（搞出來的程序還未必比得過人家，還要后面花幾十年的時間進行維護），花費幾十萬購買國外軟件的現(xiàn)成使用權(quán)他不香嗎？作為高校院所，鐵打的支撐流水的學(xué)者，與其花費大量精力搞一套國產(chǎn)自主軟件，直接用破解版的國外軟件算出自己想要的結(jié)果，不是更快更準？畢竟，論文、畢業(yè)、評職稱才是現(xiàn)在的王道啊。

總之，我國工業(yè)軟件遲遲發(fā)展不起來，其中固然和國外軟件以摧枯拉朽的勢頭長驅(qū)直入有關(guān)，然而國家不重點扶持、缺少人才和資金的投入，軟件本身研發(fā)難度大、周期長、收益慢，當前資本市場和實際應(yīng)用無價值導(dǎo)向，這其中一系列原因，都是制約國產(chǎn)軟件的瓶頸。工業(yè)軟件在現(xiàn)代工業(yè)的重要地位不言而喻——我們總把汽車比喻為現(xiàn)代工業(yè)的皇冠，將發(fā)動機和芯片比喻為皇冠上的明珠，然而，工業(yè)軟件才是實現(xiàn)皇冠和明珠從無到有的關(guān)鍵工具。

德國率先喊出智能制造和工業(yè)4.0，引起全世界國家的高潮，各大企業(yè)紛紛高喊數(shù)字化轉(zhuǎn)型，讓軟件貫穿設(shè)計、研發(fā)、制造、管理、運營、維護、反饋各領(lǐng)域，相應(yīng)也帶動了資本市場關(guān)于半導(dǎo)體、人工智能、芯片等行業(yè)的一波紅利。卻沒人愿意承認，我國至今都沒有一款真正達到世界先進水平的工業(yè)軟件，我們針對數(shù)字化和智能制造浪潮所做的所有努力，都必須基于國外軟件這一基本工具；誰又能想到，一旦從工業(yè)軟件上卡脖子，一旦停止對工業(yè)軟件的使用權(quán)，以上所有的行業(yè)前景不過是紙上談兵而已。

06  參考文獻

1.Breitenberger M, Apostolatos A, Philipp B, et al. Analysis in computer aided design: Nonlinear isogeometric B-Rep analysis of shell structures[J]. Computer Methods in Applied Mechanics and Engineering, 2015, 284: 401-457.

2.Rank E, Ruess M, Kollmannsberger S, et al. Geometric modeling, isogeometric analysis and the finite cell method[J]. Computer Methods in Applied Mechanics and Engineering, 2012, 249-252: 104-115.

3.Marussig B, Hughes T J R. A Review of Trimming in Isogeometric Analysis: Challenges, Data Exchange and Simulation Aspects[J]. Archives of Computational Methods in Engineering, 2017.

4.Kim H J, Seo Y D, Youn S K. Isogeometric analysis with trimming technique for problems of arbitrary complex topology[J]. Computer Methods in Applied Mechanics & Engineering, 2010, 199(45): 2796-2812.

5.Kim H J, Seo Y D, Youn S K. Isogeometric analysis for trimmed CAD surfaces[J]. Computer Methods in Applied Mechanics & Engineering, 2009, 198(37): 2982-2995. Piegl L A, Tiller W. Geometry-based triangulation of trimmed NURBS surfaces[J]. Computer-Aided Design, 1998, 30(1): 11-18.

6.Sederberg T W, Zheng J, Bakenov A, et al. T-splines and T-NURCCs

7. Lespiegl, Waynetiller. 非均勻有理B樣條[M]. 清華大學(xué)出版社, 2010.

8. 施法中. 計算機輔助幾何設(shè)計與非均勻有理B樣條[M]. 高等教育出版社, 2001.

9. 朱心雄等. 自由曲線曲面造型技術(shù)[M]. 科學(xué)出版社, 2000.

10. 孫家廣. 計算機圖形學(xué)(第三版)[M]. 清華大學(xué)出版社, 1998.

11. 邢譽峰. 計算固體力學(xué)原理與方法[M].

12. 林雪萍. 工業(yè)軟件黎明靜悄悄.

13. 蘇楠. 工業(yè)軟件是制造強國之重器.