有限元理論基礎(chǔ)及Abaqus內(nèi)部實(shí)現(xiàn)方式研究系列27： Abaqus內(nèi)部計(jì)算和顯示的應(yīng)變

SnowWave02

瀏覽：5732 評論：31 收藏：6

（原創(chuàng)，轉(zhuǎn)載請注明出處）

==概述==

有限元理論基礎(chǔ)及Abaqus內(nèi)部實(shí)現(xiàn)方式研究系列27： Abaqus內(nèi)部計(jì)算和顯示的應(yīng)變的圖3 本系列文章研究成熟的有限元理論基礎(chǔ)及在商用有限元軟件的實(shí)現(xiàn)方式，通過

（1）基礎(chǔ)理論

（2）商軟操作

（3）自編程序

三者結(jié)合的方式將復(fù)雜繁瑣的結(jié)構(gòu)有限元理論通過簡單直觀的方式展現(xiàn)出來，同時(shí)深層次的學(xué)習(xí)有限元理論和商業(yè)軟件的內(nèi)部實(shí)現(xiàn)原理。

有限元的理論發(fā)展了幾十年已經(jīng)相當(dāng)成熟，商用有限元軟件同樣也是采用這些成熟的有限元理論，只是在實(shí)際應(yīng)用過程中，商用CAE軟件在傳統(tǒng)的理論基礎(chǔ)上會做相應(yīng)的修正以解決工程中遇到的不同問題，且各家軟件的修正方法都不一樣，每個(gè)主流商用軟件手冊中都會注明各個(gè)單元的理論采用了哪種理論公式，但都只是提一下用什么方法修正，很多沒有具體的實(shí)現(xiàn)公式。商用軟件對外就是一個(gè)黑盒子，除了開發(fā)人員，使用人員只能在黑盒子外猜測內(nèi)部實(shí)現(xiàn)方式。

有限元理論基礎(chǔ)及Abaqus內(nèi)部實(shí)現(xiàn)方式研究系列27： Abaqus內(nèi)部計(jì)算和顯示的應(yīng)變的圖4

一方面我們查閱各個(gè)主流商用軟件的理論手冊并通過進(jìn)行大量的資料查閱猜測內(nèi)部修正方法，另一方面我們自己編程實(shí)現(xiàn)結(jié)構(gòu)有限元求解器，通過自研求解器和商軟的結(jié)果比較來驗(yàn)證我們的猜測，如同管中窺豹一般來研究的修正方法，從而猜測商用有限元軟件的內(nèi)部計(jì)算方法。我們關(guān)注CAE中的結(jié)構(gòu)有限元，所以主要選擇了商用結(jié)構(gòu)有限元軟件中文檔相對較完備的Abaqus來研究內(nèi)部實(shí)現(xiàn)方式，同時(shí)對某些問題也會涉及其它的Nastran/Ansys等商軟。為了理解方便有很多問題在數(shù)學(xué)上其實(shí)并不嚴(yán)謹(jǐn)，同時(shí)由于水平有限可能有許多的理論錯(cuò)誤，歡迎交流討論，也期待有更多的合作機(jī)會。

自主結(jié)構(gòu)有限元求解器iSolver介紹視頻：

有限元理論基礎(chǔ)及Abaqus內(nèi)部實(shí)現(xiàn)方式研究系列27： Abaqus內(nèi)部計(jì)算和顯示的應(yīng)變的圖5

http://www.yqgqt.org.cn/college/video/c12884

==第27篇：Abaqus內(nèi)部計(jì)算和顯示的應(yīng)變 ==

應(yīng)變度量在線性情況下沒有差異，但在幾何非線性情況下有不同的度量方式。系列文章18：幾何非線性的應(yīng)變中，我們提到了幾何非線性常用的三種應(yīng)變：工程應(yīng)變、Green應(yīng)變、對數(shù)應(yīng)變，其中工程應(yīng)變和Green應(yīng)變分別用于線性和小應(yīng)變情況，對數(shù)應(yīng)變用于幾何大應(yīng)變情況，那么在Abaqus的幾何大應(yīng)變中，是否真的就簡單的采用了對數(shù)應(yīng)變呢？Abaqus后處理顯示的應(yīng)變和計(jì)算的應(yīng)變是一樣的嗎？我們將舉殼單元和體單元兩個(gè)例子通過iSolver和Abaqus的結(jié)果對比來驗(yàn)證幾何大應(yīng)變Abauqs內(nèi)部計(jì)算用的應(yīng)變和后處理顯示的應(yīng)變。前面有朋友也問幾何非線性時(shí)iSolver計(jì)算用的應(yīng)變是什么？趁此文章我們也將證明iSolver計(jì)算時(shí)不同的單元用不同的應(yīng)變，每種應(yīng)變的選取方式和Abaqus一致。

