雙對稱斷面開口薄壁梁之幾何非線性動態分析 的一致性共旋轉推導法

國

立

交

通

大

學

機械工程學系碩士班

碩 士 論 文

雙對稱斷面開口薄壁梁之幾何非線性動態分析

的一致性共旋轉推導法

A consistent co-rotational formulation for Geometric Nonlinear

Dynamic Analysis of Doubly Symmetric Thin-walled open-section

Beams

研 究 生：高嘉鴻

指導教授：蕭國模 博士

雙對稱斷面開口薄壁梁之幾何非線性動態分析

的一致性共旋轉推導法

A consistent co-rotational formulation for Geometric Nonlinear Dynamic Analysis of Doubly Symmetric Thin-walled

open-section Beams

研 究 生： 高嘉鴻 Student： Chia-Hung Kao

指導教授： 蕭國模 博士 Advisor： Dr. Kuo-Mo Hsiao 國 立 交 通 大 學

機械工程學系碩士班

碩 士 論 文 A Thesis

Submitted to Department of Mechanical Engineering College of Engineering

National Chiao Tung University in Partial Fulfillment of the Requirements

for the Degree of Master of Science

in

Mechanical Engineering August 2014

Hsinchu, Taiwan, Republic of China

雙對稱斷面開口薄壁梁之幾何非線性動態分析

的一致性共旋轉推導法

A consistent co-rotational formulation for Geometric Nonlinear

Dynamic Analysis of Doubly Symmetric Thin-walled

open-section Beams

研究生：高嘉鴻 指導教授：蕭國模博士

國立交通大學機械工程學系碩士班

摘要

本研究以一致性共旋轉法推導一個雙對稱斷面開口薄壁尤拉梁元素， 探討三維梁之非線性動態反應。 本文中推導的梁元素有兩個節點十四個自由度。本研究用元素節點對 固定座標的位移及旋轉向量更新其位置及其斷面的方位，並在梁元素節點 當前的位置及斷面方位建立一個固定元素座標及與其重合的移動元素座

標。本研究用三個旋轉參數來描述元素斷面的方位與元素座標的關係。本 研究利用在當前固定元素座標的元素節點位移、旋轉向量及其擾動量、速 度、加速度、角速度、角加速度，推導出元素節點在當前固定元素座標的 擾動位移和擾動旋轉向量與元素節點旋轉參數之擾動量的關係、擾動後之 移動元素座標與當前固定元素座標的關係、移動元素座標的角速度及角加 速度、元素節點的變形參數對時間的一次及二次微分。本研究利用虛功原 理和 D’Alembert 原理，以及完整的幾何非線性梁理論的一致性二階線性化 在當前的固定元素座標推導元素節點變形力及慣性力。本研究由元素節點 變形內力的改變與擾動位移的關係推導梁元素的切線剛度矩陣。本研究將 元素節點慣性力表示成元素節點之絕對速度、加速度、角速度、角加速度 的函數，所以元素的質量矩陣與陀螺效應之關係矩陣分別由元素的節點慣 性力對元素之節點加速度向量的微分及元素之節點速度向量的微分求得。 本研究發現以一致性共旋轉推導法與共旋轉全拉格朗日推導法得到的節點 慣性力有些差異，但這些差異會在元素增加時會趨近於零。 本研究以數值例題說明本研究之梁元素的準確性，並探討本文之一致 性共旋轉法與文獻上共旋轉全拉格朗日法在動態分析時之差異性。另外本 研究比較 Newmark 積分法與中央差分法在不同例題之效率及準確性。

A consistent co-rotational formulation for Geometric Nonlinear

Dynamic Analysis of Doubly Symmetric Thin-walled

open-section Beams

Student：Chia-Hung Kao Advisor：Dr. Kuo-Mo Hsiao

Department of Mechanical Engineering National Chiao Tung University

Abstract

A consistent co-rotational (CCR) finite element formulation for geometrically nonlinear dynamic analysis of doubly symmetric thin-walled beam with large rotations but small strain is presented. The element developed here has two nodes with seven degrees of freedom per node. The element nodes are chosen to be located at the centroids of the end cross sections of the beam element and the centroid axis is chosen to be the reference axis. A rotation vector is used to represent the finite rotation of coordinate systems rigidly tied to each node of the discretized structure. The incremental nodal displacement vectors and rotation vectors in global coordinates are used to update the node locations and orientation of the element. The deformations of the beam element are described in a current moving element coordinate system constructed at the current node locations and orientation of the beam element. Three rotation parameters are defined to describe the relative orientation

unwrapped cross section and the current element coordinate system. The element equations are derived in a fixed current element coordinates which are coincident with the current moving element coordinates. The perturbed displacements and spatial rotation, velocity and acceleration, angular velocity and angular acceleration of the current moving element coordinates, and the variation of the element nodal rotation parameters corresponding to the perturbation of element nodal displacement vectors and rotation vectors in the current fixed element coordinates are consistently determined and expressed in terms of the current element nodal displacements and rotation parameters, nodal velocities and accelerations, and nodal angular velocities and angular accelerations. The element deformation and inertia nodal forces are derived using the virtual work principle, the d’Alembert principle, and the consistent second order linearization of the fully geometrically nonlinear beam theory. In element deformation nodal forces, all coupling among bending, twisting, and stretching deformations of the beam element is considered. In the element inertia nodal forces, the terms up to the second order of time derivatives of deformation parameters are retained. However, the coupling between rotation parameters and their time derivatives are not considered in the element inertia nodal forces. In this study, the element inertia nodal forces are expressed in terms of element nodal velocities and accelerations, and nodal angular velocities and accelerations. Thus, the element inertia matrices may be obtained by differentiating the element inertia nodal forces with respect to the element nodal velocities and accelerations.

There is a slight difference between the present element inertia nodal forces and that derived using corotational total Lagrangian (CRTL) formulation. However, the effect of the difference may be negligible with the decrease of

element size.

An incremental-iterative method based on the Newmark direct integration method and the Newton-Raphson method is employed for the solution of nonlinear equations of motion. Numerical examples are presented to investigate the accuracy and efficiency of the proposed method. It is found that the difference between the dynamic responses obtained using CCR formulation and CRTL formulation is negligible for all examples studied.

誌謝

衷心感謝指導教授 蕭國模博士這兩年在學業上的指導以及生活上的 照顧，使得本文能夠順利完成，在此致上最高的敬意和謝意。同時也感謝 蔡佳霖老師、尹慶中老師擔任口試委員，對本文提出了寶貴的建議。非常 感謝黃楚璋、林琮琪、沈佳鴻、莊士緯學長對於學業的細心指導及日常生 活的照顧，同學林群禮、學弟金長虹、蔡耀賢的朝夕相處以及課業上的相 互討論，使我這兩年能順順利利的度過，在我的生活上亦幫助很多，十分 感謝你們。 最後要感謝我的家人、親戚及朋友的鼓勵與支持，使我能完成這段學 業，在此以此成果獻給所有關心我的人。

