幾何非線性微分剛度的有限元分析

摘 要：本文對大位移幾何非線性問題作有限元分析，將簡單情形分析的基本結論推廣到一般情形，給出一般情形幾何非線性問題微分剛度的計算方法，并對微元體、桿單元及等截面梁單元的微分剛度進行計算，這種微分剛度計算對于大位移幾何非線性問題的有限元分析具有非常重要的理論價值及應用價值。

關鍵詞：有限元方法；大位移；廣義力；微分剛度

中圖分類號：TU32 文獻標志碼：A文章編號：1671-7953(2009)01-0103-03

Finite Element Analysis of Differential Rigidity

About Geometry Nonlinearity

HU Ke

（Guangdong Institute Of Science And Technology Zhuhai Guangdong 519090，China）

Abstract: The article has made a finite element analysis about large displacement geometry nonlinearity question，spreading simple analyzing result to an ordinary way，calculating common method has been givens out about geometry nonlinearity differential rigidity，as an example it has been set up differential rigidity of small element，bar element and constant cross section beam element，in theory and use it is of great value to that calculating differential rigidity about large displacement geometry nonlinearity question.

Key words: Finite element；Large displacement；Generailzed forces；Differential rigidity.

對幾何非線性問題的有限元分析，計算大位移引起的附加剛度效應，比較簡單的方法是微分剛度計算法，如圖1所示，剛性桿AB長為l，A端鉸鏈連接并有一個扭轉彈簧（扭轉剛度為K），B端有一水平荷載P。

圖1

當桿由AB位置轉過θ角到達AB′位置時，P力做功為W=PuB=Pl（cosθ-1），彈性勢能為V=1/2Kθ 2，以轉角θ為廣義坐標，相應于彈性勢能的廣義力為

Qθ=-Vθ=-Kθ（1）

而相應于P力做功的廣義力為

Qθ=Wθ=-Plsinθ（2）

若考慮小轉動，即θ為小量（sinθ≈θ），則有Qθ=-Plθ，由（1）式與（2）式可見，Pl與剛度K相當，這里Pl就是P力在小轉動過程中附加的剛度，稱為微分剛度Kd，通常線剛度K由（1）式得

K=-Qθθ

同樣，微分剛度由（2）式得

Kd=-Qθθ=-2Wθ2=Pl

將上述結論推廣到多自由度彈性系統，設想將一般彈性結構離散成具有n個自由度彈性系統，其上作用m個力Pi（i=1，2，…，m），沿作用力方向位移為ui，則m個力總功為

W=mi=iPiui

與第r個廣義坐標qr相應的廣義力為

Qr=Wqr=mi=iPiuiqr；（r=1，2，…，n）

式中假定Pi不隨位移而變化，但位移ui可以是廣義坐標的函數，第r個自由度對第s個自由度的微分剛度為

Kdrs=-Qrqs=-2Wqrqs；（r，s=1，2，…，n）（3）

這就是一般彈性體出現大位移時，在線性剛度陣上附加的微分剛度。顯然它是一個對稱矩陣，即Kdrs=Kdsr，式中功W本來是外荷載的功，但把彈性結構離散成n個自由度系統后，可用與外荷載等效的內力的功來代替。

1 微元體的微分剛度

微元體單元如圖2所示，它有三個移動自由度和三個轉動自由度（ωx、ωy、ωz），在小應變條件下，應力在移動位移中所做功之和與轉動位移相比可以略去不計。

圖2

幾何非線性微分剛度的有限元分析——

胡 可



第1期船 海 工 程第38卷

當微元體環繞z軸轉動ωz角度時（設為小轉動），①點向左位移1/2（1-cosωz），向上位移1/2sinωz，作用在點①上各力的功為

Wzl=-12（1-cosωz）σx+12sinωzτxy

總和①②③④各點力的功，由于剪應力的功互相抵消

Wz=-（1-cosωz）（σx+σy）≈-12ω2z（σx+σy）

根據(3)式得繞z軸轉動的微分剛度為Kdzz=σx+σy，考慮同時存在ωx、ωy、ωz三個轉動時，三個正應力做功為

Wσ=-12［（σx+σy）ω2z+（σy+σz）ω2x+（σz+σx）ω2y］

