2-1 彈性波傳理論
2-1-1 彈性波傳簡介
波傳可視為一種能量的傳遞,而能量向外傳遞必需仰賴介質,當 傳遞的介質為彈性時,此時即稱為彈性波,彈性波又因其傳遞路徑的 不同,區分為實體波(body wave)和表面波(surface wave)兩種形 式。
實體波在物體內部和表面傳播,又依波傳遞方向與介質振動方向 的關係而分為壓力波(Primary wave)及剪力波(Secondary wave);
壓力波通過物體時傳播介質的振動方向和波傳方向一致或平行,如圖 2-1。剪力波通過物體時,傳播介質的振動方向和波傳方向垂直,如 圖2-2 所示。
表 面 波 運 動 範 圍 主 要 在 自 由 表 面 附 近 , 又 稱 為 雷 利 波 (Rayleigh wave) ,其傳播介質的振動方式是在垂直面上的橢圓運 動,如圖 2-3。由於表面波的傳遞速度以及傳遞時能量隨距離的 衰減均較實體波慢,因此常被視為造成地表振動的主要來源之 一,於相關的研究中相當受到重視。
2-1-2 彈性波傳控制方程
可用兩個彈性常數(Lame’s constants)λ、μ 表示。
⎥⎥
) G
λ、μ= Lame’s constants E= Young’s Modulus G= Shear modulus υ= Poisson ratio
σ=作用於介質元素上的應力
⎥⎥
2
2
其中 ρ 解,即為雷利波。雷利波為Lord Rayleigh【2】首先提出,其後 Lamb
【3】曾詳細描述其性質,又因其傳遞範圍受限於半無限域自由表面
附近,而稱為表面波。
(2-14b)對 z 做微分,再進行相加可得
將式(2-19a、b) 代入式(2-15)、(2-17)控制方程,可得
= ,則式(2-20a)、(2-20b)可寫為 F(z)、
G(z)之二階常微分方程式
) exp(
A
G(z)= 2 −sz (2-23b) 將式(2-23a)、(2-23b)代回式(2-19a)、(2-19b),於是勢能函數Φ、Ψ最 後可表示成
整理得
數解,即可求出雷利波波速νR =K νS。
比較上述提到的三種波波速,壓力波波速最快,其次為剪力波,
雷利波則最慢;一般而言,三種波的波速為有條件之比例關係,此關 係視土壤之柏松比(Poisson ratio)而定,如圖 2-6 所示。
就傳遞能量而言,Miller 和 Pursey【4】對一均勻分佈之圓形基礎,
在一均質、等向、彈性半無限空間地表上,以垂直方向震動所得之能 量分配結果顯示,三種彈性波的能量分配為雷利波 67%、剪力波 26%、壓力波 7%,如圖 2-7。
由於雷利波的能量佔總能量的67%,而且其傳遞主要集中於介質 表面,是以在考慮地表振動時,雷利波將是一項很重要的影響因素。
2-2 有限元素法理論
2-2-1 有限元素法分析原理
有限元素法是一種新興分析結構問題的數值計算方法,這個名稱 是由美國加州大學教授克拉夫Clough,R.W.【5】所首次引用。
有限元素法原理是將原本複雜、龐大的彈性連續體結構,分割成各個 簡單的計算元素或單元(Element),這個人為地劃分元素的過程稱為彈
性連續體的離散化,若是把離散化後各個元素組合成計算模型,並假 設元素與元素之間僅在有限的節點(Node Point)上連接,來近似的替 代原來的彈性連續體,如此就可以將原來分析彈性連續體上無限個質 點變形和應力的問題,轉化成分析離散化計算模型中有限個節點變形 和應力的問題。
理論上來說,當元素劃分的越小、元素和節點的數目越多,則離 散化模型越接近真實彈性連續體結構,計算分析的結果也與實際越接 近,對於極限狀況,當將元素劃分無限小、元素和節點的數目無限多 時,則離散化模型就變成與彈性連續體完全相同,計算分析的結果就 與實際完全一樣了。但對於真實問題而言,元素和節點的數目總是有 限的,因而此種方法稱之為有限元素法。
2-2-2 靜力學的有限元素法分析
{ }
α ] [ } {q = C求解出待定係數得
{ }
q C] 1 [ }{α = − (2-33)
{ }
q = [u1v1 w1u2 v2 w2"] T為元素的節點位移將式(2-33)代回式(2-32)得
{ }
q N{ }
q CS
d} [ ][ ] [ ]
{ = −1 = (2-34) 其中[N]稱為形狀函數矩陣 式(2-34)即是用節點位移
{ }
q 來表示元素體內任意點位移函數{d}的關 係式。3. 根據幾何方程式由位移函數求應變{ε}和應變矩陣[B]
將式(2-34)以節點位移表示之u、v、w,代入幾何方程關係式(2-3) 中,計算元素的應變
{ }
ε ,把結果寫成矩陣形式,則有{ }
ε =[B(x,y,z)]{ }
q (2-35)其中[B]稱為元素的應變矩陣
4. 依據物理方程,透過應變{ε}求應力{σ}
將式(2-35)代入應力應變的關係式(2-1)中,計算元素的應力
{ }
σ 得{ }
σ =[D]{ε}=[D][B]{ }
q (2-36) 5. 以虛功方程式求出元素勁度矩陣就靜力學而言,建立元素勁度矩陣常用虛功方程代替平衡微分方
程,而三維問題的虛功方程為
{ }
q{ } { } { }
dxdydzV ε* T σ
* T F =
∫∫∫
(2-37) 式中{ }
F =作用於元素各節點的外力{ }
q* =元素各節點處的虛位移{ }
σ =外力作用下元素任意點上的應力{ }
ε* =元素任意點上與虛位移相對應的虛應變dx、dy、dz為微小元素體的邊長,它們的乘積表示微小元素體的體 積dV,由於
{ }
ε* 與{ }
q* 相協調,它們之間也可以寫出與式(2-35)類似的 關係式{ }
ε* =[B]{ }
q*{ } { }
ε* T = q* T[B]T (2-38) 將式 (2-38)代入式(2-37),虛功方程可改寫成{ }
q* T{ }
F =∫∫∫
V{ }
q* T[B]T{ }
σ dV ={ }
q* T∫∫∫
V [ B]T{ }
σ dV由於虛位移為任意,所以可以由等號兩端消去
{ }
q* T得{ }
F =∫∫∫
V [ B]T{ }
σ dV (2-39)最後將式(2-36)代入式(2-39)中,可得元素勁度矩陣的一般公式
{ }
B D B dV{ }
q K{ }
qV [ ] [ ][ ] [ ]
F =
∫∫∫
T ⋅ = (2-40)∫∫∫
= B D B dV
K] [ ] [ ][ ]
[ T
[ ]
K 稱為元素勁度矩陣 二、 整體結構分析1. 建立整體結構的靜力平衡方程式
將局部座標系中各個元素的勁度矩陣
[ ]
K ,經由座標轉換並加以 組集,可得結構的總勁度矩陣K,把各個節點外力對應組集,可得結 構的總外力行矩陣F,最後可以求出整體結構的平衡方程式Kq
F = (2-41) 2. 進行邊界條件處理
將結構節點位移的支承條件引入式(2-41),使之成為正定的線性 方程組。
3. 解方程組,求節點位移
針對結構的總勁度矩陣K的特點,選用合適的線性方程組求解方 法,對邊界條件處理後的方程組F =Kq求解,可求得結構各節點的位 移q。
4.
根據節點位移求應力若要計算應力,則在計算出節點位移q後,從中挑選出欲計算元 素對應的節點位移,組成
{ }
q ,應用{ }
ε =[D]{ }
q 和{ }
σ =[D]{ }
ε =[D][B]{ }
q , 即可求出相應節點的應力。2-2-3 結構動力學的有限元素法分析
對於複雜結構的動力學計算分析,有限元素法也是一種很有效的 工具。動力學問題的有限元素分析方法也是要把彈性連續體結構離散 為有限個元素群的組合體。首先進行每個元素的特性分析,其中包含 了元素勁度矩陣、元素質量矩陣、元素阻尼矩陣的計算,接著再把各 個元素的特性矩陣組集起來,組成整個結構的總勁度矩陣、總質量矩 陣、總阻尼矩陣,從而形成結構的動力學方程式之後再進行求解。以 下就動力學的有限元素法分析的原理加以介紹。
一、 元素特性矩陣的推導
由虛功原理的物理意義可知,外力在虛位移上所做的虛功等於彈 性體內各點的應力在虛應變上所做的虛功的總和。在動力學問題中,
元素外力項除了施加於節點的動態外力,還包括了元素的慣性力及阻 尼力。若元素節點上的激振力為
{ }
F ,則{ }
F 所做的虛功為{ }
q* T{ }
F (2-42a)元素的慣性力與加速度成正比,方向相反。設元素材料的密度為ρ, 則單位體積慣性力可寫成
{ }
f ρ =−ρ{ }
d}
{d 為元素體內任意點加速度函數
將上式乘以元素任意點上的虛位移
{ }
d* ,並對元素體積域積分,可得 慣性力所做的虛功為{ }
d{ }
d dVV
* Tρ
∫∫∫
− (2-42b)
假設阻尼力與速度成正比,方向相反,即假設為線性黏性阻尼,其阻 尼係數為υ,則單位體積阻尼力為
{ }
f υ =−υ{ }
d}
{d 為元素體內任意點速度函數 同理,阻尼力所做的虛功為
{ }
d{ }
d dVV
* Tυ
∫∫∫
− (2-42c) 由式(2-42a)、(2-42b)、(2-42c),可將動力學問題的虛功方程寫為
{ }
q{ } { }
d{ }
d dV{ }
d{ }
d dV{ } { }
dVV V
V * Tρ * Tυ ε* T σ
* T F −
∫∫∫
−∫∫∫
=∫∫∫
(2-43)因為有限元素法中施加的外力必須先轉化成作用於節點上的等值 力,因此須經由關係式{d*}=[N]
{ }
q* 、{d}=[N]{ }
q 、{d}=[N]{ }
q 等,將 式(2-43)改寫成節點位移表示式{ }
q B D B{ }
q dV{ }
q N N{ }
q dVV
V * T[ ]T[ ][ ]
∫∫∫
* T [ ]Tρ[ ]∫∫∫
+{ }
q* T [ N]T [ N]{ }
q dV{ }
q* T{ }
FV =
+
∫∫∫
υ (2-44) 其中虛位移{ }
q* 和{ }
q 、{ }
q 、{ }
q 均可看成常量,可以提到積分號外,又
{ }
q* 為任意值上均成立,可以把它從等式兩端消去,這樣式(2-44)可寫為
{ }
q N N dV{ }
q dVB D
B V
V ⋅ +
∫∫∫
⋅∫∫∫
[ ]T[ ][ ] [ ]Tρ[ ]{ } { }
F ][ ]
[ T ⋅ =
+
∫∫∫
V N υ N dV q (2-45) 式(2-45)的物理意義為動力學問題中元素節點上的彈性力、慣性力、和阻尼力與外加激振力相平衡。
其中若令元素質量矩陣
[ ]
M N N dVV[ ]T [ ]
∫∫∫
= ρ (2-46) 元素阻尼矩陣
[ ]
C =∫∫∫
V [ N]Tυ[ N]dV (2-47)元素勁度矩陣
∫∫∫
= V B D B dV
K] [ ] [ ][ ]
[ T (2-48) 則式(2-45)可改寫成
[ ]
M{ }
q +[ ]
C{ }
q +[ ]
K{ } { }
q = F (2-49) 式(2-49)稱為元素的動力學方程式。二、 整體結構分析 1. 元素矩陣的座標轉換
在建立元素的動力學方程式之後,可經由座標轉換來得到整體彈 性結構方程式。已知元素的座標轉換矩陣
[ ]
T ,可寫出元素節點虛位移{ }
q 、* 節點位移{ }
q 、節點加速度{ }
q 在統一座標系及局部座標系的轉換關係式
{ }
q =[ ]
T{ }
q (2-50){ }
q =[ ]
T{ }
q (2-51){ }
q* =[ ]
T q{ }
* (2-52) 同理,元素節點力的座標變換關係為{ }
F =[ ]
T{ }
F (2-53){ }
F =[ ]
T T{ }
F (2-54)上式中上方加橫線的量是相對於統一座標系的,而不加橫線的量是相 對於局部座標系。
首先介紹元素勁度矩陣的座標轉換。將式(2-40)代入式(2-54)可得
{ }
F =[ ] [ ]
T T K{ }
q再將式(2-50) 代入得
{ }
F =[ ] [ ][ ]
T T K T{ }
q (2-55) 式(2-55)左端{ }
F 為統一座標系下元素的節點力,右端{ }
q 為統一座標 系下元素的節點位移,因此統一座標系下元素的勁度矩陣[ ]
K 可知為[ ]
K =[ ] [ ][ ]
T T K T (2-56)也就是說,只要先求得局部座標系下的元素勁度矩陣
[ ]
K ,以及元素 的座標轉換矩陣[ ]
T ,之後再進行[ ] [ ][ ]
T T K T 的乘積計算,就可求得統一[ ]
接著介紹元素質量矩陣的座標轉換。於兩種座標系中,元素慣性 力所做的虛功應相等,即
{ } [ ]
q* T(
− M{ }
q)
={ }
q* T(
−[ ]
M{ }
q)
(2-57) 將式(2-51)、 (2-52)代入等號右項,可得{ } [ ]
q* T(
− M{ }
q)
={ }
q* T[ ] [ ][ ]
T T(
− M T{ }
q)
(2-58) 因為虛位移{ }
q* T為任意,所以可以由等號兩端消去得[ ]
M{ }
q =−[ ] [ ][ ]
T T M T{ }
q−
比較上式兩端,就可以得出統一座標系下的質量矩陣
[ ]
M =[ ] [ ][ ]
T T M T (2-59) 若將式(2-46)代入得[ ]
M[ ]
T N N[ ]
T dV N N dVV
V T[ ]T [ ]
∫∫∫
[ ]T [ ]∫∫∫
== ρ ρ (2-60)
式中
[ ]
N 稱為統一座標系下的形狀函數矩陣,若已知[ ]
N ,則統一座標系下的質量矩陣
[ ]
M 可由式(2-60)求得。最後由相同原理,可得兩種座標系中元素阻尼矩陣的變換公式
[ ]
C =[ ] [ ][ ]
T T C T (2-61) 2. 整體結構動力學方程式的建立與靜力學問題相同,於動力學中建立元素特性矩陣之後,也要把 元素群綜合成為結構的計算模型,再建立起該計算模型的動力學方程 式,並以求解計算模型的動力學方程式來近似替代求解原來的彈性連
續體結構的動力學方程。先將各元素座標變換後的
[ ]
K 、[ ]
M 、[ ]
C 進行組集,則得總勁度矩陣 K、總質量矩陣 M 和總阻尼矩陣 C,再將各 元素的
{ }
F 進行組集,則得總外力行矩陣 F,因此對於整個結構而言,就有
F Kq q C q
M+ + = (2-62) 式(2-62)即為整個彈性結構的動力學方程式,也就是用有限元素法來 求解動力學問題的基本方程式。與靜力學分析相同,式(2-62)尚未考 慮結構的支承條件,因此也要先進行支承條件處理後再加以求解。