有限元學習必須知道：隱式與顯式有限元算法

2017-02-27 by:CAE仿真在線來源:互聯(lián)網(wǎng)

導語:對于有限元的學習,最最關鍵的其實不是本構方程,也不是屈服準則,而在于對求解算法的理解。本人根據(jù)長期學習經(jīng)驗,在美國作訪問學者學習經(jīng)歷,純屬興趣,在這里做一點膚淺的總結,完全原創(chuàng)。

隱式與顯式有限元最大的區(qū)別在于是否迭代,是否所有的物理量在同一時刻獲得。采用隱式迭代求解平衡方程(位移、速度和加速度)、而不管是否用隱式與顯式的方法(前向或者后向歐拉求解方法)求解本構方程(應力和應變)叫做隱式有限元;用顯式時間積分的方法求解叫做顯式有限元。

首先,對于本構方程的求解,通常分為前向和后向歐拉算法。對于后向歐拉算法求解彈塑性問題,所有的物理量(包括等效塑性應變增量、N+1迭代步的應變和應力以及相關依賴于solution的狀態(tài)變量)均是同時求解獲得,因為涉及到多個物理量,而通常情況下他們是相互依賴、相互成為函數(shù),所以必須通過牛頓迭代同時求解幾個方程(如采用試應力方程、屈服函數(shù)徑向返回算法(對于各向異性,也叫回映算法,最近點的投射算法)聯(lián)合求解等效塑性應變增量)。對于前向歐拉,直接由N時刻的應力和應變求出N+1時刻的應力和應變,無需迭代。

其次,對于平衡方程的求解,通常分為隱式和顯式有限元算法。對于隱式有限元算法,由應力平衡方程+邊界條件變分之后獲得的剛度方程KU=F,隱式求解必須引入雅可比矩陣(二次收斂、只影響計算速率、不影響數(shù)值精度;K又稱為雅可比),其是實時更新的,是N+1時刻的應力、應變以及狀態(tài)變量(如損傷內(nèi)變量)的函數(shù),隱式求解是很robust的,確保了計算精度,但是不足之處在于計算非常expensive,每次迭代都要計算K的逆矩陣,也容易產(chǎn)生數(shù)值收斂性問題,目前解決的方法有弧長法、粘性阻尼法等,個人認為粘性阻尼法效果最好。

對于顯示算法,采用時間積分,用t+1時刻的積分點應力、應變,獲得t+1時刻的節(jié)點位移,無需迭代求解,也不需要雅可比矩陣(應力對應變偏導數(shù));如果硬是要有,連續(xù)雅可比,基于本構模型而不是剛度方程推導近似的連續(xù)雅可比。對于顯示算法,單元高斯積分點應力、應變的求解可用前向或者后向歐拉方法,然后通過時間積分求取節(jié)點位移。本質(zhì)上,平衡方程中位移的迭代求解與本構方程中的應力、應變求解沒有關聯(lián),這點很容易造成誤解,很多時候?qū)⑶啊⒑髿W拉算法視為顯式和隱式的區(qū)別,大大錯誤。通常應用較廣的顯示算法紐馬克法、威爾遜-sita法,其中改變紐馬克法中的兩個參數(shù),可以實現(xiàn)隱式與顯式求解,其中alpha=0.5和beta=0是中心差分法(二階精度)。目前一個大的誤區(qū)認為只有顯示算法可以求解動力學問題,隱式只能求解準靜態(tài)問題(如低速沖擊),alpha=0.5和beta=0.25就是隱式,所有的物理量在t+1時刻同時求解,通常ABAQUS軟件中所說的隱式動力學求解采用了斯坦福大學Hilber、HUGHES院士(現(xiàn)在德克薩斯大學奧斯丁分校)和加州大學伯克利分校Taylor院士提出的無條件穩(wěn)定隱式差分算法,可以求解低速動力學問題,缺點是不適合含阻尼的求解、計算效率不高;alpha=0.5和beta=0時的紐馬克法更適合求解動力學問題,主要原因在于比隱式求解計算效率更高,不足之處在于其是條件穩(wěn)定,時間增量過大位移解容易震蕩,根本原因是差分算法的條件穩(wěn)定導致的,時間增量必須非常小(其值越大,一方面不穩(wěn)定、另一方面計算誤差也更大),其依賴于波速、彈性模量和最小單元網(wǎng)格尺寸,這是顯式算法計算最耗時的地方。相對于隱式算法,顯式算法的功能更強大,通常計算依賴于率的變形和應力,也可以求解穩(wěn)態(tài)問題,如alpha=0.5和beta=0時,對于剛度方程中引入阻尼矩陣后,叫做動態(tài)松弛法,可以解決靜力學問題的一些穩(wěn)態(tài)問題(如重力、預應力引起的初始應力)。此外,一些準靜態(tài)的剪切自鎖問題,本質(zhì)上有解,但是用牛頓法求解失效,中心差分引入質(zhì)量矩陣后,可以獲得正常的解。需要注意的是,時間積分算法通常采用Lumped集中對角質(zhì)量矩陣而不是一致質(zhì)量矩陣,以提高計算效率。總體來說,由于計算效率的問題,隱式時間積分算法ABAQUS-Standard特別適合于低速沖擊問題;對于高速沖擊問題,由于存在不連續(xù)非線性接觸的動響應過程,隱式算法解決不好,使用顯式時間算法ABAQUS-Explicit更好。此外,對于瞬態(tài)和穩(wěn)態(tài)熱傳導問題,半離散的拋物線方程,中心差分法可以較好獲得溫度分布。

對于依賴于率的粘塑性問題(對于本質(zhì)上的粘性材料),與彈塑性材料的根本區(qū)別在于,一般來說是一致性條件不滿足(排除彈塑性材料在高溫下的軟化問題,對于這種問題,屈服條件也可以滿足),即屈服條件不滿足,N+1時刻的物理量不用回映到N+1時刻的屈服面上,粘塑性模型成為過應力模型,顯示和隱式算法都可以求解。對于依賴于率的本構模型,其可解決模擬高速沖擊、爆炸、彈道射擊問題時存在的動態(tài)應變局部化問題(對于動態(tài)問題,平衡方程喪失雙曲線特性;對于靜態(tài)問題,平衡方程失去橢圓性),解決網(wǎng)格尺寸效應,其實質(zhì)上是引入了適當?shù)淖枘徇t滯效應。需要注意的是,對于大變形(又稱為有限變形)問題,Cauchy應力率和速度梯度(包括客觀和對稱的扭曲張量率D、不客觀和反對稱的spin旋轉(zhuǎn)張量W兩個部分)均是不客觀的,為解釋剛體旋轉(zhuǎn)(如純剪切變形就包含剛體旋轉(zhuǎn)),在共旋坐標系下面求解真實應力和應變,應力和應變積分求解的時候應首先求解客觀性的Jaumann應力率(相對于真實應力,空間坐標系),相對于Second Piola-Kirchhoff應力是Truesdell率(材料坐標系)。ABAQUS軟件對于大變形問題已經(jīng)做了旋轉(zhuǎn)。

對于一些耦合場問題,由于計算量非常大,同時要求解太多物理量,如熱-流-固耦合,要求解位移、壓力、溫度,采用純隱式算法或顯式算法基本不太獲得收斂或準確的解,這時候可采用混合的隱式與顯式有限元格式mixed implicit-explicit partitoning方法,將剛度矩陣和阻尼矩陣分成兩個部分,在同一區(qū)域采取不同算法,提高計算效率和精度、穩(wěn)定性和收斂性。

開放分享：優(yōu)質(zhì)有限元技術文章,助你自學成才