1、1.1.彈性動力學(xué)有限元基本解法結(jié)構(gòu)系統(tǒng)的通用運動學(xué)方程為：MU CU KU 二 Rt（ 1）求解該動力學(xué)振動響應(yīng)主要有三類方法：（1）時域法（2）頻域法（3）響應(yīng)譜法 時域法又可分為：（1）直接積分法，（2）模態(tài)疊加法。直接積分法又可分為中 心差分法（顯式）， Wils on '（隱式）法以及Newmark （隱式）法等。本文介紹中心差分法（顯式）與 Newmark （隱式）法。1中心差分法（顯式）假定0, t1, t2,，tn時刻的節(jié)點位移，速度與加速度均為已知，現(xiàn)求解 tn（t *：t）時刻的結(jié)構(gòu)響應(yīng)。中心差分法對加速度，速度的導(dǎo)數(shù)采用中心差分代替，即 為：UtUtRt1 2（U

2、t_t -2Ut Ut T Lt1*（UtUt J將M?U式M?（2）式代入（1）式后整理得到t R（3）（3）中112 M C ：t22：t2 1 1二 Rt _（k - M）Ut -（M -石 C）U 分別稱為有效質(zhì)量矩陣，有效載荷矢量。R，M，C，K為結(jié)構(gòu)載荷，質(zhì)量，阻尼，剛度矩陣。求解線性方程組（3），即可獲得r . ：t時刻的節(jié)點位移向量Ut.、.t，將U,*代回幾 何方程與物理方程，可得t時刻的單元應(yīng)力和應(yīng)變。中心差分法在求解t .迸瞬時的位移Ut t時，只需t 氏時刻以前的狀態(tài)變量Ut和 Utj：，然后計算出有效質(zhì)量矩陣M?，有效載荷矢量R，即可求出U,t，故稱此解法 為顯式算法

3、。中心差分法，在開始計算時，需要仔細(xì)處理。 t=0時，要計算U t，需要知道u J 的值。因此應(yīng)該有一個起始技術(shù)，因而該算法不是自動起步的。由于 Uo， Uo， Uo是 已知的，由t=0時的（2）式可知： : t2 "Uj 二 U。- ：tU° U。2中心差分法中時間步長：t的選擇涉及兩個方面的約束：數(shù)值算法的穩(wěn)定性和計算 時間。中心差分法的實質(zhì)是用差分代替微分，并且對位移和加速度的導(dǎo)數(shù)采用線性外 插，這限制了氏的取值不可過大，否則結(jié)果可能失真過大。可以證明：中心差分法是條件穩(wěn)定的。即當(dāng)時間步長-：t必須小于由該問題求解方程性質(zhì)所決定的一個時間步長的臨界值。LS-DYNA中

4、，采用“變時間步長法”，即每一時刻的步長 氏由當(dāng)前結(jié)構(gòu)的穩(wěn)定性條件來控制。具體算法為：計算每一個單元的極 限時間步長"：tei，i=1,2，取t = mi n（Utei）為下一個時刻的時間步長。各種單元的 八tei計算方法如下。(1) 1D桿，梁單元e C其中為時間步長因子，系統(tǒng)默認(rèn)為0.9。L為桿，梁單元的長度。E材料聲速。(2)2D板，殼單元觸LminC其中為時間步長因子，系統(tǒng)默認(rèn)為0.9。Lmin為殼單元的最小單元邊長度E2 為材料的聲速。(1()(3) 3D單元Let e22 1/2Q (Q2 C2)其中GC +C°Le ®kk ( % £ 0)

5、0( ；kk - 0)Co和Ci為無量綱常數(shù)，默認(rèn)Co =1.5,G =0.06.V (對8節(jié)點體單兀)、Le =二AemaxLe為單元等效長度，Ve為單元體積，Aemax為單元Lmin(對4節(jié)點體單元)最大側(cè)面積。C=J E (1-v)為材料聲速。V(1+v)(1-2v) P時間步因子可由用戶設(shè)置，減小:相當(dāng)于減少時間步長。設(shè)置時間步長因子: 的關(guān)鍵字為*CONTROL_TIMESTEP，控制參數(shù)為TSSFAC。另外，質(zhì)量縮放可以人為控制時間步長。即調(diào)整單元密度，來改變時間步長。以殼單元為例說明質(zhì)量縮放改變時間步長。氏廠血(假如=1)C4 specified、2 由(飛) min(1 -)

6、： i 得至卩-("tspecified ) E得"Lu2)LS-DYNA中有2種質(zhì)量縮放方案，修改*CONTROL_TIMESTEP中的參數(shù)。方案1: DT2MS為正的時間步通過調(diào)整單元密度，使所有單元時間步相同，只用于慣性效應(yīng)不重要的情況。方案2: DT2MS為負(fù)的時間步質(zhì)量縮放只用于小于指定時間步長 DT的單元。慣性效應(yīng)應(yīng)該通過變形體的動能與內(nèi)能的比例進(jìn)行衡量(一般應(yīng)該小于 10%等)2 Newmark 法(隱式)Newmark假定在時間間隔t,t .4內(nèi)，加速度線性變化，即采用如下的加速度，速 度公式：Ut t 二Ut 心一、)5：t 1 " “ 2Ut&