1.1 對數(shù)應(yīng)變的計(jì)算方法

當(dāng)存在幾何大變形的情況，變形比較大時(shí)，此時(shí)終止時(shí)刻位移x-初始時(shí)刻X已經(jīng)和X比擬了。

對數(shù)應(yīng)變?nèi)某跏紩r(shí)刻到最終時(shí)刻點(diǎn)這段時(shí)間的累加，與材料點(diǎn)在整個(gè)時(shí)間段內(nèi)的變形路徑，是一種積分行為。

（1）在一維情況下，設(shè)拉伸率為r：

這個(gè)等式其實(shí)隱含了應(yīng)變是沿x方向的，如果在三維空間可以認(rèn)為有一個(gè)主方向n=（1,0,0）’，因?yàn)閼?yīng)變是二維張量，那么上式就是：

（2）在三維情況下，存在三個(gè)主方向n1,n2,n3，拉伸率分別為r1，r2，r3：

剩下的問題就是怎么求這三個(gè)主方向，這個(gè)主要是求一個(gè)與變形F相關(guān)的3X3矩陣的特征值和特征向量，具體可參考其它論文書，其實(shí)和Abaqus后處理中顯示應(yīng)變時(shí)有三個(gè)Principal的值是一樣的求法。

1.2 變形率積分的計(jì)算方法

由上面計(jì)算方法發(fā)現(xiàn)每次計(jì)算對數(shù)應(yīng)變都需要求一個(gè)特征值和特征向量，在數(shù)值計(jì)算中，特征問題的求解耗時(shí)相對較多，且計(jì)算相對復(fù)雜（一般人都是認(rèn)為計(jì)算復(fù)雜才采用別的應(yīng)變，個(gè)人不太認(rèn)可），而實(shí)際許多非線性材料中，都有這樣一個(gè)規(guī)律，就是彈性應(yīng)變都相對較小，譬如典型的鋼材料，楊氏模量為2.1e11Pa，屈服應(yīng)力為235Mpa，那么達(dá)到屈服時(shí)的應(yīng)變?yōu)?35e6/2.1e11=1e-3，同時(shí)，典型的應(yīng)力應(yīng)變本構(gòu)曲線如下圖，在塑性段譬如C點(diǎn)的彈性應(yīng)變和屈服應(yīng)變差異并不大。

因此，Abaqus中假定內(nèi)置的所有材料都滿足彈性應(yīng)變相對較小，此時(shí)，理論可以證明，對數(shù)應(yīng)變可以簡單的取為變形率D的積分：

上式無法直接得到數(shù)學(xué)表達(dá)式，但在有限元中，可以通過增量形式累加。

1.3 調(diào)試Abaqus內(nèi)部計(jì)算應(yīng)變的方法

由于對數(shù)應(yīng)變和試驗(yàn)最接近，因此Abaqus后處理中的E都是用對數(shù)應(yīng)變來顯示的，Abaqus為了進(jìn)一步提示對數(shù)應(yīng)變，直接在后處理中如果選了NLGeom=On，應(yīng)變的顯示從E變味了LE，但Abaqus幾何非線性實(shí)際計(jì)算的應(yīng)變并不完全一致。

Abaqus真正計(jì)算的變量度量可以通過它的子程序的輸入?yún)?shù)獲取。在Abaqus中，增量步即代表時(shí)刻點(diǎn)，可以查看增量點(diǎn)時(shí)刻的子程序輸入來猜測Abaqus的內(nèi)部量描述方式。UMAT子程序中，在材料本構(gòu)函數(shù)中要利用應(yīng)變增量和當(dāng)前應(yīng)力等物理量更新應(yīng)力，查看UMAT等子程序的接口：

可知其中STRAN和DSTRAN分別表示當(dāng)前增量步最后時(shí)刻的應(yīng)變?nèi)亢驮隽俊＞唧w的介紹也可參考下面視頻講解：

http://www.yqgqt.org.cn/college/video/c13034

變形率D在一維上代表對數(shù)應(yīng)變的導(dǎo)數(shù)，但三維上并不是對數(shù)應(yīng)變的導(dǎo)數(shù)，這是有很大區(qū)別的，同時(shí)，可以利用iSolver分別采用上述兩種應(yīng)變度量和Abaqus子程序接口的結(jié)果比對來確認(rèn)Abaqus計(jì)算的應(yīng)變是哪種度量。所以下面我們將找一個(gè)體單元和一個(gè)殼單元的例子來驗(yàn)證到底Abaqus計(jì)算和顯示的應(yīng)變是什么。

