緒論

用有限元素法分析平面應力或平面應變的問題時，最常使用的元素為平 面三角形元素或平面四邊形元素，對工程結構進行分割。對任何不規則形 狀的平面結構，我們都可以輕易地將其切割成有限的三角形組合，但不一 定適合將其切割成四邊形的組合，故三角形元素在文獻上被廣泛的探討及 使用。

最早出現的三角形元素為 CST 元素(constant strain triangle element)[1]，

其節點數目與節點參數如圖 1.1(a)所示。CST 元素是利用拉格朗日插值法 (Lagrange interpolation)將位移場假設為x與y的一次完整多項式，故在 CST 元素內部的應變及應力都是常數，這通常與實際情況不符，因此使用 CST 元素時，如需較精確的結果通常需要很密的網格。同樣地利用拉格朗日插 值法將位移場假設為x與y的二次完整多項式，可以得出 LST 元素(linear strain triangle element)[2]，其節點數目與節點參數如圖 1.1(b)所示。LST 元 素內部的應變與應力都是線性變化，所以精度較 CST 元素佳。同樣地利用 拉格朗日插值法將位移場假設為x與y的三次完整多項式，可以得到 QST 元素(quadratic strain triangle element)[3]，其節點數目與節點參數如圖 1.1(c) 所示。根據其節點數目與節點參數，文獻[4]中將其命名為 QST-10/20C 元 素。在拉格朗日元素(Lagrange elements)中，隨著位移場多項式的次數增高 可以提升精度，但是其元素節點不會單純的位於三角形的頂點上，且節點 參數缺少物理上對應的旋轉自由度(drilling degree of freedom)。

節點參數為何需要旋轉自由度？文獻[4]對於加入旋轉自由度的原因提 到：改善元素的性能並避免使用到三角形邊上的節點，因為邊上的節點會 影響到網格生成，而且在模擬非線性分析與動態分析時較為困難；當三角 平面元素與三角板元素疊加時，能滿足物理上一個節點有 3 個旋轉自由度

2

的要求；三角形元素與殼元素、板元素或是梁元素同時使用時，能使接合 簡單化。在有限元素法發展的過程中，有一個積極的目標是如何去創造高 性能(high-performance)的元素。文獻[5]中對高性能元素的定義為『簡單的 元素在粗糙的網格下卻可以擁有工程上的精度』。至於『簡單的元素』，文 獻[5]中的定義為『元素只有位於頂點的節點以及物理上定義的自由度』，一 般物理上定義的自由度是指位移和旋轉。

基於上述理由，在 1964 至 1983 年期間，許多人在研究如何在 3 節點 的平面三角形元素上加入節點旋轉自由度，希望能得到一個 3 節點 9 個自 由度且具節點旋轉自由度的平面三角平面元素，但是都沒有得出可用的元 素。一直到 1984 年文獻[6]與 1985 年文獻[7]才成功的得出帶有旋轉自由度 3 個節點 9 個自由度的三角形元素。

從 1984 年至今，有許多具旋轉自由度三角形元素被提出[4,6-20]。其中 除了文獻[19]的 DLST 元素是 12 個自由度外，其餘的元素都是 3 個節點 9 個自由度的元素。文獻上用來測試具旋轉自由度之三角形元素性能的例題 通常有：懸臂梁受彎矩[7]、懸臂梁受端點剪力[3]及 Cook 例題[21](一懸臂 梯形平板受端點剪力)。從文獻上的結果可以發現當元素網格很密時，所有 元素的結果都差不多，但使用 2×2 的網格分析懸臂梁受彎矩的例題時，僅 有文獻[4]中之 OPT 元素有很好的結果。使用 2×2 的網格分析 Cook 例題時，

除了文獻[19]的 DLST 元素外，其他元素的結果都不太好。文獻[19]稱 DLST 元素的性能和 LST 元素的性能差不多。由文獻的結果可發現在使用同樣數 目的元素網格時，LST 元素的結果比具旋轉自由度之 3 節點 9 個自由度之 元素的結果好，這應是合理的，因 LST 元素有 12 個自由度。

QST 元素在很多有限元素法的書上[22-25]都有提到，文獻[4]中將 QST 元素分成 20 個自由度和 18 個自由度兩大類。20 個自由度的 QST 元素之位 移場是完整的三次多項式[22]，18 個自由度的 QST 元素之位移場是不完整

