本文解(2.69)式的非線性平衡方程式所使用平衡迭代的數值計算方 法是採用文獻[38]中所提出基於牛頓-拉福森(Newton-Raphson)法和弧 長控制(arc length control)法的增量迭代法。為了本文的完整性,以下將 簡單介紹文獻[38]中提出的增量迭代數值計算方法與程序。
3.1 增量迭代法
若第I 個增量的平衡位置為已知,令其位移向量為QI、負荷參數 為λI,則第I +1個增量的初始增量位移向量Δ ,可以利用尤拉預測Q 值(Euler predictor)求得[39]
rT
Q=Δλ
Δ (3.1)
( )
K PrT = TI −1 (3.2)
其中Δ 為初始增量負荷參數,λ KTI 為第I 個平衡位置的系統切線剛度矩 陣。(3.1)式中的Δ 可利用下式求出[38] λ
( )
T 1/2 tl rTr Δ
±
=
Δλ (3.3)
其中正負符號決定方法為,當第I −1與I 個增量收歛時,其系統切線 剛度矩陣之行列式值同號,則Δ 的正負符號和第λ I 個增量時相同;若 異號則符號相反。 lΔ 表示第I +1個增量的增量弧長,其值可以如下決 定[38]
(
D I)
12I J J
l l =Δ
Δ (3.4)
其中JD為給定的期望迭代次數,JI為第I 個增量迭代至平衡所使用的 迭代次數,ΔlI為第I 個增量的增量弧長。本文第一個增量的增量弧長
l1
Δ 是由下式決定
17
18
其中r =δ Q為增量位移修正量,整理(3.10)式,可得到位移修正量
(
I)
TT Ψ P r r
K
r= − −1 +1 −δλ = 0 +δλ (3.11) 其中r0 =−KT−1ΨI+1為標準牛頓法的位移修正量,rT在(3.2)式中已定義。
3.2 弧長控制法
(3.11)式中的負荷參數修正量δλ可利用文獻[38]中所提出的定弧長 控制法決定,其方法在每一個增量中固定其增量位移向量的 Euclidean norm 為一定值Δl,由新的增量位移向量
(
ΔQ+r)
可以得到(
Q+r) (
Q+r)
= Δ Δ
Δl2 t (3.12)
將(3.11)式代入(3.12)式可以得到
(
T) (
t T)
l2 = ΔQ +r0 +δλr ΔQ+r0 +δλr
Δ (3.13)
(3.13)式經過整理後可以得到δλ的二次方程式
3 0
2 2
1 + a +a =
aδλ δλ (3.14)
其中
T t
a1 =rTr (3.15)
( )
t Ta2 = 2 ΔQ+r0 r (3.16)
(
0) (
0)
23 Δ Δ Δl
a = Q+r t Q+r − (3.17) 將(3.11)式求得之增量位移修正量r以及(3.14)式求得之負荷參數修正 量δλ加入上次迭代之ΔQ與Δ 中,可以得到新的增量位移向量與增量λ 負荷參數,再進行下一次迭代,迭代過程將一直重複至(3.6)式中的系 統不平衡力向量滿足(2.70)式的收斂準則為止。
3.3 數值程序
本文使用的增量迭代法之數值程序可以分為三個主要部分:
1. 輸入並計算開始分析所需的資料
19
(a)輸入結構資料、邊界條件與負荷參數。
(b)選擇一個參考自由度,並給定期望此自由度應達到的位移。
(c)給定最大增量次數、最大迭代次數、期望迭代次數與容許誤差值。
(d)形成系統剛度矩陣並求得(3.5)式中的R0。
(e)利用(3.1)式、(3.3)式與(3.5)式計算初始增量位移向量、初始增量 負荷參數與第一次增量弧長。
2. 使用迭代法求增量的收歛解
(a)由 2.6 節的方法,利用已知的增量位移求得當前元素的變形向量 qθ,並計算(2.41)式中元素節點內力fθ,接著由2.5 節座標轉換將fθ 利用(2.57)式轉換至fB,然後將元素節點內力fB組合成結構系統節點 內力F 。
(b)計算(2.69)式的不平衡力Ψ 。
(c)檢查(2.70)式的收斂準則,若滿足則進入第 3 部分;否則檢查迭代 數,如果小於給定之最大迭代次數,則進行步驟(d);否則減少增量 弧長並以(3.1)式與(3.3)式計算新的增量位移向量與增量負荷參數,
回到步驟(a)重新計算。
(d)利用(3.11)式與(3.14)式計算增量位移修正量與增量負荷參數修正 量,回到步驟(a)重新計算。
3. 計算下一次增量所需的資料
(a)檢查參考自由度的位移及進行的增量次數是否已達給定值,若已 達到給定值則停止分析工作;否則進行步驟(b)。
(b)計算(3.11)式中的切線剛度矩陣KT。本文中KT的計算方法是利 用當前變形位置的元素座標重新計算(2.42)式中的元素剛度矩陣 kθ,本文中忽略元素其他的幾何變形與內部應力對元素剛度矩陣造 成的影響,再利用 2.5 節座標轉換將kθ利用(2.61)式轉換至kB,然
20
後將元素剛度矩陣kB疊加成結構系統剛度矩陣KT。
(c)利用(3.1)式、(3.3)式與(3.4)式計算下一次增量的增量位移向量、
增量負荷參數與增量弧長。
(d)回到第 2 部分進行迭代工作。