1.4 體單元的例子

1.4.1 算例介紹

體單元算例參數(shù)如下：

尺寸：5X1X0.1。

材料：Young’s Modulus 1e8， Poisson Ratio 0.3。

左側(cè)四個(gè)節(jié)點(diǎn)固支。

右側(cè)四個(gè)節(jié)點(diǎn)約束位移為5,1,1。

劃分為一個(gè)殼單元C3D8R。

幾何非線性開關(guān)NLGeom=On，且控制只迭代一次。

1.4.2 Abaqus的應(yīng)變

Abaqus中采用殼的UMAT子程序進(jìn)行計(jì)算。

（1）顯示應(yīng)變：Abaqus計(jì)算完畢后得到導(dǎo)入結(jié)果，在后處理中查看，應(yīng)變E11=0.6819,E22=5.621e-3如下：

（2）計(jì)算應(yīng)變：Abaqus中采用UMAT子程序，利用我們的子程序調(diào)試插件DUS調(diào)試UMAT，在Visual Studio中查看dStran的值，發(fā)現(xiàn)在計(jì)算完應(yīng)變后，進(jìn)入U(xiǎn)MAT時(shí)，顯示如下，可得E11=0.6667， E22=0：

1.4.3 iSolver的應(yīng)變

猜測體單元中，abaqus顯示的是對數(shù)應(yīng)變，而計(jì)算還是用的變形率。為了證明這點(diǎn)，在iSolver以同樣方式建模，只是iSolver中采用自帶材料進(jìn)行計(jì)算，材料參數(shù)和UMAT的輸入完全一致。

計(jì)算中，iSolver也采用變形率計(jì)算方式，得到的應(yīng)變顯示E11=0.6667， E22=0，如下，可發(fā)現(xiàn)和Abaqus完全一致。

1.5 殼單元的例子

1.5.1 算例介紹

殼單元算例參數(shù)如下：

尺寸：5X1，厚度0.1。

材料：Young’s Modulus 1e8， Poisson Ratio 0.3。

左側(cè)兩個(gè)節(jié)點(diǎn)固支。

右側(cè)兩個(gè)節(jié)點(diǎn)每個(gè)加集中力6.73128E+006，x方向。

劃分為一個(gè)殼單元S4R。

幾何非線性開關(guān)NLGeom=On，且控制只迭代一次。

1.5.2 Abaqus的應(yīng)變

Abaqus中采用殼的UMAT子程序進(jìn)行計(jì)算。

（1）顯示應(yīng)變：Abaqus計(jì)算完畢后得到導(dǎo)入結(jié)果，在后處理中查看，應(yīng)變E11=8.528e-1,E22=-5.173e-1如下：

（2）計(jì)算應(yīng)變：Abaqus中采用UMAT子程序，利用我們的子程序調(diào)試插件DUS調(diào)試UMAT，在Visual Studio中查看dStran的值，發(fā)現(xiàn)在計(jì)算完應(yīng)變后，進(jìn)入U(xiǎn)MAT時(shí)，E11=8.528e-1,E22=-5.173e-1，調(diào)試如下：

可以發(fā)現(xiàn)殼單元Abaqus的計(jì)算應(yīng)變和顯示應(yīng)變一樣，猜測都是對數(shù)應(yīng)變。

1.5.3 iSolver的應(yīng)變

iSolver中采用自帶材料進(jìn)行計(jì)算，材料參數(shù)和UMAT的輸入完全一致。

為了計(jì)算和Abaqus完全一致，iSolver也采用對數(shù)應(yīng)變計(jì)算方式，得到的應(yīng)變顯示如下，可發(fā)現(xiàn)和Abaqus完全一致。

==總結(jié)==

由上可以看到，在實(shí)際計(jì)算中，對體單元，Abaqus和iSolver都采用變形率積分方式來計(jì)算應(yīng)變，對殼單元，Abaqus和iSolver都采用對數(shù)應(yīng)變。一般理論書都認(rèn)為Abaqus是因?yàn)閷?shù)應(yīng)變計(jì)算復(fù)雜才采用別的應(yīng)變，但個(gè)人認(rèn)為應(yīng)該不是這個(gè)原因，因?yàn)锳baqus對體單元為了顯示對數(shù)應(yīng)變，依然重新計(jì)算了一遍，說明Abaqus體單元采用變形率是有其它原因的，具體什么原因我也沒研究清楚，歡迎探討。

如果有任何其它疑問或者項(xiàng)目合作意向，也歡迎聯(lián)系我們：

snowwave02 From www.yqgqt.org.cn

email: snowwave02@qq.com

以往的系列文章：