微分方程數(shù)值計算方法簡介
求解微分方程的常用數(shù)值方法有有限差分法、有限元法、有限有限體積法、邊界元法、離散元法、光滑粒子流體動力學(xué)法等。
1.有限差分法
有限差分法(FDM)是一種直接將微分問題化為代數(shù)問題的數(shù)值方法。有限差分法具有簡單、靈活以及通用性強(qiáng)的特點,容易在計算機(jī)上實現(xiàn),是發(fā)展比較早且比較成熟的數(shù)值方法。
有限差分法的基本思想是先把問題的定義域進(jìn)行網(wǎng)格劃分,然后在網(wǎng)格點上,按適當(dāng)?shù)臄?shù)值微分公式把定解問題中的微商換成差商,從而把表示變量連續(xù)變化關(guān)系的偏微分方程離散成為有限個線性代數(shù)方程,進(jìn)而求出數(shù)值解。此外,還要研究差分格式的解的存在性和唯一性、解的求法、解的數(shù)值穩(wěn)定性、差分格式的解與原定解問題的真解的誤差估計、當(dāng)網(wǎng)格尺寸趨于零時差分格式的解是否趨于真解等。
目前主要采用泰勒級數(shù)展開方法來構(gòu)造差分格式,有一階向前差分、一階向后差分、一階中心差分和二階中心差分等,其中前兩種差分格式為一階計算精度,后兩種格式為二階計算精度。
有限差分法已經(jīng)成為求解各類數(shù)學(xué)物理問題的主要數(shù)值方法之一。在固體力學(xué)中,有限元法出現(xiàn)以前,主要采用差分方法;?在流體力學(xué)中,差分方法仍是較為常用的數(shù)值方法。代表性的商業(yè)軟件有爆炸沖擊分析軟件?AUTODYN、巖土工程專業(yè)分析軟件?FLAC。有限差分法適于求解邊界較為規(guī)則的問題,如果邊界條件比較復(fù)雜,就不太適用了。
2.有限元法
有限元法(FEM)是一種求解偏微分方程邊值問題近似解的數(shù)值技術(shù),求解時對整個問題域進(jìn)行分解,每個子區(qū)域都成為簡單的部分。
有限元法通過變分法,使得誤差函數(shù)達(dá)到最小值并產(chǎn)生穩(wěn)定解。其基本思想是將連續(xù)的求解區(qū)域分解成一組有限個小的互聯(lián)單元,對每一單元假定一個合適的(較簡單的)近解,然后推導(dǎo)求解這個域總的滿足條件,從而得到問題的解。有限元法的實質(zhì)是最小勢能原理的應(yīng)用,對每一有限元假定一個合適的、較簡單的近似函數(shù)來分片地表示全求解域上待求的未知函數(shù),單元內(nèi)的近似函數(shù)通常由未知場函數(shù)或其導(dǎo)數(shù)在各個節(jié)點的數(shù)值及其插值函數(shù)來表達(dá)。如此一來,在一個問題的有限元分析中,未知場函數(shù)或其導(dǎo)數(shù)在各個節(jié)點上的數(shù)值就成為新的未知量(即自由度),從而使一個連續(xù)的無限自由度問題變成離散的有限自由度問題。一經(jīng)求出這些未知量,就可以通過插值函數(shù)計算出各個單元內(nèi)場函數(shù)的近似解,從而得到整個求解域上的近似解,達(dá)到用較簡單的問題代替復(fù)雜的實際問題的目的。
有限元法的解不是準(zhǔn)確解,而是近似解,因為實際問題被較簡單的問題所代替。而且.隨著單元數(shù)目的增加,即單元尺寸的減小,或隨著單元自由度的增加及插值函數(shù)精度的提高,解的近似程度將不斷改進(jìn)。如果單元是滿足收斂要求的,近似解將收斂于精確解。由于大多數(shù)實際問題難以得到準(zhǔn)確解,而有限元法不僅計算精度高,而且能適應(yīng)各種復(fù)雜形狀因此有限元法成為一種廣泛應(yīng)用、行之有效的數(shù)值分析手段。
有限元法的缺點是難以應(yīng)用于具有極大單元變形的情況,求解大型問題時需要的內(nèi)存和計算量比其他數(shù)值算法要大得多、預(yù)估計算產(chǎn)生的誤差比較困難。
采用有限元法編制的軟件非常多,代表性的商業(yè)軟件有?ANSYS Mechanical、LS-DYNA、MSC.Nastran?等。
3.?有限體積法
有限體積法(finite volume method,FVM)又稱有限容積法、控制體積法,是計算流體力學(xué)中常用的一種數(shù)值計算方法,有限體積法基于積分形式的守恒方程而不是微分方程,該積分形式的守恒方程描述的是計算網(wǎng)格定義的每個控制體。有限體積法著重從物理觀點來構(gòu)造離散方程,每一個離散方程都是較小體積上某種物理量守恒的表達(dá)式,推導(dǎo)過程物理概念清晰,離散方程系數(shù)具有一定的物理意義,并可保證離散方程具有守恒特性。
有限體積法的基本思路是:將計算區(qū)域劃分為一系列不重疊的控制體積,并使每個網(wǎng)格點周圍有一個控制體積,將待解的微分方程對每個控制體積積分,便得出一組離散方程其中的未知數(shù)是網(wǎng)格點上的因變量的數(shù)值。為了求出控制體積的積分,必須假定數(shù)值在網(wǎng)格點之間的變化規(guī)律,即假設(shè)值的分段的分布剖面。
有限體積法的優(yōu)點是:
(1)?具有很好的守恒性
(2)?更加靈活的假設(shè),可以克服泰勒展開離散的缺點。
(3)對網(wǎng)格的適應(yīng)性好,可以很好地解決復(fù)雜的工程問題
(4)在進(jìn)行流固耦合分析時,能夠和有限元法完美融合。
大多數(shù)計算流體力學(xué)軟件如?Fluent,STAR-CD、Flotherm?等采用的都是有限體積法此外有限體積法也是爆炸沖擊計算軟件?MSC.Dytran?的基本算法之一。
4 .邊界元法
邊界元法(boundary element method,BEM)是繼有限差分法有限元法之后發(fā)展起來的又一種重要的數(shù)值方法。邊界元法只在定義域的邊界上劃分單元,用滿足控制方程的函數(shù)去逼近邊界條件。與有限元法相比,邊界元法降低了問題的維數(shù)和自由度數(shù),具有單元數(shù)量少、數(shù)據(jù)準(zhǔn)備簡單、計算時間少等優(yōu)點,已經(jīng)成為現(xiàn)代科學(xué)和工程數(shù)值分析的有效工具。邊界元法的基本思想是應(yīng)用格林定理等,將問題的控制方程轉(zhuǎn)換成邊界上的積分方程然后引人位于邊界上的有限個單元將積分方程離散求解。借鑒有限元法劃分單元的離散技術(shù),通過對表面邊界進(jìn)行離散,得到邊界單元,將邊界積分方程離散成線性方程、經(jīng)過離散后的方程組只含有邊界上的節(jié)點未知量,因而降低了問題的維數(shù),最后求解方程的階數(shù)降低,數(shù)據(jù)準(zhǔn)備方便,計算時間縮短。另外,邊界元法通過引用問題的基本解而具有解析與離散相結(jié)合的特點,使得計算精度較高。由于積分方程可以用加權(quán)余量法得到,這就避免了尋找泛函的麻煩。
邊界元法的主要缺點是它的應(yīng)用范圍以存在相應(yīng)微分算子的基本解為前提,對于非均勻介質(zhì)等問題難以應(yīng)用,故其適用范圍遠(yuǎn)不如有限元法廣泛,而且通常由它建立的求解代數(shù)方程組的系數(shù)矩陣具有非對稱、稠密,甚至在某些情況下病態(tài)的特性。應(yīng)用邊界元法分析大規(guī)模問題,系統(tǒng)方程求解過程消耗大量計算資源和時間,對解題規(guī)模產(chǎn)生較大限制。對一般的非線性問題,由于在方程中會出現(xiàn)域內(nèi)積分項,從而部分抵消了邊界元法只離散邊界的優(yōu)點。
邊界元法代表性的商業(yè)軟件有三維邊界元流體力學(xué)分析軟件?ANSYS LINFLOW等在?LS-DYNA?軟件中也包含了邊界元聲學(xué)計算功能。
5.?離散元法
離散法(discrete element method?或者?distinct element method,DEM)是?20?世紀(jì)70?年代由?Cundall?首先提出來的,起源于分子動力學(xué),是為研究巖體等非連續(xù)介質(zhì)的力學(xué)行為而發(fā)展起來的一種數(shù)值方法。
離散元法的基本原理是牛頓第二定律,其基本思想是把求解域分離為離散單元的集合使每個離散單元滿足牛頓第二定律,離散單元本身一般為剛體,單元間的相對位移等變形行為一般由聯(lián)接于節(jié)點間的變形元件來實現(xiàn)。離散單元可以平移、轉(zhuǎn)動或變形。各個離散單元在外界的干擾下就會產(chǎn)生力和力矩的作用,由牛頓第二定律可以得到各個離散單元的加速度,然后對時間進(jìn)行積分,就可以依次求出離散單元的速度、位移,最后得到離散單元的變形量。離散單元在位移向量的方向會發(fā)生調(diào)整,這樣又會產(chǎn)生力和力矩的作用。如此循環(huán)直到所有離散單元達(dá)到一種平衡狀態(tài)或者處于某種運動狀態(tài)之下,使各個離散單元滿足運動方程,用時間步長迭代的方法求解各個離散單元的運動方程,繼而求得不連續(xù)體的整體運動形態(tài)。離散單元法允許單元間的相對運動,不一定要滿足位移連續(xù)和變形協(xié)調(diào)條件,計算速度快,所需存儲空間小,尤其適合求解大位移和非線性的問題。
離散元法比較適合于模擬節(jié)理系統(tǒng)或者離散顆粒組合體在準(zhǔn)靜態(tài)或動態(tài)下的變形過程,如顆粒材料的混合、分離、注裝、倉儲和運輸過程,以及連續(xù)體結(jié)構(gòu)在準(zhǔn)靜態(tài)或動態(tài)條件下的變形及破壞過程。
目前,離散元的理論體系并不完善,其局限性表現(xiàn)在接觸模型和參數(shù)難以確定、穩(wěn)定性難以保證等方面。
基于離散元法的商業(yè)軟件有?UDEC、PFC、GRANULE、EDEM?等,離散元也是?LSDYNA?軟件的重要功能之一。
6 .光滑粒子流體動力學(xué)法
光滑粒子流體動力學(xué)(smoothed particle hydrodynamics,SPH)法是一種為解決天體物理中涉及的流體質(zhì)團(tuán)在三維空間無邊界情況下的任意流動問題的純拉格朗日方法。后來發(fā)現(xiàn)解決如連續(xù)體結(jié)構(gòu)的解體、碎裂、固體的層裂、脆性斷裂等物理問題,以及產(chǎn)生大變形的流體動力學(xué)問題、考慮材料強(qiáng)度的固體動力學(xué)問題也非常有效??傊?SPH?法適合用來求解具有大變形、復(fù)雜邊界和物質(zhì)交界面等復(fù)雜的單相或多相流體動力學(xué)問題,是最早的無網(wǎng)格方法。
光滑粒子流體動力學(xué)法的理論基礎(chǔ)是插值理論,通過一種稱之為核函數(shù)的積分核進(jìn)行插值近似.從而將流體力學(xué)方程轉(zhuǎn)換為數(shù)值計算用的?SPH方程組。光滑粒子流體動力學(xué)法的基本思想是用一系列的粒子來表示求解域,這些粒子具有流體的密度、速度、熱能等物理量,粒子之間不需要任何連接,即具有無網(wǎng)格特性,粒子的密度、位移、速度、壓力等物理量的更新只與時間有關(guān):用積分近似和粒子近似離散流體動力學(xué)控制方程,產(chǎn)生帶狀或離散化的稀疏系數(shù)矩陣,某一粒子場函數(shù)的值通過支持域內(nèi)相鄰粒子的疊加求和計算得到;由于所使用的用于求和的局部粒子為當(dāng)前時間步的粒子,所以?SPH?法具有較強(qiáng)的適應(yīng)性;將粒子近似法應(yīng)用于所有偏徽分方程的場函數(shù)相關(guān)項中,可得到一系列只與時間相關(guān)的離散化形式的富徽分方程,應(yīng)用顯式積分法來求解常微分方程以獲得最快的時間積分,并可得到所有粒子的場變量隨時間的變化值。
SPH法相對于傳統(tǒng)的基于網(wǎng)格的方法具有一些特別的優(yōu)點:具有自適應(yīng)性、無網(wǎng)格性以及拉格朗日公式與粒子近似法的結(jié)合。一方面,由于?SPH?法的近似過程不受網(wǎng)格的限制,因此可避免網(wǎng)格極度大變形造成的精度下降問題;另一方面,SPH?法是一種拉格朗日方法,每個粒子點代表一種獨立的物質(zhì),可以自然且直觀地跟蹤粒子所代表的物質(zhì)屬性,能夠方便地指述多物質(zhì)流動問題。
SPH法也存在一些數(shù)值方面的困難,如拉伸不穩(wěn)定性會引起數(shù)值斷裂、形函數(shù)不一致性、引入本質(zhì)邊界條件存在困難,需要比較復(fù)雜的接觸算法等。
基于SPH法的商業(yè)分析軟件有美國?CentroidLab?公司開發(fā)的?NEUTRINO,開源?SPH軟件有?DualSPHysics,SPHinXsys、SPHERAL?等,在?LS-DYNA、AUTODYN、ABAQUS軟件中均包含?SPH計算功能。
摘自:《Ls-dyna有限元建模、分析和優(yōu)化設(shè)計》辛春亮