的三次多項式[22,26]，QST 元素可以使用不同數目的節點及不同的節點參 數[4,23]。文獻上使用 QST 元素分析例題的結果甚為少見，這可能是因為 QST 元素不滿足簡單元素的要求，且現在的趨勢是使用低階元素。但由文 獻[26]的例題，可以發現 QST 元素有很好的性能。3 節點 18 個自由度的 QST 元素中有一種的節點自由度為 2 個位移、1 個旋轉及 3 個應變，在文獻[4]

中稱其為 QST-3/18RS 元素。廣義的來說，應變也是物理上的自由度，所以 該 QST 元素應可稱為高階的簡單元素。因文獻上缺乏用 QST 元素分析例題 的結果，故本研究的目的之一為用 QST-3/18RS 元素分析文獻上的例題，以 探討其性能。

造成幾何非線性的原因主要為剛體的大位移及旋轉，若用共旋轉法[27]

將剛體位移及旋轉從總位移及旋轉中扣除，則剩下的變形位移及旋轉仍為 小位移及旋轉。所以若使用共旋轉法，則在線性分析使用的元素可以用在 大位移小應變問題之幾何非線性分析。共旋轉法在梁與殼結構的幾何非線 性分析已經被廣泛的使用[27-34]，但在平面應變的問題，則很少見到被使 用。據本人所知，僅有文獻[35]提出一個 4 節點平面元素的共旋轉推導法，

並成功地使用在平面問題的幾何非線性分析。文獻[35]使用的 4 節點元素並 不具旋轉自由度。據本人所知，文獻上並沒有具旋轉自由度之三角平面元 素的共旋轉法被提出，所以本研究的主要目的是提出一具旋轉自由度之三 角平面元素的共旋轉推導法，將 QST-3/18RS 元素應用於平面應變問題的幾 何 非 線 性 分 析 ， 本 研 究 擬 採 用 類 似 文 獻 [27] 中 的 共 旋 轉 法 來 描 述 QST-3/18RS 元素的變形，並提出一運動過程來決定當前元素座標的元素節 點位移、節點變形角及節點應變。

因文獻上並無現成 QST-3/18RS 元素的形狀函數可以用，本文將在第二 章利用文獻[22]中以節點位移及位移的一次微分為節點參數之 QST 元素(文 獻[4]中稱其為 QST-3/18G 元素)的形狀函數推導 QST-3/18G 元素的剛度矩

4

陣，然後推導不同元素節點參數間的轉換關係，再利用不同等效節點參數 的轉換關係求出 QST-3/18RS 元素的剛度矩陣及元素的節點內力。本文並將 在第二章說明如何使用共旋轉法來決定元素變形參數。在第三章中將介紹 本文的數值計算方法及程序。在第四章中將先以線性例題探討 QST-3/18RS 元素的性能，再以非線性例題說明本文提出的共旋轉法及決定元素節點變 形參數的方法是可行的。本文也將測試不同的邊界條件、元素網格及元素 形狀對於數值結果的影響。

二、理論推導

本章中將提出一個具旋轉自由度(drilling degree of freedom)之三角平面 元素(triangular plane element)的共旋轉推導法(co-rotational formulation)

，本章中採用文獻[22]中 3 節點 18 個自由度埃爾米特元素(Hermite elements) 的形狀函數(shape functions)，但轉換其節點參數，使其具旋轉自由度。在 本章將描述元素變形的位移場與節點參數，推導元素節點內力與元素剛度 矩陣，亦將提出一個決定元素節點變形參數的方法。

2.1 基本假設

本文中對三角平面元素的變形做如下的假設：

元素的變形位移及旋轉為小位移及小旋轉。

2.2 座標系統

為了描述系統的運動及元素的變形，本文定義了三組座標系統：

(a) 固定總體座標系統(global coordinate system)：X_i^G

(

)

如圖2.1 所示，結構體所有節點的座標、系統的邊界條件與其他座標系統的 基底，以及結構的平衡方程式，均在此座標系統中定義。在X_i^G座標系統內 之座標值以

(

)

(b) 元素座標系統(element coordinate system)：x_i^E

(

)