還有剪應力τxy對交叉乘機ωxωy也做功，因為轉動ωx時，引起②點向外和④點向里的位移均為1/2ωx，接著再轉動ωy，使②、④點均沿著剪應力方向位移1/2ωxωy，兩者做功之和為τxyωxωy，同理，加了另外兩對剪應力，共做功

Wτ=τxyωxωy+τyzωyωz+τzxωzωx

現在取無限小的微元體積ΔV，所有位移自由度都存在時，全部應力所做的總功為

ΔW=-ΔV2［（σy+σz）ω2x+（σz+σx）ω2y+（σx+σy）ω2z

-2τxyωxωy-2τyzωyωz-2τzxωzωx］

根據（3）式，微元體ΔV的微分剛度系數矩陣為

ΔKdωω=ΔV σy+σz-τxyσz+σx-τzx-τyzσx+σy（4）

有了微元體微分剛度陣，就可把每個單元再分成許多的微元體，按（4）式經積分后，便得各單元附加的微分剛度陣，對于小應變幾何非線性問題，先按線性剛度陣算出各單元的內力，接著便能算出各單元的微分剛度陣，把它與線性剛度陣加在一起，既得考慮幾何非線性問題的單元剛度陣，然后按通常的辦法，將單元剛度陣作坐標變換并組裝成總剛度陣，即可求解幾何非線性問題，為了提高計算精度還可以進行迭代運算。

2 桿單元的微分剛度

桿單元如圖3所示，桿ab長為l，截面積為A，指定x軸沿ab方向，則只有σx是非零應力，軸力Fx=σxA，能產生微分剛度的轉動角有

ωy=（uza-uzb）/l

ωz=（uyb-uya）/l

圖3

整個桿中，內力Fx做功為

W=∫ΔW=∫V-ΔV2σx（ω2y+ω2z）=-Fx2l［（uza-uzb）2+（uyb-uya）2］

根據（3）式，相應q=［uyauzauybuzb］T四個自由度的微分剛度是下式q前的系數陣

PyaPzaPybPzb=Fxl101-1 0 10 -1 0 1 uyauzauybuzb

3 等截面梁單元的微分剛度

等截面梁單元如圖4所示，在等截面梁ab中，除軸力Fx外，還有彎矩My、Mz及剪力Vy、Vz也產生微分剛度。

圖4

在梁單元x、y、z處取微元體ΔV，其上有應力σx、σxy、σzx。微元體的元功為

ΔW=-ΔV2［σx（ω2y+ω2z）-2τxyωxωy-2τzxωzωx］（5）

將微元體轉角用梁中性軸參數表示，即ωx=θx、ωy=-u′z-yθ′x、ωz=-u′y-zθ′x，代入（5）式，并積分得

W=-12∫V｛σx［（-u′z-yθ′x）2+（u′y-zθ′x）2］+2θx［τxy（u′z+yθ′x）-τxz（u′y-zθ′x）］｝dAdx

首先設q=［u′y u′z θ′x θx］T為廣義坐標，由于廣義坐標為梁中性軸的參數，上式積分應先在橫截面A上進行（即假設dx=1），于是得到單位長度梁的微分剛度為

［Kdij］=Fx0Fx-My -Mz IFxA-VzVy00

式中I為截面對形心的極慣性矩，再作梁的軸向積分，得梁單元的總功為

W=12∫baqTi［Kdij］qjdx

由(3)式得相對于梁端點位移q［uyauzaθxaθyaθzauybuzbθxbθybθzb］T為廣義坐標的微分剛度為

3 結束語

微分剛度計算對于大位移幾何非線性問題的有限元分析是一個比較適用的方法，對于這類問題，首先按上述方法得到微分剛度，然后用微分剛度修正有限元方法的求解方程，于是大位移幾何非線性問題的有限元分析得到簡化。

參考文獻

［1］ 胡 可.連續梁結構電算方法研究［J］.科技咨詢導報，2007(16).

［2］ 胡 可.結構電算方法中非節點荷載處理方法的研究［J］.科技與企業，2006(8).

［3］ 胡 可.平面組合結構電算方法研究［J］.科技咨詢，2006(25).

［4］ 龍馭球.有限元法概論［M］.人民教育出版社，1981.