7、#174; 小()5 5小 t式中:為按積分的精度和穩(wěn)定性要求可以調(diào)整的參數(shù)。根據(jù)(4)式可給出Ut . t和Ut. t用Ut. t，Ut，Ut表示的表達(dá)式，代入(1)式中整 理得到KUt t 朋 t (5)其中K 2 M C Kr, =t2，J1 1 1Rt t =Rt t M 2Ut Ut (1)Ut CUt (1)Ut (1).：tUt'':八t:汀t2t2:稱之為有效剛度矩陣和有效載荷矢量。由上式可以看出求解當(dāng)前Ut.t，需要用到當(dāng)前時刻的R t，因此該算法為隱式算法。當(dāng)載荷歷史全部已知時，F(xiàn)t . t為已知量，求解需要迭代實現(xiàn)。可以證明，當(dāng)參數(shù)-0.5，：_0.25(

8、0.5，)2時，Newmark法是無條件穩(wěn)定的，即 t的大小不影響數(shù)值穩(wěn)定性。此時時間步長:t的選擇主要根據(jù)解得精度確定。一般，Newmark法可以比中心差分法的時間步長大得多。3.結(jié)論比較兩種算法，顯式中心差分法非常適合研究波的傳播問題，如碰撞、高速沖 擊、爆炸等。分析式(3)發(fā)現(xiàn)，顯式中心差分法的 M與C矩陣是對角陣，如給定某 些有限元節(jié)點以初始擾動，在經(jīng)過一個時間步長后，和它相關(guān)的節(jié)點進(jìn)入運動，即U中這些節(jié)點對應(yīng)的分量成為非零量，此特點正好和波的傳播特點相一致。另一方面， 研究波傳播的過程需要微小的時間步長，這也正是中心差分法的特點。而Newmark法更加適合于計算低頻占主導(dǎo)的動力問題，

9、從計算精度考慮，允許采 用較大的時間步長以節(jié)省計算時間，同時較大的時間步長還可以過濾掉高階不精確特 征值對系統(tǒng)響應(yīng)的影響。隱式方法要轉(zhuǎn)置剛度矩陣，增量迭代，通過一系列線性逼近(Newton-Raphsor)來求解。正因為隱式算法要對剛度矩陣求逆，所以計算時要求整 體剛度矩陣不能奇異，對于一些接觸高度非線性問題，有時無法保證收斂。1.2.1顯式算法顯式算法基本假定為：在一微小時間段內(nèi)，模型任意點速度、加速度為常 數(shù)°AEAQUS軟件Explicit模塊應(yīng)用中心差分法對運動方程進(jìn) 行顯式時間積分， 運動方程的解為'u(i)=M-1-(F(i)-I(i)(1)式中：M為集中質(zhì)量矩陣

10、；F 為外荷載向量；I為單元內(nèi)力向量。 由于顯式算法中不需要對剛度矩陣求逆，集中質(zhì)量矩陣為對角矩陣，求逆簡 便，使顯式算法并行計算數(shù)據(jù)傳輸量較小；且顯式算法剛度矩陣大小與自由 度數(shù)成線性關(guān)系，因此顯式算法用于自由度數(shù)龐大的數(shù)值計算時具有很大優(yōu) 勢。1.3.2隱式算法隱式算法含義為：t+ t時刻狀態(tài)不僅與t時刻狀態(tài)有關(guān)，且與t+At時刻某些量有關(guān)°AEAQUS軟件Standard應(yīng)用H ilber Hughes-Taylor 隱式算法、NcwtonRaphson迭 代法進(jìn)行動力方程求解10。運動方程解為 u(i + l) = Au(i)+K t1(F(i)-I(i)(2)式中：Kt為當(dāng)

11、前切線剛度矩陣：Au為位移增量。求解方程位移增量Au(i + 1)時，必須對剛度矩陣K求逆。當(dāng)自由度數(shù) 非常龐大時，這項計算消耗資源。對K矩陣的求逆計算，計算機之間數(shù)據(jù)傳輸量非常大，隨著自由度數(shù)增加，剛度矩陣K大小成指數(shù)增長。因此，隱式算法用于自由度數(shù)龐大的數(shù)值計算時，優(yōu)勢不明顯，甚至?xí)档陀嬎阈?率。根據(jù)地鐵地下結(jié)構(gòu)抗震研究需要，對有限元并行計算顯式算法和隱式算法計 算精度和效率進(jìn)行比較，主要結(jié)論如下：(1) 有限元并行計算中心差分顯式算法與H ilberHughes Taylor隱式算法計算精度相當(dāng)，顯式算法計算效率遠(yuǎn)高于隱式算法。(2) 黏彈性人工邊界在顯式算法和隱式算法中，都能起到很好的模擬效 果。(3) 顯式算法和隱式算法計算相對位移、相對速度時程基本一