目錄

中文摘要 ...I 英文摘要 ... III 誌謝 ...VI 目錄 ... VII 表目錄 ... X 圖目錄 ...XI 第一章 緒論 ... 1 第二章 理論推導 ... 5 2.1 基本假設... 5 2.2 座標系統描述... 5 2.3 旋轉向量... 7 2.4 梁之剪心軸的位移及其斷面的旋轉... 9 2.5 節點參數與節點力... 13 2.6 元素斷面座標、移動元素座標、元素節點位移及旋轉 參數之決定... 16 2.7 移動元素座標與固定元素座標的關係 ... 17 2.7.1 擾動後之移動元素座標及元素節點參數 ... 18 2.7.2 移動元素座標的角速度、角加速度與元素節點

參數對時間的微分... 22 2.8 梁元素之變形描述... 27 2.8.1 梁元素之位置向量... 27 2.8.2 梁元素之應變... 28 2.8.3 梁元素應變之變分... 31 2.9 元素節點力之推導... 33 2.9.1 元素節點變形力向量之推導 ... 35 2.9.2 元素節點慣性力向量之推導 ... 39 2.10 元素剛度矩陣及質量矩陣 ... 44 2.11 與變形位置相關之節點作用力與負荷剛度矩陣 ... 50 2.12 系統平衡方程式與收斂準則 ... 52 第三章 數值計算方法與程序... 54 3.1 Newmark 直接積分法... 54 31.1 Newmark 直接積分法...54 31.2 Hilber-Hughes-Taylor scheme ...56 3.2 中央差分法... 57 第四章 數值例題與結果... 60 第五章 結論 ... 70 參考文獻 ... 72 附表 ... 76

附圖 ... 84 附錄 A ... 139 附錄 B ... 148 附錄 C ... 161 附錄 D ... 163 附錄 E... 172 附錄 F... 177 附錄 G ... 184

表目錄

表一 例題二之 I 型斷面性質... 76 表二 例題三之 I 型斷面性質... 77 表三 例題四之矩型斷面性質... 78 表四 例題五之 I 型斷面性質... 79 表五 例題六之矩型斷面性質... 80 表六 例題六質心位置分析... 81 表七 例題七之 I 型斷面性質... 82 表八 例題八之十字斷面性質... 83

圖目錄

圖 2.1 元素座標以及元素截面座標... 84 圖 2.2 固定元素座標與移動元素座標... 85 圖 2.3 旋轉向量 ... 86 圖 2.4 元素節點參數與節點力... 86 圖 2.5 決定x₂軸與x₃軸之第一步驟 ... 87 圖 2.6 決定x₂軸與x₃軸之第二步驟 ... 88 圖 2.7 剛性桿受力作用機制圖... 89 圖 4.1 例題一懸臂直角梁之幾何及受力圖... 90 圖 4.2 例題一 Newmark 積分法位移-時間曲線圖... 91 圖 4.3 例題一數值方法不同之 A、B 點X₃G方向位移比較 ... 92 圖 4.4 例題二 I 型斷面簡支梁之幾何及受力圖 ... 93 圖 4.5 例題二簡支梁受偏心軸力之負荷-位移曲線圖 ... 94 圖 4.6 例題二數值方法不同端點 B 之位移-時間曲線圖... 95 圖 4.7 例題二數值方法不同端點 C 之位移-時間曲線圖... 96 圖 4.8 例題二簡支梁受偏心軸力之位移-時間曲線圖(P₀ 50kip).... 97 圖 4.9 例題二簡支梁受偏心軸力之位移-時間曲線圖(P₀ 100kip) .. 98 圖 4.10 例題二簡支梁受偏心軸力之位移-時間曲線圖 (P 150kip)... 99

圖 4.11例題三I型斷面簡支梁之幾何及受力圖... 100 圖 4.12例題三簡支梁受力之負荷-位移曲線圖... 101 圖 4.13例題三簡支梁受力之位移-時間曲線圖... 102 圖 4.14例題四矩形斷面固端梁之幾何及受力圖... 103 圖 4.15例題四數值方法不同中點C之位移(Case (a)) ... 104 圖 4.16例題四數值方法不同中點C在XG₂方向的位移(Case (b)) .. 105 圖 4.17例題四數值方法不同中點C在XG₃ 方向的位移(Case (b)) .. 106 圖 4.18例題四數值方法不同中點C的轉角(Case (b))... 107 圖 4.19例題五直角懸臂梁之幾何及受力圖(Case (a))... 108 圖 4.20例題五直角懸臂梁之幾何及偏心受力圖(Case (b)) ... 109 圖 4.21例題五C點之負荷-位移曲線圖(Case(a))... 110 圖 4.22例題五C點之負荷-位移曲線圖(Case(b))... 111 圖 4.23例題五C點之位移(Case(a), P₀ 10kip, Q 0) ... 112 圖 4.24例題五C點之位移(Case(a), P₀ 10kip, Q0.001P)... 113 圖 4.25例題五C點之位移(Case(a), P₀ 10kip, Q 0.01P)... 114 圖 4.26例題五C點之位移(Case(a), P₀ 15kip, Q0) ... 115 圖 4.27例題五C點之位移(Case(a), P₀ 15kip, Q0.001P) ... 116 圖 4.28例題五C點之位移(Case(a), P₀ 15kip, Q 0.01P)... 117 圖 4.29例題五C點之位移(Case(b), P₀ 10kip)... 118

圖 4.30例題五C點之位移(Case(b), P₀ 15kip)... 119 圖 4.31 例題六三維矩形斷面梁之幾何及受力圖... 120 圖 4.32例題六E點之位移... 121 圖 4.33例題六質心在X₂G位移... 122 圖 4.34例題六梁在三維空間中之變形位置圖(t=0s-t=20s)... 123 圖 4.35例題七I型斷面梁中點承受一偏心力之結構圖... 124 圖 4.36例題七I型斷面梁中點承受偏心力之負荷-位移 曲線圖(e 0) ... 125 圖 4.37 例題七 I 型斷面梁中點承受偏心力之負荷-位移 曲線圖(e b/2) ... 126 圖 4.38 例題七中點 C 在X₂G方向的位移-時間曲線圖(P₀ 105lb). 127 圖 4.39 例題七中點 C 在X₃G方向的位移-時間曲線圖(P₀ 105lb). 128 圖 4.40 例題七中點 C 的轉角-時間曲線圖(P₀ 105lb)... 129 圖 4.41 例題七中點 C 在X₂G方向的位移-時間曲線圖 (P₀ 2105lb) ... 130 圖 4.42 例題七中點 C 在X₃G方向的位移-時間曲線圖 (P₀ 2105lb) ... 131 圖 4.43 例題七中點 C 的轉角-時間曲線圖(P₀ 2105lb)... 132 圖 4.44 例題八十字斷面梁之幾何及受力圖... 133