如圖2.1 所示，此座標系統是建立在每一個元素變形後的最新位置上，其座 標原點為元素節點 1，x₁^E軸為元素節點 1 與元素節點 2 在元素平面上的連 線，x₂^E軸是在元素平面上垂直於x₁^E軸，且朝著元素節點3 的方向。元素的 位移、元素變形、元素節點內力與元素剛度矩陣是在此座標系統中定義，

然後經由座標轉換，將其轉換至總體座標系統及基礎座標系統。本文中⁰x_i^E

6

表示初始未變形時的元素座標，^Ix_i^E表示第I 個增量迭代收斂後的元素座 標，x_i^E表示當前變形位置的元素座標。在x_i^E座標系統內之座標值以

( )

(c) 節點基礎座標系統(base coordinate system)：x_ij^B

(

)(

)

如圖2.1 所示，此座標系統的原點是剛接在結構離散後的每一個節點，並與 對應的節點一起移動及旋轉。本文中其初始位置與固定總體座標系統一 致，節點的應變是在此座標系統中定義。本文中⁰x_ij^B表示元素節點 j在初始 未變形時的節點基礎座標，^Ix_ij^B表示元素節點 j在第I個增量迭代收斂後的 節點基礎座標，x_ij^B表示元素節點 j在當前變形位置的節點基礎座標。

2.3 三角平面元素的變形描述

本 文 採 用 的 三 角 平 面 元 素 為 文 獻[22] 中 不 完 整 三 階 埃 爾 米 特 元 素 (incomplete cubic Hermite element)，其位移場為三次變化，應變場則為二次 變化，因此又稱為QST 元素(quadratic strain triangle element)。此元素有 3 個節點，每個節點有 6 個自由度，此元素可使用彼此間能互相轉換的不同 節點參數。本文中採用的節點參數為節點 j在元素座標x₁^E、x₂^E軸的位移分 量u 、_j v ，應變分量_j ε_xj、ε_yj、γ_xyj及逆時鐘方向的旋轉θ_j。但為了方便推 導，本文在元素推導時使用的節點參數是節點 j在x₁^E、x₂^E軸的位移分量 u 、j v 以及_j x₁^E、x₂^E軸的位移分量分別對面積座標(area coordinates)ξ、η的 微分u,_ξj、u,_ηj、v,_ξj、v,_ηj，面積座標的介紹詳見附錄 A。本章中將推導不 同節點參數間的轉換關係，再推導元素節點內力及元素剛度矩陣。此元素 的位移場可表示為[22]：

qξ

N^t_u

u = (2.1)

qξ

8

θ

將(2.20)式與(2.23)式代入(2.1)式與(2.2)式可得：

θ

將(2.26)式、(2.27)式代入(2.7)至(2.9)式可得：

θ

10

其中E是楊氏模數(Young’s modulus)，ν 是蒲松比(Poisson’s ratio)。將(2.30) 式代入(2.33)式可得： 廣義節點力矩，而非應力。將(2.30)式、(2.36)式代入(2.37)式可得

θ

θ

12 轉。將(2.44)式、(2.46)式、(2.51)式代入(2.17)式可得：

B

由反梯度法則(contragradient law)[23]及(2.52)式可得：

fθ

應於θ^B_j 的節點傳統力矩。將(2.41)式、(2.52)式代入(2.57)式可得：

B B

B k q

f = (2.60)

EB t

EB

B T k T

k = _θ (2.61)

其中k_B是對應於q_B的元素剛度矩陣。

2.6 元素節點變形參數的決定方法

本文中採用增量迭代法解非線性平衡方程式，假設第I 個位置已知，此 處的第I 個位置，是指第I 個增量迭代收斂後的平衡位置。^Ix_j為元素在第I 個位置的元素座標，^I X_j、^Ix_ij^B、ΔU_j、Δε^B_j 以及Δθ_j^B分別為元素節點 j在 固定總體座標中第I 個位置的位置向量、節點基礎座標軸、增量位移向量、

增量應變向量以及增量旋轉。本文中假設元素節點受ΔU_j、Δε^B_j 以及Δθ^B_j 作

增量應變向量以及增量旋轉。本文中假設元素節點受ΔU_j、Δε^B_j 以及Δθ^B_j 作