本文為了模擬 2.1 節所敘述之問題，將採兩階段的位移負荷來分 析。在第一階段的模擬中，同時施加三組位移負荷，這三組位移負荷 分別為施加 向的位移X λd於X =0的左夾持端，和施加 向的轉角Y λθ 於X=W的右夾持端，以及施加 向的轉角-Y λθ於X=0的左夾持端，
如圖 3.1(a)所示。為了使施加的位移與轉角這三組位移負荷同時到達 期望值，位移與轉角三組負荷參數間必須維持一固定比率。待位移與 轉角達到給定期望值，則形成圓柱薄板受側向集中負荷之前的初始狀 態。
第一階段進行平衡迭代時的數值方法大部分採用文獻[42]中所提 出基於牛頓-拉福森(Newton-Raphson)法和弧長控制(arc length control) 法的增量迭代法，而為了使第一階段的邊界位移與轉角精確地達到期 望值，將於第一階段的最後一步增量使用標準牛頓法來達成。第二階
段一開始先重新設定兩夾持端的邊界條件，再施加一個 向的集中
位移負荷
−Z
λE於圓柱板的中心(E 點)，如圖 3.1(b)所示。第二階段則全 程採用文獻[42]中所提的增量迭代法。
本文中為了求得分歧點，以系統切線剛度矩陣之行列式值為零來 判斷。本文採用二分法決定增量位移向量的長度，以求得系統切線剛 度矩陣之行列式值為零的平衡位置。為了求得次要平衡路徑，本文在 平衡路徑的第一個挫曲負荷分歧點加入一個與第一挫曲模態向量方 向一致的擾動位移。
3.1 增量迭代法
若第 I 個增量的平衡位置為已知，令其位移向量為QI、負荷參
數為λI，則在此位置的系統切線剛度矩陣 可以求得，且第 I+1 個 增量的初始增量位移向量
KT
∆ ，可以利用尤拉預測值(Euler predictor)Q 求得[40]
2.5 節的方法，則可求得系統中各殼元素在當前的元素座標、節點位 移、節點變形角。再利用 2.4 節的方法求得當前元素座標上的節點內 力及剛度矩陣，再由座標轉換轉到總體座標的節點內力及剛度矩陣。
當系統內力及外力求得後，不平衡力向量 可由(2.54)求得。若(2.55) 的收斂準則不能滿足，則利用定弧長控制法[41]求得一位移修正量
Ψ
δQ與負荷參數修正量δλ，並加入前一次迭代的∆ 與Q ∆λ中，而得到 一新的增量位移向量與增量負荷參數，再進行下一次迭代。δQ與δλ 可由下列兩式決定
δQ=−KT−1(ΨI+1+δλFP) (3.5) δl2 =(∆Q+δQ)t(∆Q+δQ) (3.6) 以上迭代過程一直重複至不平衡力滿足(2.55)式的收斂準則。
3.2二分法
利用 3.1 節的增量迭代法可以求得結構之主要平衡路徑。在每個
增量的迭代收歛時，可以得到該增量在其平衡位置的負荷參數λ及結
構剛度矩陣的行列式值D(λ)。令λI 及D(λI)分別表示第 I 個增量在其 平衡位置的λ及D(λ)值。λI+1及D(λI+1)分別表示第 I+1 個增量在其 平衡位置的λ及D(λ)值。∆lI+1表示第 I+1 個增量的增量位移向量之 弧長。若D(λI)大於零且D(λI+1)小於零則可利用以下二分法求得挫曲 負荷參數λcr：
(1) 令∆lL =0，∆lR =∆lI+1，λL = ，λI λR =λI+1。其中下標 L 及 R 表示 左界及右界。
(2) 取
1 2
R L I
l l = ∆l +∆
∆ + ，重做第 I+1 個增量迭代，並求得新的λI+1及 )
( I+1
D λ 。
(3) 若D(λI+1)大於零，則令λL =λI+1，∆lL =∆lI+1。若D(λI+1)小於零，
則令λR =λI+1，∆lR =∆lI+1。
(4) 若下列二式挫曲誤差準則同時滿足，則λI+1即為系統挫曲負荷，否 則回到步驟(2)，重新展開下一次二分法迭代。
I eD D
D + <
) 0 (
) (λ 1
λ
λ λ
λ e
I L R − <
+1
其中 及 為給定的容許誤差值，本文例題之計算給定
及 。
eD eλ eD =10−8
10−5 λ = e
3.3數值程序
本文針對此問題之模擬與分析，採用的數值程序說明如下：
1. 輸入初始的結構資料、材料常數，與分析所需的階段數目、每個
階段的位移負荷數，以及轉換至下一階段所需變換的自由度數目。
2. 設定該階段結構的邊界狀態，選擇該階段的參考自由度，並給定
期望此自由度於該階段應達到的位移，以及設定該階段收斂時的 容許誤差。
3. 設定該階段的位移負荷，若該階段有兩個以上之位移負荷同時作
用，則輸入位移負荷間的比例大小。
4. 使用迭代法求增量的收歛解
(a)利用已知的增量位移求當前元素變形並計算元素節點內力。
(b)計算不平衡力Ψ，檢查是否滿足收斂準則，若滿足，則進行 5.；
否則減少增量弧長並計算新的增量位移與增量負荷，回到(a)重 新此增量。
5. 計算下一次增量所需的資料
(a)檢查該階段之參考自由度的位移是否已達的給定值。若已達到
給定值則進行6.；否則進行下一步工作。
(b)計算切線剛度矩陣 與參考負荷向量 ，並計算下一個增量
的增量位移、增量弧長與增量負荷參數。
KT FP
(c)回到2.執行迭代工作。
6. 更換邊界狀態，並把邊界條件更動的節點位移儲存，供下一階段
的分析工作所使用。