圖 4.45 例題八十字斷面梁中點之負荷-位移曲線圖 ... 134 圖 4.46 例題八十字斷面梁中點承受一偏心力之 負荷-位移曲線圖... 135 圖 4.47 例題八中點 C 在X₂G方向的位移-時間曲線圖... 136 圖 4.48 例題八中點 C 在X₃G方向的位移-時間曲線圖... 137 圖 4.49 例題八中點 C 的轉角-時間曲線圖... 138

第一章 緒論

三維梁在結構工程系統中，長久以來一直扮演著非常重要的角色，在 機械、航空太空、建築、車輛及土木工程中皆有很廣泛的應用。這些結構 在使用中，常伴隨著高速旋轉及大變形，所以在結構設計及分析時，必須 考慮其幾何非線性的動態反應。 常 見 的 幾 何 非 線 性 分 析 的 推 導 方 法 有 三 種 ： 全 拉 格 朗 日 法 (total

Lagrangian formulation)、更新拉格朗日法(updated Lagrangian formulation)和

共旋轉法(Co-rotational formulation)。全拉格朗日法是用初始狀態為參考位 置來表示總位移和旋轉；更新拉格朗日法是以結構上一個平衡狀態為參考 位置來表示增量位移和旋轉；共旋轉法是利用建立在元素當前變形位置的 元素座標將剛體位移及旋轉從總位移及旋轉中扣除，剩下的位移和旋轉即 為小位移和小旋轉，因此若使用共旋轉法，原本在線性分析的元素也可以 應用在大位移、小應變的幾何非線性分析，共旋轉法在梁結構的幾何非線 性分析已經被廣泛的使用[1-13]。 共旋轉法很少應用在非線性動態分析[7]，因為將元素運動分解成剛體 運動以及變形，導致慣性項變得相當複雜[7]，所以在動態分析，大部分的 文獻採用全拉格朗日法[14-18]，慣性項較為簡潔。

理慣性力項。文獻[1]使用 Consistent 2D Co-rotational formulation，該方法在 動 態 與 靜 態 皆 是 採 用 共 旋 轉 法 。 文 獻 [9] 將 文 獻 [1] 的 Consistent 2D Co-rotational formulation 延伸至 3D 梁結構，因為慣性項相當複雜，所以在 慣性矩陣省略一些次要的項[7]，以加強迭代的計算效率。 文獻[2,6]在梁元素當前位置建立固定元素座標系統，元素受擾動前及 擾動後變形皆在該固定元素座標系統上描述，再利用非線性梁理論、 d’Alembert 原理、虛功原理推導梁元素的節點變形力、節點慣性力，此方法 在推導上結合共旋轉法及全拉格朗日法，文獻[2,6]將此方法稱為共旋轉全

拉格蘭日推導法(Co-rotational Total Lagrangian(CRTL) formulation)。文獻[2]

的梁元素有兩個節點，十二個自由度，並無考慮翹曲的自由度，該梁元素 保留變形力至節點參數二次項，忽略變形對慣性力之影響，且不考慮轉換 矩陣對慣性力之影響，慣性矩陣僅考慮質量矩陣。文獻[6]推導一個雙對稱 斷面的開口薄壁梁元素，該梁元素有兩個節點，十四個自由度，考慮翹曲 的自由度。該梁元素保留變形力至節點參數二次項，忽略變形對慣性力之 影響，但考慮轉換矩陣對慣性力之影響，文獻[6]推導出的慣性力項相當簡 單。 文 獻 [5] 利 用 一 致 性 共 旋 轉 推 導 法 (Consistent Co-rotational(CCR) formulation) 推導一個不對稱斷面的開口薄壁梁元素，在移動元素座標上描 述元素的變形，在一個與當前的移動元素座標重合的固定元素座標上，推

導出元素節點在當前固定元素座標的擾動位移和擾動旋轉與元素節點旋轉 參數的擾動量的關係。文中利用元素節點在當前固定元素座標的位移和旋 轉及其擾動量、速度、加速度、角速度、角加速度，推導出移動元素座標 的角速度及角加速度。最後利用非線性梁理論、d’Alembert 原理和虛功原 理，推導出元素的節點變形力、節點慣性力，文獻[5]忽略變形及轉換矩陣 對慣性力之影響，並在推導慣性力時，將節點的虛位移直接表示成定義在 當前的固定元素座標上之節點參數的變分，以簡化慣性力的表示式，其慣 性矩陣僅考慮質量矩陣。 上述的文獻[5]與文獻[6]在推導三維梁元素，方法上主要的差異為一致 性共旋轉推導法(CCR) [5,9]在移動元素座標上描述元素的變形，共旋轉全拉 格蘭日推導法(CRTL) [2,6]則在當前變形位置建立的固定元素座標上描述元 素的變形。由於文獻[5]在移動元素座標上描述元素的變形，所以其絕對加 速度的表示式包含了定義在當前的固定元素座標及移動元素座標之節點參 數對時間的一次及二次微分，導致推導出之慣性力的表示式相當複雜，又 因為文獻[5]推導的梁元素為不對稱斷面，使其慣性力更為複雜，而文獻[6] 所推導出的慣性力則相對簡單。文獻[2,6]在推導上，使用之共旋轉全拉格 蘭日推導法，僅需保留到變形參數的二次項，不過其精神接近全拉格朗日 法。因為在文獻上非線性動態分析很少使用一致性共旋轉推導法，因此本 研究利用文獻[5]之一致性共旋轉推導法(CCR)推導一雙對稱斷面之開口薄

壁三維梁元素，探討三維梁之非線性動態行為，並與文獻[6]的梁元素比較 兩種推導法的慣性項之差異與動態分析時之效率及準確性。為了方便與文 獻[6]比較，本研究以文獻[5]的慣性力表示式為基礎，將元素的節點慣性力 表示成定義在當前的固定元素座標上之節點參數對時間微分的函數。由於 完整的慣性矩陣相當複雜，文獻[5]僅考慮對應於加速度的質量矩陣，本研 究還考慮了對應於速度耦合項的慣性矩陣，即陀螺矩陣，因文獻[6]中發現 考慮陀螺矩陣會改善平衡迭代的收斂速率。 第二章為理論推導，在推導時採用一致性共旋轉法，為了描述系統的 運動、元素的變形、邊界條件、以及與結構變形位置相關的外力，定義四 套座標系統，一致性共旋轉法在當前的移動元素座標定義元素的變形、節 點變形參數，在當前的固定元素座標定義元素的節點位移向量及旋轉向 量、節點位移向量及旋轉向量的擾動量，節點速度及加速度、節點角速度 及角加速度、節點力、剛度矩陣、質量矩陣、陀螺矩陣。第三章介紹本研 究求解非線性運動方程式的數值計算方法與程序，隱積分法使用 Newmark 直接積分法和 Newton-Raphson 增量迭代法來求解非線性動態平衡方程式， 顯積分法使用中央差分法。第四章以數值例題探討本文所推導之慣性力、 剛度矩陣、質量矩陣、陀螺矩陣之效率及準確性，並與文獻[6]比較結果之 差異，以及比較 Newmark 積分法與中央差分法效率及準確性。

第 二 章 理 論 推 導 本章採用一致性共旋轉法推導一個雙對稱斷面開口之薄壁梁元素， 梁斷面之形心與剪心重合，本章推導梁元素之節點變形力、慣性力、剛度 矩陣及慣性矩陣的推導方法和文獻[5]中所採用的方法大都相同。為了本文 的完整性，本章將會重複與文獻[5]中相同的部分。 2.1 基本假設 本文對非線性梁元素的推導，作如下的假設： (1)梁為細長的等斷面及雙對稱的開口薄壁梁，且Euler-Bernoulli 假說成 立。 (2)梁元素的形心軸之單位長度伸長量(unit extension)為均勻的伸長。 (3)梁元素在斷面平面上沒有變形。 (4)梁元素斷面的翹曲為梁元素的軸向扭轉率與該梁的聖維南 (Saint Venant)翹曲函數的乘積。 (5)梁元素的變形與應變皆為小變形與小應變。 2.2 座標系統

本文採用一致性共旋轉法(Consistent Co-rotational formulation)。為了 描述系統的運動、元素的變形、邊界條件、以及與結構變形位置相關的外 力(configuration dependent load)，本文中共定義了四套直角座標系統，並說 明如下：

(1) 固定總體座標系統(圖2.1)， ；系統的節點座標、方位、 位移、旋轉、速度、加速度、角速度、角加速度，系統的運動方程 式，其他座標系統之座標軸的方向餘弦，皆在此座標系統中定義。 3) 2, , 1 (i X_iG (2) 元素座標系統(圖2.1，圖2.2)，x_i(i1,2,3)、x_i(i1,2,3)；本研究採 用兩組元素座標， 一組為與元素一起剛體運動，但不一起變形的移 動元素座標 ，一組為固定在元素當前的變形位置之固定 元素座標 i x (i1,2,3) i x (i1,2,3)。兩組元素座標在當前梁元素變形的位置上是 重合的。此座標系統附屬在每一梁元素上，其原點位於該元素的節 點1上， 軸通過該元素的兩端節點(1、2，即兩端斷面的剪心)， 軸 與 軸在元素變形前與斷面的主軸方向一致，而元素變形後的 軸 與 軸，可由該元素未翹曲的兩端斷面的方位來決定[13]，即分別將 位於節點 1、2 變形後的斷面繞一個與該斷面之法線及 軸垂直的 旋轉軸旋轉一角度使斷面之法線方向與 軸方向一致(此時並不考慮 斷面之翹曲變形，否則斷面的法線方向無法定義)，然後再以兩斷面 主軸方向的角平分線作為 軸及 軸的方向。本研究在當前的移動 元素座標定義元素的變形、節點變形參數，在當前的固定元素座標 定義元素的節點位移向量及旋轉向量、節點位移向量及旋轉向量的 擾動量，節點速度及加速度、節點角速度及角加速度、節點力、剛 度矩陣、質量矩陣、陀螺矩陣。當梁元素在當前的變形位置受到擾 動時，擾動後的移動元素座標及對應的元素節點變形參數是由元素 當前的元素座標及節點變形參數、節點位移向量與旋轉向量、擾動 節點位移向量與旋轉向量決定。元素節點變形參數的擾動量是指元 素在擾動前及擾動後，移動元素座標上之節點變形參數的差。移動 元素座標的速度、加速度、角速度、角加速度是由當前的元素節點 速度、加速度、角速度、角加速度、節點變形參數決定。 1 x x₂ 2 x 3 x 3 x 1 x 1 x 3 x 2 x

(3) 元素斷面座標系統(圖2.1)， ；此座標系統與元素的斷 面一起平移和旋轉，其原點剛接於未翹曲斷面的形心上， 軸為未 翹曲斷面的法線方向， 軸與 軸分別與未翹曲斷面的主軸重合。 元素的變形是由斷面座標相對於元素座標的旋轉來決定。 3) 2, , 1 (i x_iS S x₃ S x₁ S x₂ (4) 負荷基底座標系統， ；此座標系統是用來描述與結構 變形位置相關的作用力機制。該作用力機制造成的系統節點外力及 力矩，和負荷剛度矩陣(load stiffness matrix)，皆建立在此座標系統 中。此座標系統的原點剛接於與結構變形位置相關的作用力作用的 節點上。 3) 2, , 1 (i X_iP 本文中用符號{}代表行矩陣。總體座標系統 與元 素座標 、元素斷面座標 、負荷基底座標 之關係可表示如下， } X , X , X { ₁G G_{2 3}G G _ X } x , x₂S ₃S } 3 x , x , {x = _{1 2} x } X , X , X₁P ₂P ₃P , x { ₁S S _ x { P _ X , x A XG _{ GE} , S GS G _{A x} X  XG A_GPXP, (2.1) 其中 、 、 分別代表元素座標、元素斷面座標、負荷基底座標 相對於固定總體座標系統的方向餘弦矩陣。 GE A A_GS A_GP 2.3 旋轉向量及其對時間的微分 本文中使用旋轉向量來表示一個有限旋轉，如圖2.3所示，一向量b 受 到一旋轉向量



e的作用而轉到一個新的位置b~，向量b~與b 之間的關係 可表示成[19]

) ( sin ) )( cos 1 ( cos ~ _{b e b e e b} b



 



 



 b I e e I e I sin ( ) 1 cos ( )] [     ₂   









=_{R b} (2.2)  其中I 為33的單位矩陣，表向量外積，．表向量內積，



表繞旋轉軸的 旋轉角，e 表旋轉軸的單位向量，R_稱為旋轉矩陣。 1  當



時，sin 1





， 2 1 cos 1 2  





，R_的近似式可表示成 2 2 1 Φ Φ I R_    (2.3)                 0 0 0 1 2 1 3 2 3













I Φ  (2.4) 其中Φ 為對應於旋轉向量 的反對稱矩陣 由(2.2)式 b~對時間的微分可表示成 b R R b R b ~ ~ t dt d        (2.5) 因 b~的長度固定，所以其對時間的微分可表示成 b ω b ~ ~   dt d (2.6) 其中ω 為角速度向量。 從(2.5)式和(2.6)式可得 ωIR_R_t (2.7) 因ωI為一反對稱矩陣，所以 t 也為反對稱矩陣。  R R 由(2.2)及(2.6)式可得ω和旋轉向量 對時間的一次微分有以下關係[3]：

ωΓ() )] ( ) ( [ ) ( I ₁ I ₁ I Γ   a  b  (2.8) 其中 2 1 cos 1





  a ， 1 (1 sin ) 2 1

_





  b ，當



0， 2 1 1  a ， 6 1 1  b 。 當旋轉向量 有一微小變量 



時會使向量 b~繞x_i (i 1, 2,3)軸做微小旋轉 i



，



與



的關係和 與 的關係相同[3]，即  ω



Γ()



 (2.9) 將(2.8)式對時間微分可得 ω Γ()Γ() (2.10) ) ) ( ) ( ) ( ) ( ₁ I ₁ I ₁ I Γ  a_  b  a  ( ) ( ₁ 1  I   I b  b  其中ω 為角加速度， 為旋轉向量對時間的兩次微分。 _  由(2.8)、(2.10)式可知當 0時，Γ(0)I，Γ )( 0所以 ω  (2.11) ω  (2.12) 即在旋轉向量0時，其對時間的一次、二次微分之值，等於角速度、角 加速度。 2.4 梁之剪心軸的位移及其斷面的旋轉 本文是在當前的移動元素座標x_i(i1,2,3)上，描述梁元素當前的的 變形，由2.1節中的基本假設可知，梁元素的變形可由其剪心軸在移動元素 座標上的位移及其斷面繞剪心軸的旋轉決定。本文以梁元素之剪心軸為參 考軸、以梁元素兩端斷面的剪心為節點，節點自由度及廣義力都在剪心上

定義，所以本文在描述梁元素的變形前，先描述其剪心軸的位移及其斷面 的旋轉。 本文中若( )表示定義於移動元素座標中的變數，則( )表示定義於固定 移動元素座標中的對應變數，在當前梁元素變形的位置上，固定元素座標 與移動元素座標是重合的，所以定義在兩座標系統的對應變數有相同的 值，但其擾動量及對時間的微分並不相同。 令梁剪心軸上的任一點P(見圖2.1) ，在不計扭轉翹曲時，變形前後的 位置向量分別為{x,0,0} 及{x_p(x,t),v(x,t),w(x,t)} ，其中v( tx, )與w( tx, )為 剪心軸的側向位移。 本文中符號( ) 代表( )_,x () x。變形後剪心軸的單位切線向 量可表示為 } , , {cos





₃ 



₂  _n t (2.13) 2 1 2 3 2 2 ) 1 ( ) , ( cos



 







   s t x x_p n (2.14) o w s t x w





        1 ) , ( 2 (2.15) o v s t x v





      1 ) , ( 3 (2.16) x w w     ， x v v     ， 1    x s o



(2.17) 其中



_n為元素座標之x₁軸和 的夾角，t s_{為節點1至點P間的剪心軸在不計} 扭轉翹曲時，變形後的弧長，



_o為剪心軸的單位伸長量(unit extension)。

由梁剪心軸上的任一點P在變形後的位置向量及(2.14)至(2.17)式，在



o、 與 都遠小於1時， 可以表示成 x v_, x w_, x_p( tx, ) dx w v t x x x(1



_o  _{x x} ) , 0 ( t p p 



₀  2 , 2 , ) 2 1 2 1 ) , ( (2.18) 由(2.18)式可知 為節點1在 方向上的位移，由移動元素座標系 統的定義方式可知 、其擾動量 1 x u  1 u 1 x 1 u



及對時間的微分都為零。 由梁之形心軸單位長度的伸長量為均勻的伸長量之假設及(2.18)式， 可以得到形心軸單位長度伸長量



_o可表示如下 



   _{x x} o w dx L v L L L 2 , 2 , ) 



0 ( 2 1 (2.19)  x_p L t 



L  _o  v_x  w_x dx 0 2 , 2 , ) 2 1 2 1 1 ( ) , (



 (2.20) 其中 L 為梁元素變形前的長度， 為梁元素變形後之形心軸的弦長， 由移動元素座標的定義方式節點2在 方向的位移x₁ u₂可以表示成 u₂ L (2.21) 本 文 用 元 素 斷 面 座 標 軸 的 旋 轉 表 示 梁 之 斷 面 的 旋 轉 。 令e 與 分別代表移動元素座標的 i 3) 2, , 1 i e (S i x_i (i1,2,3) , 1 (i S i e 3) 2, 方向的單位向量與元素 斷面座標的 軸方向的單位向量。由座標系統的定義方式可 知，在變形前 軸與 軸的方向是一致的，而且變形後 與(2.13)式的t 方向一樣。在本文中假設變形後的單位向量 的方向是由以下兩 個旋轉向量連續作用於單位向量 S i x x ) 3 , 2 , 1 (i i xiS e1S 3) 2, , 1 (i i e 來決定： θn 



nn (2.22) θ_t 



₁t, (2.23)

n{0,



₂ (



₂2 



₃2)1/2,



₃ (



₂2 



₃2)1/2}{0,n₂,n₃}, (2.24) 其中n 為垂直於e₁與 之單位向量，t



_n定義於(2.14)式， 定義於(2.13)式，t 1



為斷面繞 旋轉的角度。 t 旋轉向量θ_n作用在 上，將其轉至一中繼位置e_i e_i (i1,2,3)，此時e₁與 重合，再將θ 作用在 ，將其轉到 。若 、 、以及 已知，則元素 斷面座標 就唯一決定；反之，若e_i S i e 知，則旋轉向量θ_n θ 亦唯一_t 決定。 t _t e_i S i e 與 i e θ_n θ_t 與 S i e 已 由(2.2)、(2.13)、(2.14)式與(2.22)－(2.24)式， 與 之關 係可表示如下 S i e e_i (i1,2,3) (2.25) i i S i t R R e R e e [ _{1 2}]  _ 2 1 1 1 1 cos r sin r R 







2 1 1 1 2 sin r cos r R 







} ) cos 1 ( , ) cos 1 ( cos , { _{3 2}2 _{2 3} 1 





n  



n n 



n n n r } ) cos 1 ( cos , ) cos 1 ( , { _{2 2 3 3}2 2 







n n n



n  



n n r 其中R_稱為旋轉矩陣。因R_為



_i(i1,2,3)的函數，所以本文中稱



_i為旋 轉參數。

本文中假設梁元素變形後的剪心軸，其側向位移v( tx, )與w( tx, )及剪 力中心軸的扭轉角



₁(x,t)皆為 x 的三次Hermitian氏多項式。因此(2.17)式之 ) , ( tx v 與w( tx, )可表成 (2.26) b t b t _{v v v v} N N N N t x v( , ){ ₁, ₂, ₃, ₄}{ ₁, ₁, ₂, ₂}N u (2.27) c t c t w w w w N N N N t x w( , ){ ₁, ₂, ₃, ₄}{ ₁, ₁, ₂, ₂}N u



₁(x,t){N₁,N₂,N₃,N₄}t





₁₁,



₁,



₁₂,



₂



Nt_du_d (2.28) ), 1 )( 1 ( 8 L ), 2 ( ) 1 ( 4 1 ), 1 )( 1 ( 8 L ), 2 ( ) 1 ( 4 1 2 4 2 3 2 2 2 1                      N N N N (2.29) L x 2 1  



(2.30) 其中v_j與w_j (j1,2)分別是 v 與 在節點w j 的節點值， 及 則是(2.17) 式中v 及 在節點 j v w_j  w j ( j 1,2)之節點值，



_1j(j1,2)是(2.23)式之



₁在節點 j 的節點值，



_j(j1,2)是



_1,x在節點 j 的節點值。N_i (i14) w  代表形狀函 數(shape function)。由移動座標的定義方式可知v_j與 _j(j 1,2)，其變分 及其對時間的微分皆為零。 2.5 節點參數與節點力 本文中用旋轉向量描述梁元素兩端節點之斷面的有限旋轉，但用旋 轉參數描述梁元素在兩端節點及內部之斷面的有限旋轉，因對應於旋轉向

量及旋轉參數的廣義力矩不一樣且非向量，所以不同元素在共同節點的廣 義力矩不能以向量的方式相加，傳統力矩為向量，所以本文將元素節點的 廣義力矩轉換成等效的傳統力矩，使其能以向量的方式相加，因對應於傳 統力矩向量的廣義位移為繞該力矩向量的微小旋轉，所以本文在推導梁元 素時需使用以下七類元素節點參數： （1）u_ij (i1,2,3; j1,2)，如圖2.4所示，u_ij (u_1j u_j，u_2j  ，v_j j j w u₃  ) 為 元 素 節 點 j 的 位 移 向 量 u 在 其 當 前 的 固 定 元 素 座 標 軸_j i x (i1,2,3)方向的分量，由固定元素座標系統的定義方式可知u 中除了_ij 12 u 外，其餘的值皆為零，但u 的增量_ij  、u_ij u 的擾動量_ij



u_ij、及u 對時_ij 間的微分u_ij、u_ij並不為零。對應於



u_ij的廣義節點力 f_ij，為在x_i軸方向的 力。 （2）



_ij (i1,2,3; j1,2)，



_ij是元素節點 j 繞其當前的固定元素 座標軸x_i (i1,2,3)的擾動旋轉，對應於



_ij的廣義節點力m ，為繞_ij x_i軸的 傳統力矩(見圖 2.4)。 （3）



_ij (i1,2,3; j1,2)，



_ij是元素節點 j 的旋轉向量 在其固定 元素座標軸 j  i x (i1,2,3) j  方向的分量，本研究中，在任何時刻及位置都將節 點的旋轉向量 的值重新設定為零，但



_ij的增量 、擾動量



_ij



_ij及對時間 的微分並不為零。對應於



_ij的廣義節點力m 為一廣義力矩。因_ij 的值重 新設定為零，由(2.8)及(2.9)式可知 j  ij



和



_ij的值相同，所以廣義力矩m 和_ij 傳統力矩m 的值相同，且可視為是向量，但對應於_ij



_ij和



_ij的切線剛度 矩陣並不相同。 （4）u_ij (i1,2,3; j1,2)，u_ij (u_1j u_j，u_2j v_j， )為元素 節點 j j w u₃  j 的位移向量u_j在其當前的移動元素座標軸x_i (i1,2,3)方向的分量，

由移動元素座標系統的定義方式可知除了 外， _ij、 _ij的增量 _ij、u_ij 的擾動量 u_ij 12 u u u u



及u 對時間的微分_ij u_ij、u_ij的值皆為零。對應於



u_ij的廣義節 點力 f ，為在_ij x 軸方向的力。_i （5）



_ij (i1,2,3; j1,2)，



_ij是元素旋轉參數



_i在節點 的值，j



₁定 義於(2.23)式，



₂與



₃定義於(2.15)與(2.16)式，



_ij是用來描述梁元素在兩端 節點及內部之斷面的有限旋轉。因元素變形後



_ij的值不為零，所以



_ij並 不是繞x_i軸的無限小旋轉。 (j （6）



_j 1,2)，如圖 2.4 所示，



_j是梁元素之剪心軸的扭轉率 x   /₁ 





在元素節點 j 的值。對應於



_j之擾動量



_j的廣義節點力 為雙 力矩(Bimoment)。 j B j



及



_j與標系統無關。 （7）



_ij* (i1,2,3; j 1,2)，



₁*_j 



_1j、



₂*_j w_j、



₃*_j v_j，其中



_1j 在節點 j 的節點值， ₁ 是元素旋轉參數



₁



定義於(2.23)式， 與 為(2.17) 式中v 與 在節點 j v w_j  w j 的節點值。本文中採用



_ij*及



_j來決定梁元素剪心軸的 側向位移及軸向扭轉角。對應於



_ij*之擾動量



_ij*的節點力為廣義力矩m 。_ij 因為 在變形後不為零，所以 並不是繞固定元素座標 軸的無限小旋 轉，所以廣義力矩 * ij



* ij



x_i  ij m 並非繞固定元素座標x_i軸的傳統力矩。 本文中在元素節點 j (j1,2)對應於元素節點參數u_ij、



_ij、



_ij及 j



(i1,2,3)的系統節點參數為u_ijG、



_ijG、



_ijG及



_j，其中u 是節點_ijG j 的 位移向量uG_j 在固定總體座標軸X_iG(i1,2,3)方向的分量；



_ijG是元素節點 j 的旋轉向量 在固定總體座標軸 方向的分量。本研究中，在 任何時刻及位置都將節點的旋轉向量 的值重新設定為零，但 的增量 、擾動量 及對時間的微分並不為零。本文用元素節點 G j G  G(i1,2,3) G j  X_i G ij



G ij







_ij j ( j1,2) 的

位移增量 及 決定當前的元素斷面座標、移動元素座標、元素節點 變形位移 及旋轉參數 G ij u  ij u G ij



 ij



(i1,2,3; j1,2) G ij



  j (詳見2.6節)，因本文中 的 值重新設定為零，所以由(2.9)式可知 和 的值相同。對應於 的 廣義節點力，為繞 軸的傳統力矩；對應於 的擾動量 的廣義節點 力為在 方向的力；因 G j  G ij



G ij



G ij u G i X



u_ijG G i X



與座標系統無關，對應於



_j的廣義節點力， 亦為廣義雙力矩 (Bimoment)B_j。u_ijG、



_ijG與 、u_ij  的關係可以由



_ij 標準的座標轉換求得。 2.6 元素斷面座標、固定元素座標、元素節點位移及旋轉參數之決定 本文是使用增量迭代法(incremental iterative method)解非線性平衡方 程式。假設第 I 個位置為已知，此處所謂的第 I 個位置，是指第 I 個增量的 平衡位置，或是指某一增量中第 I 次迭代後的位置。本節中所有的向量都是 表示成在總體座標的分量。令IXG_j 、 S ij I x 、Ie (_ijS i1,2,3; j1,2)分別是在 第 I 個變形位置時，元素節點j_{在總體座標的位置向量、元素斷面座標與元} 素斷面座標軸的單位向量。令uG_j 與G_j (j 1,2)_{分別代表元素節點 在總}j 體座標系統中的增量(或改正)位移向量與增量(或改正)旋轉向量。當元素節 點 j (j 1,2)受到 及 作用後，其當前的座標 可由 加上 得到，當前的元素弦長_(即兩元素節點的距離)及當前的固定元素座標之 G j u  G j   G j X I G j X G j u  1 x 軸可以由XG_j 求得，u₁₂ u₂的值可以由(2.21)式求得。當前的元素斷面座標 S ij x 之座標軸的單位向量e 是將旋轉向量_ijS G_j 作用在Ie 得到。 _ijS 本文以下列方法決定在當前變形位置的固定元素座標[13]。

由 通 過 當 前 元 素 兩 節 點 的 軸 決 定 當 前 固 定 元 素 座 標 的 x₁軸 。 令 nj θ (j 1,2)_{為(2.22)式中旋轉向量}θn在節點 j 之節點值，由 的定義方式， 可以表示成 n θ nj θ j nj nj n θ 



_(2.31) S j j nj 1 1 sin



n e e 其中e₁為x₁軸的單位向量，e 為₁S_j S j x₁ (j 1,2)軸的單位向量。 當前固定元素座標的x₂及x₃軸可由下列兩個步驟決定： 步驟1：將旋轉向量θnj作用在x 軸上(圖2.5)，使其旋轉到標示為ijS S ij x (圖2.6) 的新位置。此時，x₁S_j軸與x1軸重合，且 S j x₂和x₃S_j軸與x1軸垂直。 步驟2：定義固定元素座標之xi(i2,3)軸的單位向量ei (i2,3)如下 ' 2 ' 1 ' 2 ' 1 S i S i S i S i i e e e e e    _(2.32) 其中e_ijS'(i2,3; j 1,2)為x 軸的單位向量。 _ijS' 在當前的變形位置及固定元素座標上，元素節點的旋轉參數，可依下 列的方法決定。令



1je1(j 1, 2)表示將x2軸旋轉到 S j x₂軸的旋轉向量，則由 (2.23)式可知



1j為節點旋轉參數，本文中以下式計算



1j ) ( sin 1 _{2 2}' ₁ 1  e e e  S j j



(2.33) 由(2.22)、(2.24)、以及(2.31)式，節點旋轉參數



_ij(i2,3; j=1,2)可用 下式求得 i S j ij e1e1 e



(2.34)

移動元素座標與固定元素座標在梁元素當前的變形位置上是重合的， 但定義在固定元素座標的元素節點參數受到擾動時，移動元素座標就不再 與固定元素座標重合，且定義在其上的元素節點參數亦有一對應的擾動， 本節中將探討擾動後之移動元素座標與固定元素座標的關係及元素節點參 數的改變量。本文中移動元素座標對固定元素座標的速度及加速度是由定 義在固定元素座標的節點速度及加速度決定，本節中亦將推導其關係。 2.7.1 擾動後之移動元素座標及元素節點參數 令 t 表示當前的時間，xi與xi (i1,2,3)表示在當前變形位置重合的固 定元素座標與移動元素座標，在當前變形位置的固定元素座標與移動元素 座標是重合的。本節中所有向量的分量都是在元素當前的變形位置之固定 元 素 座 標 的 分 量 。 令 u_j  u_j(t){u_1j,u_2j,u_3j}{u_j,v_j, w_j} 、 } , , { ) ( _{1j 2j 3j} j j   t 







 (j 1,2)分別表示元素節點 j 當前的位移向量、旋 轉向量，



_ij 



_ij(t) (i1,2,3; j1,2)及(t)分別表示元素節點 j 當前的 旋轉參數及元素當前的弦長。本文中旋轉向量 的值在當前變形位置重新j 設定為零，u 的分量除了j u12(即u2)外其餘的值皆為零，但因推導移動元素 座標的速度、角速度、加速度及角加速度過程的需要，故仍保留全部的 及j j u 。 令



u_j與



_j (j 1,2)分別代表元素節點 j 的位移向量u 與旋轉向量_j _j 的擾動量， t



代表元素節點 j 受擾動的時間。當元素節點 j 受到



u_j 與 j 



(j1,2)作用後，可由 2.5 節的方法及一致性一階線性化(consistent first order linearization) 決 定 受 擾 動 後 的 移 動 元 素 座 標 、 元 素 節 點 旋 轉 參 數 ) t (t ij





 (i1,2,3; j1,2)、元素弦長_(t 



t)，其推導過程詳見附錄 A。 令r(t)、r(t)是為梁元素變形後，梁元素中任意點在當前的固定及移動

元素座標上的位置向量。梁元素受擾動後，r(t



t)與r(t



t)為該點在當 前的固定元素座標及擾動後之移動元素座標的移動元素上的位置向量。由 附錄 A 可知位置向量r(t)與r(t)及r(t



t)與r(t



t)有如下的關係 ) ( ) ( ) ( ) (t 0A_xx t r t u₁ t r   (2.35) ) (t



t ) ( ) ( ) ( ) (t



t 1A_xx t



t 0A_xx t r t



t u₁ r (2.36) x x x x x x A A A 10 11 1 _{ } (2.37)                           1 2 2 1 1 12 11 12 11 0









    w v w v x x A (2.38)                           1 2 2 1 1 12 11 12 11 10

















    w v w v x x A (2.39)                                0 0 0 1 1 2 1 1 2 1 2 2 11 A v u w A w u v v u w u v m m m x x                A      1 w m  (2.40) 其中0A_xx0A_xx(t)，1A_xx A1 _xx(t



t)分別為擾動前的移動元素座標相對於 固定元素座標的方向餘弦矩陣及擾動後的移動元素座標相對於擾動前的元 素座標的方向餘弦矩陣， ) ( 4 1 ) 4 1 ) ( 4 1 ) ( 4 1 32 31 22 32 31 22 21 1     (₂₁                     w w v A  v

) ( 4 1 ) ( 4 1 ) ( 4 1 ) ( 4 1 32 22 31 21 22 32 21 31                             v v w w 2 2 ₂ 2_{ } v w w v  _   





， 2 12 11 1







_m   ， 2 12 11 1







_m   。 擾動後元素弦長可表示成(詳見附錄A) w w v v u t t



















   (  )    (2.41) 1 2 u u u    ，v v₂ v₁，w w₂ w₁ (2.42) 1 2 u u u







   ，



v 



v₂ 



v₁，



w 



w₂ 



w₁ (2.43) 擾動後元素節點 j 的旋轉參數



ˆ_ij(t



t) (i1,2,3; j1,2)在擾動後的 移動元素座標上可表示成(詳見附錄A) 11 12 ˆ ˆ





 (2.44) ( ) 4 1 ) ( 4 1 2 2 2 2 2 22 21 32 31 12 11 12 11 11 12





























 _{       }    w v ) ( 4 1 4 1 4 1 4 1 4 1 21 22 32 22 31 21 22 32 21 31





















        v ( ) 4 1 ) ( 4 1 ) ( 4 1 31 32 21 22 31 32

















            w v w j j j j j j j w w 1 3 3 12 11 2 2 2 2 2 ˆ

























          ( ) 2 ) ( 2 3 12 11 3 12 11 j j v v

_









_







         j j j j u w 3 1 1 3 2 ₂ 1 2 1

_

_

_

_



      j j j j j j j 3 3 3 11 12 2 2 1 3 2 ˆ





v







v















         

) ( 2 ) ( 2 2 12 11 2 12 11 j j w w

_









_







 _{ }  _    j j j j u v 2 1 1 2 2 ₂ 1 2 1

_

_

_

_



 _{ }    (2.44)式中加底線之項可視為在擾動前之元素座標上之節點旋轉參數 。本文中元素節點旋轉參數的擾動量 ) ( ˆ t_ij





_ij(i 1,2,3)(j1,2)為擾動後 元素節點旋轉參數



ˆ_ij(t



t)與擾動前元素節點旋轉參數



ˆ_ij(t)的差，即 ) ( ˆ ) ( ˆ _{t t t} ij ij ij









   (2.45) 詳見附錄A之(A. 38)－(A. 43)式。 本文中元素的弦長及元素節點旋轉參數的擾動量指擾動前後之元素 弦長的差及節點旋轉參數的差。由(2.41)式可得 w w v v u t t





















     (  )    (2.46) 由移動元素座標 的定義可知元素節點位移在 方向的擾 動可表示成 i x (i1,2,3) x_i 0 1  vj  wj  u







_(2.47) 由 (2.21)、(2.41)及(2.46)式可得 u   u  v v  ww   











₂ (2.48) 由 (2.15)、(2.16)式可得 o j j o j j



w











* _{2 2} 2   (1 )  (2.49) o j j o j j



v











* _{3 3} 3   (1 )  其中



_o為剪心軸的單位伸長量，定義於(2.19)式，



ij (i2,3; j 1,2)為 (2.45)式之元素旋轉參數



_i在元素節點 j的擾動量。

2.7.2 移動元素座標的角速度、角加速度與元素節點參數對時間的微分 本文中(˙)代表()對時間的微分。令uj {u1j,u2j,u3j}{uj,vj,w j}、 } , , } , , { _{1j 2j 3j j j} j u u u v w    u {u_j 、  _j {



_1j,



_2j,



_3j} 與  _j {



_1j,



_2j,



_3j} ) 2 , 1 (j  _{分別代表在時間 t，元素節點} j在當前的固定元素座標中的絕對速 度、絕對加速度、絕對角速度、絕對角加速度。 令r(t



t)r(t)



r及u₁(t



t)u₁(t)



u₁，由(2.35)－(2.40)式，梁元 素中任意點在當前的固定元素座標上的絕對速度可表示成 ) ( ) ( ) ( ) ( 1 0 0 ₁ 0

lim

t t t t t t _{xx xx xx} t u r A r A A r r r         





 (2.50) x x x x x x x x t x x t t t I _{A A Ω A} A A 3 10  11 11 1 0 1 ( )

lim

       





 (2.51)                                                0 2 2 0 0 0 0 0 12 11 1 2 12 11 1 2 1 2 1 2 10









                 w w v v w w v v x y x z y z x x Ω A (2.52)                                         0 0 0 1 2 1 2 2 2 11 A v u w A w u v v u w w u v x x x x               _xx A (2.53) 2 12 11





    _x ,   ₂ w₁ w y     ,   ₂ v₁ v z    (2.54) 1 2 u u u      (2.55) ) ( 4 1 ) ( 4 1 ) ( 4 1 ) ( 4 1 32 31 22 21 32 31 22 21 1











 







 







     v w _{z y} A       

) ( 4 1 ) ( 4 1 ) ( 4 1 ) ( 4 1 32 22 31 21 22 32 21 31                         _{y y z z}   2 2 z y w v       其中(2. 50)式中u 為元素節點1 1的絕對速度，亦為移動元素座標的速度，(2.

51)式中之1A_xx為移動元素座標的角速度矩陣(Angular velocity matrix)，因本

文中不考慮隨元素減小時會趨近於零的變形參數，角速度矩陣僅考慮(2.52) 式中之Ω，故(2.51)式中之x、y及z可視為移動元素座標對固定元素座 標之x_i軸的角速度。因推導移動元素座標的角加速度的需要，(2.53)式中之 x x A 11 _{仍需保留， 梁元素中任意點在當前的移動元素座標上的速度。r} 由移動元素座標 的定義，可知元素節點在 方向對移動元 素座標的速度中 i x (i1,2,3) x_i 0 1 vj wj  u   (j1,2) (2.56) 由(2.21) 及(2.41)式可得 w v u t t t t u _{z y} t                  



_

 ) ( ) (

lim

0 2 (2.57) 由(2.45)、(A.38)－(A.43)式及 t ij t ij

_







lim

0    (i1,2,3; j1,2)可得 ) ( 4 1 ) ( 4 1 2 2 22 21 32 31 12 11 11 12

















                 v w (2.58) ) ( 4 1 ) ( 4 1 4 1 4 1 4 1 4 1 31 32 21 22 32 22 31 21 22 32 21 31

























                    v w ( ) 4 1 ) ( 4 1 31 32 21 22









          w v j j j j j w 1 3 3 12 11 2 2 2

















            j j j j j j u w v v 3 1 1 3 2 3 12 11 3 12 11 2 1 2 1 ) ( 2 ) ( 2





















_{ }                       

j j j j j v 1 2 2 12 11 3 3 2

















            j j j j j j u v w w 2 1 1 2 2 2 12 11 2 12 11 2 1 2 1 ) ( 2 ) ( 2





















_{ }                      11 31 31 12 11 21 21 2

















        w     31 11 11 31 2 31 12 11 31 12 11 2 1 2 1 ) ( 2 ) ( 2





















_{ }                        v v w u 12 32 32 12 11 22 22 2

















       _ w _{ }  _ 32 12 12 32 2 32 12 11 32 12 11 2 1 2 1 ) ( 2 ) ( 2





















_{ }                        v v w u 11 21 21 12 11 31 31 2

















        v     21 11 11 21 2 21 12 11 21 12 11 2 1 2 1 ) ( 2 ) ( 2





















_{ }                      w w v u 12 22 22 12 11 32 32 2

















        v     22 12 12 22 2 22 12 11 22 12 11 2 1 2 1 ) ( 2 ) ( 2





















_{ }                      w w v u 將(2.15)、(2.16)式對時間微分並取其在節點 j (j1,2)之值可得 o j j o j j w











  2 2  * 2   (1 )  (2.59) o j j o j j v











* _ _{3 3 } 3   (1 )  其中



_ij (i1,2,3; j1,2)已在(2.58)式中定義。 梁元素中任意點在當前的固定元素座標上的絕對加速度可由(2.50)式 對時間微分求得，若不考慮隨元素減小時會趨近於零的變形參數，(2.50)式 對時間微分可表示成