- PDF 免费下载

Size: px

Start display at page:

Download ""

芸弘
7 years ago
Views:

1 科組別 : 物理科別 : 高中組作品名稱 : 非線性振動 - 單擺運動方程式之數值分析與研究關鍵詞 : 非線性單擺數值分析編號 :416

3 作品名稱 : 非線性振動 - 單擺運動方程之數值分析與研究壹. 摘要 : 在我們的這主題中, 以單擺為主要研究對象, 及探討其中在於我們所學之外, 所無法經由更精確的運算, 求出其運動方程式之解析解但如何利用數值方法寫成計算程式, 並藉以整理運用繪圖軟體, 描繪其運動的模式之物理意義, 便是我們此主題的重點所在, 並且我們更深入於實際物理情況中, 消耗力或摩擦力終將阻滯運動以迄振動不在發生就我們而言一開始利用 RK4 寫成 Visual Basic 6. 計算程式來解其運動方程得到等時距之角度及角速度, 再藉由統整後使用 Matlab 5.3 繪出我們想要的圖形 - 角度與時間, 角速度與時間以及角速度與角度之相圖所得的圖形中, 也分為有或無阻尼, 每項中又有不同的討論, 如 : 在已知的任何一個初始狀態下, 其擺動的情況當然也有互相之對照比較, 最終得出精采且史無前例的成果 : 且將單擺作了完整且一般性的分析初角度無阻尼週期性封閉相徑 ( 力學能守恆 ) 初角速度次阻尼無固定週期性振動 ( 非封閉臨界阻尼無週期性非振動最快到達終點相徑, 力超阻尼無週期性非振動較慢到達終點學能耗散 ) 貳. 研究動機 : 在上學期物理課中, 我們學到了單擺與簡諧運動 ; 但是所接觸到的僅是小角度 ( 擺幅 5 度以內 ) 的擺動, 因為小角度擺動時, 才近似於簡諧運動 ; 它的角度與時間的關係是正餘弦函數而且其相圖是橢圓函數但是非小角度的時候呢? 經詢問老師及查閱參考資料後發現了它的運動方程式並不是線性的, 也無法找到 exactly solution, 所以只能夠利用各種數學方法, 去尋找它的近似解, 而所可利用的數學方法中, 有幾種是蠻簡單而易懂的, 包括 Runge-Kutta method 四階 (RK4, 即在泰勒展開中可準確至第四階 ) 及 Verlet method 及 Euler-Cromer method, 我們於是使用上述之方法, 去找出非小角度單擺之數值解希望對單擺之運動有一個完整的認識 : 究竟在各種擺角下它的週期如何變化? 角度與時間的關係如何呢? 相圖又是長的怎樣呢? 參. 研究目的 : 一. 已知無阻尼單擺之運動方程式為 : g θ = - Sin θ = -Ω Sin θ l 二. 有阻尼單擺之運動方程式為 : θ = -a* θ - l g Sin θ = -a*θ - Ω Sin θ 1

便是我們此主題的重點所在, 並且我們更深入於實際物理情況中, 消耗力或摩擦力終將阻滯運動以迄振動不在發生就我們而言一開始利用 RK4 寫成 Visual Basic 6.

4 上兩式為非線性微分方程式, 在給定初值條件下我們要使用微分方程式的數值方法找出 θ (t) 及 θ (t) 在 t = h,h,3h,4h nh 各等時距所對應之近似解 θ 1, θ, θ 3,., θ n 及 θ 1, θ, θ 3 θ n ( 其中 h 為近似法中所選取之時間間隔 ), 並且畫出 θ n -t 及 θ n -t 之關係圖以及 θ n- θ n 之相圖肆. 研究設備及器材 : 工程用電子計算機,Visual Basic 6., Matlab 5.3 伍. 研究過程及方法 : 一. 所使用的方法及其原理 : ( 一 ) 關於 RK4: RK4 是實用上極為廣泛的方法之一, 其精神在於這是一個精確度與計算量妥協後的最佳選擇, 我們在此不擬推導此方法之由來, 僅將其原理列於下表 : R-K-4( f, t, θ, θ, h, N ) 本演算法計算出值問題 θ = f(t, θ, θ ), θ (t) = θ, θ (t) = θ 在等時距 t1 = t + h, t = t + h,, tn = t+nh 之解 ; 此處 f 為使得本問題在區間 t, tn 有唯一解 INPUT : 初值 t, θ, θ, 間距 h, 步驟數 n OUTPUT: 解 θ (tn+1) 在 tn+1 = t + ( n+1 )h 的近似 θ n+1, 其中 n =, 1,,, N-1 對 n =, 1,,, N-1 做 : k1=h f(tn, θ n, θ n )/ k=h f(tn+h/, θ n +k, θ n +k1)/ 其中 k= h(θ n +1/k1)/ k3=h f(tn+h/, θ n+k, θ n +k)/ k4=h f(tn+h, θ n+l, θ n +k3)/ 其中 L=h(θ n +k3) tn+1= tn+h θ n+1= θ n+h( θ n + (k1+k+k3)/3) OUTPUT θ n+1, θ n+1 { 解在 θ n+1 的近似值 } θ n+1= θ n +1/3(k1+k+k3+k4) [ 下一步要用到的輔助值 ] STOP END R-K-4 得到 k1 = - h(a θ n + Ω Sin θ n)/ k = - h[a( θ n + k1)+ Ω Sin( θ n+k)]/ k=h(θ n + k1/)/

所使用的方法及其原理 : ( 一 ) 關於 RK4: RK4 是實用上極為廣泛的方法之一, 其精神在於這是一個精確度與計算量妥協後的最佳選擇, 我們在此不擬推導此方法之由來, 僅將其原理列於下表 : R-K-4( f, t, θ, θ, h, N ) 本演算法

5 k3 = - h[a( θ n + k)+ k4 = - h[a( θ n + k3)+ Ω Sin( θ n+k)]/ Ω Sin( θ n+l)]/ L= h(θ n + k3) θ n+1 = θ n+h[ θ n +( k1+ k+ k3)/3] θ n+1 = θ n +( k1+ k+ k3+ k4)/3 ( 二 ) 使用 Microsoft Visual Basic6. 將上述 RK4 的方法寫成計算程式 : Private Sub Command1_Click() h = Text1.Text a = Text.Text w = Text3.Text x = Text4.Text y = Text5.Text n = Text6.Text For i = 1 To n xo = x yo = y k1 = -.5* h *(a* yo + w * Sin(xo)) k = -.5* h *(a*(yo + k1) + w * Sin(xo +.5 * h * yo +.5 * h * k1)) k3 = -.5* h *(a*(yo + k) + w * Sin(xo +.5 * h * yo +.5 * h * k1)) k4 = -.5* h *(a*(yo + * k3) + w * Sin(xo + h * yo + h * k3)) x = xo + h *(yo + (k1 + k + k3) / 3) y = yo + (k1 + * k + * k3 + k4) / 3 Text7.Text = Text7.Text &" x(" & Str(i) &")= "& Str(x) &" y(" & Str(i) &")= "& Str(y) Text8.Text = Text8.Text &" x(" & Str(i) &")= " & Str(x) Text9.Text = Text9.Text &" y(" & Str(i) &")= " & Str(y) Next End Sub Private Sub Command_Click() Text7.Text = "" Text8.Text = "" Text9.Text = "" End Sub Private Sub Command3_Click() h = Text1.Text a = Text.Text w = Text3.Text x = Text4.Text y = Text5.Text n = Text6.Text 3

6 For i = 1 To n xo = x yo = y k1 = -.5* h *(a* yo + w * Sin(xo)) k = -.5* h *(a*(yo + k1) + w * Sin(xo +.5 * h * yo +.5 * h * k1)) k3 = -.5* h *(a*(yo + k) + w * Sin(xo +.5 * h * yo +.5 * h * k1)) k4 = -.5* h *(a*(yo + * k3) + w * Sin(xo + h * yo + h * k3)) x = xo + h *(yo + (k1 + k + k3) / 3) y = yo + (k1 + * k + * k3 + k4) / 3 Text7.Text = Text7.Text &" "& Str(x) &" "& Str(y) Text8.Text = Text8.Text &" "& Str(x) Text9.Text = Text9.Text &" "& Str(y) Next End Sub Private Sub Command4_Click() If Command4.Caption = "x,y,x,y" Then Command4.Caption = "x,x,y,y" Text7.Visible = False Text8.Visible = True Text9.Visible = True Exit Sub Else Command4.Caption = "x,y,x,y" Text7.Visible = True Text8.Visible = False Text9.Visible = False End If End Sub Private Sub Form_Load() Dim h, a, w, x, y, xo, yo, n, test, i, k1, k, k3, k4 WindowState = End Sub ( 三 ) 如何選擇適當的時間間距 h? 以及步驟數 N? 因為 RK4 的方法中每一步驟的截斷誤差為 h 5 階 ( 此證明可參考 Collatz,L.,The Numerical Treatment of Differential Equations. 3rd ed.),h 不可過大否則每一步驟之截斷誤差將會增大 ( 例如以 h 來代替 h, 每一步驟之截斷誤差會增大約 3 倍 ), 因為我們將在每個已知的 Ω 下考慮不同的初始條件時皆固定使用相同的步驟數, 似乎 h 之選擇當然是越小越好, 然而步驟數 N 與 h 之乘積也不可過小否則我們將看不到一個週期 T 的循環即 h*n 最好是超過一個週期, 所以 h 太小時步驟數 N 必須相當大, 這會我們的作業時間延誤許多 ( 電腦處理速度 ) 所以 h 值需適當小即可 4

Text &" "& Str(x) Text9.Text = Text9.Text &" "& Str(y) Next End Sub Private Sub Command4_Click() If Command4.Caption = "x,y,x,y" Then Command4.Caption = "x,x,y,y" Text7.Visible = False Text8.

7 ( 四 ) 如何可知所做出之結果其誤差大或小呢? 因為無阻力之單擺運動為力學能守恆的, 也就是其相圖在力學能 E < mgl( 取最低位置為零位能 ) 時必須是封閉的, 所以我們可以在一開始時嘗試各種不同的 h 值, 去分別作出其相圖並觀察其封閉性與唯一性, 從嘗試中去體驗 h 的控制二. 過程及步驟 ; ( 一 ) 先考慮無阻力時 : 設 Ω =1, a = θ ()= θ = π/9, π/3, π/18, π/9, π/6, π/4, π/3, 5π/1, π/, π/3, 5π/6,179π/18 θ ()= θ 步驟數 N=1 1. 先分別用 h=1.,.8,.6,.4,.,.1 來計算, 用 MatLeb5.3 將所得到的 θ n 及 θ n 去作 θn -tn, θ n - tn, θn - θ n 圖並特別留意相圖的部分. 重複上述步驟將 N 改為及 4 並做同樣之觀察 3. 由前兩步驟中體驗出 h 與 N 之取法之最佳搭配並做為爾後步驟之參考 ( 二 ) 改 Ω =.1,.6,.1,.6, 1., 6., 1, 6, 1, 6, 1 重複 ( 一 ) 中之步驟, 但 h 必須根據 ( 一 ) 之經驗及原理中之敘述作適當的更換同樣找出相圖, θ n -tn, θ n - tn 圖 ( 三 ) 由 ( 一 ) 及 ( 二 ) 所得之較佳結果中, 觀察 θ = 上述初角且 θ = 之 θn -tn, θ n - tn 圖中 : 1. 由數據中利用內插法尋找 θ n = 及 θ m = 所對應之時間 ntn+1 及 mtm+1 ; 即 θ (ntn+1) =, ntn+1 是位於 nh 及 (n+1)h 之間, 又 θ (mtm+1)=, mtm+1 是位於 mh 及 (m+1)h 之間 ; 利用牛頓向前差分內插公式 Newton s forward difference interpolation formula 如下 : n r f(x) pn(x) = ( s ) s f (x = x+h, r = (x-x)/h) s= r( r 1) r ( r 1)( r ) ( r n + 1) = f+r f+ f+ +! n! n f 去逼近找出所找出之 ntn+1 即相當於 1/4 週期, 3/4 週期,5/4 週期, ; 所找出之 mtm+1 即相當於 1/ 週期,1 週期,. 將 ( 一 ) 及 ( 二 ) 中 θ 改為 π, θ = +.15, +.5, +.4, 同樣找出相圖, θ n -tn, θ n - tn 圖 ( 四 ) 由力學能守恆知 :(θ ) = Ω (Cos θ -Cos θ )+( θ ) 將式中 θ 取 ( 一 ) 中之諸初角度 θ 取 π 時, θ = 所對應之 θ 之值叫作 ( θ )min, 再取比 ( θ )min 為大之值當作新的 θ 並配合 θ =( 一 ) 中之諸初角度之初始條件同樣找出相圖, θ n -tn, 5 θ n - tn 圖 ( 五 ) 再考慮有阻力時 : 同樣使用無阻力時之 Ω 值以及 θ, 因為阻尼擺動分為次阻尼, 臨界阻尼以及超阻尼所以必須嘗試配合各種 θ : 1. 要區分這三種運動, 首先要找出臨界阻尼所對應的 a = a, 也就是說 (1) a< a 時為次阻尼 () a=a 時為臨界阻尼

制二. 過程及步驟 ; ( 一 ) 先考慮無阻力時 : 設 Ω =1, a = θ ()= θ = π/9, π/3, π/18, π/9, π/6, π/4, π/3, 5π/1, π/, π/3, 5π/6,179π/18 θ ()= θ 步驟數 N=1 1. 先分別用 h=1.,.8,.6,.4,.,.1 來計算, 用 MatLeb5.

8 (3) a>a 時為超阻尼. 如何尋找在已知 Ω 及 θ 時之 a : 因為所謂臨界阻尼即是在已知 θ 之下其運動最快回到 θ= 的位置, 所以在使用 Microsoft Visual Basic 6. 來跑程式時, 就必須多次觀察 a 取何值時 θ n 最快接近 3. 當找到 a 之後, 就分為次阻尼與超阻尼來分頭進行 4. a< a, 使用各種 Ω ; 各種 θ 及 θ 找出相圖, θ n -tn, θ n - tn 圖, 並利用內插法求出 θn -tn 圖中曲線與 t 軸之交點並與相同初始條件之無阻尼時作比較 5. a>a, 使用各種 Ω ; 各種 θ 及 θ 找出相圖, θ n -tn, θ n - tn 圖陸. 研究結果與討論 : ( 一到七為無阻尼, 八到九為有阻尼 ) 一. 在某初角 θ ()= θ o < π, 初角速度 θ ()= θ 'o = 在 (, θ o ) 之切線斜率, 所得之 θ -t 圖 : 此圖中諸曲線分佈可被代表著不穩定振動 (E= mgl) 之兩條分界虛線 : θ 'o = ± Ω Cos( θ o ) 分為三個區域 : ' ( 一 ) θ o < ± Ω Cos( θ o ) 此區內表示為穩定振動 (E<mgl), 且各有不同週期, 圖中紅和紫之 θ ' o 大小相等, 方向相反 ( 即具有相同能量 E), 有相同週期 T, 而因其 θ ' 甚小故 T 與 o θ ' o = 之 T 很接近 ; 綠和黃亦有等大之 θ ', T 也同 ( 見圖中粗箭頭 o 所指 ), 已可發現 : E 愈大則 T 愈大之定性關係. 後面再詳述 T 與 E 及 Ω 之關係. θ o ) 及 θ ' θ o ) 此區內曲線以 θ = θ o 之虛 ' ( 二 ) θ o > Ω Cos( o < Ω Cos( 點線為對稱, 表示為非振動的, 仍是週期性繞支軸作全周旋轉的擺 (E>mgl). 在 θ = θ o + π 所對應之 t 即為 T ' ( 三 ) θ o = Ω Cos( θ o ) 表示擺可恰好衝達最高點瞬間停下, 但再來如何? 這是不穩定的平衡 ; 有可能的情形是 (1) 受到任何極微小地擾動, 造成 E=mgl ± δ ; 擺必以確切為零之角速到達 θ = nπ 各點之一, 有確切 T, () 否則只在無限長時間後 (T ) 6

a>a, 使用各種 Ω ; 各種 θ 及 θ 找出相圖, θ n -tn, θ n - tn 圖陸. 研究結果與討論 : ( 一到七為無阻尼, 八到九為有阻尼 ) 一.

9 二. 此圖為 θ =π ; θ ' > 及 θ ' <; ( 亦即 E >mgl) 條件下, θ -t 圖, 當然以 θ =π 之水平線為對稱的, 圖中 θ ' 為綠 > 紅 > 藍, E 綠 > E 紅 > E 藍, 這些全都是繞支軸作 7

10 全周旋轉的擺, 可由虛線對應之 t 值看出 E 愈大者旋轉週期愈小. 這些情形之旋轉周期可以由 θ ' = g Sin ( θ θ l )-Sin ( )+ l g θ ' 1 dt= 1 l g 1-Sin ( θ )+ g l θ ' 1 d θ,( θ =π ) T= l g π π 1-Sin ( θ )+ g l θ ' 1 d θ T= 4k Ω 1 dx ( 1 k x )(1 x ) T= 4k Ω ' = θ Ω g 其中 Ω =, 1 k [( ) + 1], x = Sin( θ ) l 因為 k<1, 故 (1-k x ) -1/ 可展開為級數後 : dx k x 3k x 5k x [ ]. 5 (1 x ) 8 16 πk = Ω k 9k 5k 35 k ( +...) 如果 k 近於 1 則需許多項方能產生準確結果, 不過對於甚小的 k( 即 θ ' 甚大 ;E>>mgl) 上式迅速收斂 8

)+ g l θ ' 1 d θ T= 4k Ω 1 dx ( 1 k x )(1 x ) T= 4k Ω ' = θ Ω g 其中 Ω =, 1 k [( ) + 1], x = Sin( θ ) l 因為 k<1, 故 (1-k x ) -1/ 可展開

11 三. 這兩張圖是週期 T 與 E*( 小於 1) 之關係圖 : θ ' 同樣地, 我們將各種給定之 Ω, 使用 θ = 及 = π π, 9 45 π π π π π π π 5π 35π,,,,,,,, 去計算出 θ 與 t 之關係圖, 在 n 利用內插法得出 T 與 θ 之關係, 因為 E*=Sin ( θ ), 故也就可以得到 T 與 E* 之關係了, 上圖為全圖, 下圖為局部份放大圖, 對各曲線而言, T 皆隨 E* 之增加而增大, 並且當 E* 愈接近 1, T 增加特別快, 可以用左邊板第 πk k 9k 5k 35 k 二圖之分析中 : T= ( ) for E*>1, 令 k=1 Ω π 則 T ( ) 或是使用前述之 E*<1 之週期 Ω π k 9k 5k 35 T= ( ), 令 k=1 Ω π 則 T ( ) Ω π 故當 E* 1 時, 週期 T 會近似於 ( ) Ω 不論從 E*<1 或是 E*>1 去趨近都是相同的, 9

4 6 8 k 9k 5k 35 k 二圖之分析中 : T= ( 1+ + + + +...) for E*>1, 令 k=1 Ω 4 64 56 18 π 1 9 5 35 則 T ( 1+ + + + +.

12 1

13 四. 此部份為 ( E*<1) 週期之討論 : 這個圖為各種給定之 E* 時 ( E*<1), 經由所算出 θ n 與 t 之關係再利用內插法得出週期 T 與 Ω 關係於上圖是全圖, 下圖是局部放大圖, 現在說明之 : ( 一 ) 圖中 E*1< E*<..< E*1 < E*11 分別代表著 θ ' = 且 θ = π π, 9 45 π π π π π π π 5π 35π,,,,,,,, ; 可以理解的隨著 E* 之增加而週期迅速增大, 且因 E*1 及 E* 近似簡諧運動, 故週期幾乎重合, 可說明再很小擺角下單擺近似等時性的對於大角度之擺動討論如下 ( 二 ) 對於任意角度之擺動週期討論如下 : E 因為 E* = = [ 1 ' θ θ θ ( ) + Sin ( ) ] = ( ' Sin ) mgl 4 Ω θ = θ 1 ' g θ l θ θ = [ Sin ( ) Sin ( )] T = g l θ 1 θ [ Sin ( ) Sin ( )] dθ T = [(1 x )(1 k x )] dx Ω, g Ω =, l sin( θ ) x =, k = Sin( ) sin( ) θ θ 欲產生振動須 θ <π 即 k<1 如此上式可展開 (1 x ) k x 3k x k = dx k x 3k x π k 9k 5k 35 T = [ ]. 5 Ω = ( ) (1 x ) 8 Ω π = Γ( n +.5) n [ ( E*) ] Ω = π Γ( n + 1) n, E*=k 當然, 若 E* 為某定值, 則 T 與 Ω 成反比, 所以此兩圖為各種 E* 之下之反比等週期曲線, 11

等時性的對於大角度之擺動討論如下 ( 二 ) 對於任意角度之擺動週期討論如下 : E 因為 E* = = [ 1 ' θ θ θ ( ) + Sin ( ) ] = ( ' Sin ) mgl 4 Ω θ = θ 1 ' g θ l θ θ = [ Sin ( ) Sin ( )] T = g l θ 1 θ [ Sin ( ) Sin

14 五. 此圖為將所得 θ n 及 θ 繪成之相圖, 因力學能 ( 最低點為零位能 ) ' θ ' n 1 θ E = mgl [ + Sin ( ) ] 4 Ω E E* = = [ 1 ' θ θ ( ) + Sin ( ) ] 圖中 E=mgl 即 mgl 4 Ω ' θ θ E*=1 即 = Cos( ) 曲線為分界相徑, 在物理中一個分 Ω 界相徑永遠經過一個不穩定的平衡點, 此分界相徑分隔局部束縛運動 (E*<1, 穩定擺動 ) 與局部非束縛運動 (E*>1, 繞支軸旋轉 ) : 前者 (E*<1) 之情況等於一個質點束縛於位井 U=mgl(1-Cosθ ) 之中, 所以這區域的相徑為封閉的曲線, 由於位能對 θ 是對稱且週期性的, 在. 3 π < θ < π, π < θ < π, π θ < 3π 穩定平衡位置如預料中, 對於甚小之 E* 而言, 相徑幾近乎圓圈, 再此時運動大略為簡諧性後者 (E*>1) 之情況, 運動不是振動而是繞支軸作週期性旋轉的擺 <... 諸區域內存在著相同的相徑, 在 θ =., -π,, π,.. 各點為如果 E*=1, 其相徑不實際代表擺的可連續運動, 如擺 θ =π ( 此乃在 E*=1 相徑上一點 ), 1

U=mgl(1-Cosθ ) 之中, 所以這區域的相徑為封閉的曲線, 由於位能對 θ 是對稱且週期性的, 在.

15 則任何小的擾動便會使運動密切但不確切地沿循由 θ = π 分歧的相徑之一, 此因 E=mgl+δ 如果運動是沿 E*=1 相徑之一, 擺必以確切為零之角速到達 θ = nπ 各點之一, 但只在無限長時間後! θ ' 六. 此圖為將 θn, n, 和 E * 三者繪成 3D 之圖, 可以更明顯且容易地看出三者之關 Ω 係 : E 因為 E* = = [ 1 ' θ θ ( ) + Sin ( ) ], 每一條等能量曲線為同一相徑, 非封 mgl 4 Ω 閉之相徑能量較高代表繞支軸之旋轉運動, 而封閉之相徑能量低 ; 代表來回之擺動, 而 E*=1 之相徑代表不穩定之運動 ; 將空間隔成上下兩區, 其觀念如同平面相圖中所述 3D 相空間幾何式地代表簡單振動系統的動力學是極有用的, 在往後研究中預計朝向多質點系統之相空間探討並驗證 Liouville s Theorem 13

此圖為將 θn, n, 和 E * 三者繪成 3D 之圖, 可以更明顯且容易地看出三者之關 Ω 係 : E 因為 E* = = [ 1 ' θ θ ( ) + Sin ( ) ], 每一條等能量曲線為同一相徑, 非封 mgl

16 14

17 七. 此圖為 T 與 E*, Ω 之關係 3D 圖, 由前述中 : π T= Γ( n +.5) n [ ( E*) ] Ω = π Γ( n + 1) n 可以驗證此 3D 圖之長相由 Ω 處看去即可得 T- Ω 圖, 由 E* 處看去即可得 T-E* 圖至此, 研究單擺之力學架構之結語與心得 : ( 一 )E*(θ, θ ' ) : 能量守衡觀點與運動分析 ( 二 )T (E*, Ω ) : 週期性運動之週期探索方法 ( 三 ) 使用了很好的微分方程數值分析對於非線性力學有了初步之認識 ( 四 ) 數學軟體在數據分析時提供了適當之協助八. 次阻尼擺動 : - 雖然不像簡諧阻尼振動 ( 參考質點和系統的古典動力學一書 :(p94~p11) 有著解析的分類阻尼參數去分隔次阻尼, 臨界阻尼, 超阻尼三種情形, 我們還是可以 15

心得 : ( 一 )E*(θ, θ ' ) : 能量守衡觀點與運動分析 ( 二 )T (E*, Ω ) : 週期性運動之週期探索方法 ( 三 ) 使用了很好的微分方程數值分析對於非線性

18 利用數值分析的方法配合前面所述寫下的程式, 再將對一組已知的初始條件下臨界阻尼接近平衡的時間較次阻尼或超阻尼都要來得短的特性, 也可以找出在已知的 Ω 下所對應的臨界阻尼參數 a, 所找出的結果如下表 : Ω a 在質點和系統的古典動力學一書中對簡諧阻尼振動有詳細求出臨界的阻尼參數 a = Ω 今以 Ω =1 對應臨界阻尼參數 a = 為例來討論各種阻尼擺動 : 這六張圖是阻尼參數 a=.1 和.5 之次阻尼擺動的 θ- t,θ -t 以及相圖, 我們發現 ( 一 ) 因為 θ = -a*θ - l g Sinθ = -a*θ - Ω Sinθ 故所有可能的擺動過程中平衡點 ( 合力矩 =) 必出現在 θ =, 即曲線 Ω θ = - a Sinθ = -1* Sinθ 代表的點集合上面. 當其相徑通過該平衡曲線時, θ ' 恰為極大值或極小值 ( 在圖中以黑逗點標示 ), 可以理解平衡點並不一定是在最高或最低點上 : (1) 若通過最高點之速度不為, 則平衡點出現在重力的切線分量與阻力抵銷之處 ( 如右圖 ) () 若通過最高點之速度恰為, 則平衡點出現在最高點但卻為不穩定的平衡點, 如相圖中所示 : 深藍色的平衡曲線通過在 -3π,-π,π,3π 這 16

19 些位置上, 水藍色箭頭可表示出其為不穩定的平衡點 ( 二 ) 運動不是週期性的, 或者說單擺永不會以相同的速度通過一點兩次, 與前所討論的無阻力單擺對比, 阻尼單擺的能量不是持久不變的, 它的能量是持續地付給阻尼性介質並消耗為熱 ( 或者可能以流體波的形式化為輻射 ). 能量消耗率與 θ 平方成正比, 當單擺於平衡點上而到達最大速度時, 耗損率將為最大 ( 三 ) 相圖中所標示的 A(pi),A(4pi),A(6pi) 代表各種初狀態之相徑終點分別為 pi,4pi,6pi, 也就是說給定任何一組 (θ, θ ' ) 我們就可以根據相圖中相徑終點的分隔區域來決定它的相徑路線和終點這和之前所討論的無阻力情形是截然不同的 17

能量消耗率與 θ 平方成正比, 當單擺於平衡點上而到達最大速度時, 耗損率將為最大 ( 三 ) 相圖中所標示的 A(pi),A(4pi),A(6pi) 代表各種初狀態之相徑

20 18

21 19

22 九. 超阻尼擺動 : 依 θn -tn 圖可分為三類現在是有關於超阻尼擺動的討論我們以 Ω =1, 阻尼參數 a=.5 來作例子 : ( 一 ) 首先如同在次阻尼所討論的一樣, 仍然是有著一條以平衡點為集合的曲線 : Ω θ = - a Sinθ =-.4* Sinθ 也就是所有的相徑和該曲線之交點皆為平衡點 ( 二 ) 在這種超阻尼的擺動所有的相徑最後都會極趨近於一條直線 ( 相圖中黑色虛線 ), 且終點必在 -π,,π 這些位置上那這些黑色虛線是如何得的到呢? 因為當任何相徑極靠近終點的時候角度和角速度都非常的小, 那就可以將運動視為簡諧阻尼振動, 所以運動方程式就可以近似為 at θ =-a*θ - Ω t t θ, 此時方程式之解為 θ(t)= ( ω + ω ) e A 1 e A e

23 其中 = a ω Ω, 再代入初始條件後就可得到 θ = 1 + A 4 a = A ( ) A ( ω + a A ' 及 θ 1 ω ) 可解出 A1 及 A, 另令當 t = T 時 θ(t) = at + = A A1 A θ e 1 A1 A1 ' a θ + + ω θ 即生出所需的線方程式因此 ( 三 ) 超阻尼情形產生一種非振動性而漸近於零的幅度減小但隨 θ ' 不同, 在 θ 接近於零之前,θ 可能改變符號對於任一個初始條件 (θ, θ ' ) 所產生的各種相徑可分為三類 : Ⅰ. 擺動過程中會出現兩個平衡點,θ 一成不變趨近於終點 : Ω 即相徑會與平衡點集合曲線 θ = - a Sinθ =-.4* Sinθ 有兩個交點 ( 放大相圖中黑色逗點處 ), 其中一個代表速率極大處, 另一個代表速率極小處, 此種情況是初狀態 (θ, θ ' ) 位於相平面分界圖中 A(1),A(4),B(1), B(4),C(1), C(4) 區塊內換句話說當初狀態 (θ, θ ' ) 1

24 例如 (-4.8,5),(4.8,-5),(1.5,5),(-1.5,-5) 是在這些區塊內則其擺動過程中必會有一個極大速率與一個極小速率的位置現在用 (-4.8,5) 作例說明之 : 如右圖從出發點到第一個平衡點的過程中, 阻力大於重力分量,θ 持續變小, 且阻力減少得比重力分量減少來得快, 當兩力相等時合力矩為零即為第一個平衡點, 此時 θ 為極小 ( 約.85) 由第一個平衡點到第二個平衡點的過程中, 重力分量大於阻力,θ 逐漸變大, 且重力分量增加比阻力增加來得快, 當兩力再度相等時合力矩為零, 即到達第二個平衡點, 此時 θ 為極大 ( 約.39) 由第二個平衡點到終點的過程中, 阻力大於重力分量,θ 迅速地減少, 重力分量和阻力均同時迅速減少至零終歸於原點 Ⅱ. 擺動過程中只出現一個平衡點,θ(t) 在趨近於零之前某 t > 時到達一個 Ω 極值 : 即相徑會與平衡點集合曲線 θ = - a Sinθ =-.4* Sinθ 只有一個交點 ( 放大相圖中黑色逗點處 ), 此平衡點處必出現速率極大, 此種情況是初狀態 (θ, θ ' ) 位於相平面分界圖中 A(3),A(6),B(3),B(6), C(3), C(6) 區塊內, 換句話說當初狀態 (θ, θ ' ) 例如 (-1.5,6),(1,4), (1.5,-6),(-1,-4) 是在這些區塊內則其擺動過程中必只有一個速率極大之處現在用 (1,4) 來為例說明 ( 如右圖 ): 由出發點到最大角度之過程中, 重力分量與阻力同方向使得速率持續減少至零由最大角度到平衡點的過程中, 重力分量大於阻力使得速率漸

25 增且阻力漸增, 當阻力等於重力分量時即到達平衡點, 此時速率出現極大值由平衡點到終點的過程中, 阻力恆大於重力分量且反方向故速率持續減少至零 Ⅲ. 擺動過程中不會出現平衡點,θ 一成不變趨近於終點 : 即相徑在到達終點前決不會與平衡點集合的曲線 Ω θ = - a Sinθ =-.4* Sinθ 有任何相交, 此種情況是初狀態 (θ, θ ' ) 位於相平面分界圖中 A(), A(5),B(),B(5),C(), C(5) 區塊內, 換句話說當初狀態 (θ, θ ' ) 例如 (-4, 6),(-3., 6),(3., -6),(4, -6) 是在這些區塊內則其擺動過程中 θ 和 θ 都一成不變得減少至零現以 (4, -6) 為例說明 : 由出發點到最高點的過程中, 重力分量與阻力同方向且皆與速度反向, 故速率漸減由最高點到終點的過程中, 重力分量與阻力反方向且阻力恆大於重力, 故速率仍然持續減少至零 3

26 4

27 5

28 柒. 結論 : 一. 其實所有無阻尼的實驗結果可以分為三部份 : ( 一 ) 力學能 E<mgl( 以最低點為零位能 ): 在此條件之下, 不論初角度與初角速度如何其相圖皆為封閉的, 而且 E 越小時相圖越接近橢圓也就是其運動越接近簡諧運動, 當 E 越大時期相圖會以類似菱形放大但會到達一個極限 ( 二 ) 力學能 E=mgl( 以最低點為零位能 ): 在此條件之下, 不論初角度與初角速度如何其相圖為餘弦函數此乃代表不穩定之情況 ( 在最高點為不穩定的平衡點 ) ( 三 ) 力學能 E>mgl( 以最低點為零位能 ): 在此條件之下, 不論初角度與初角速度如何其相圖皆為非封閉的, 這代表著它不是一個振動而是同方向的迴旋轉動二. 有阻尼時其相圖皆為非封閉的, 代表著力學能為遞減, 且角速度出現極值的地方並不是在最高或最低點, 由數據發現次阻尼時不可能界定出一個固定的週期, 而且也不會以相同的速度通過同一點兩次, 而臨界阻尼與超阻尼並不是振動的同時, 在超阻尼的情況相徑最終會趨近於一條直線, 該直線方程式為 Ω ' a a 4 θ = ( ) *θ 4 該方程式可由 θ 和 θ 都趨近於時振動近似為線性簡諧阻尼振動, 找出其方程式的解, 便可求出三. 關於實際情況中之單擺擺動, 我們歸納出 : 6

29 捌. 參考資料 : 一.(1) 編譯者 : 冉長壽 () 書名 : 質點和系統的古典動力學 (3) 版次 : 增定八版 (4) 出版地 : 新店市 (5) 出版社 : 徐氏基金會 (6) 頁數 :55 頁 (7) 出版年 : 民國 77 年月 9 日二.(1) 作者 :ERWIN KREYSZIG () 書名 : 高等工程數學 (3) 版次 : 第 6 版 (4) 出版地 : 臺北 (5) 出版社 : 松崗圖書公司 (6) 頁數 :3 頁 (7) 出版年 : 民國 79 年 7

30 ( 第二名 ) 本作品是用 RK4 方法模擬非小角度單擺的運動, 除了說明有阻尼與無阻尼的兩種狀況下的擺角隨時間之變化外, 也能解釋反折點附近的能量及其消耗變化提供單擺運動相當完整的認識, 作品嚴謹完整, 有學術性 8

6-1-1極限的概念

6-1-1極限的概念選修數學 (I-4 多項式函數的極限與導數 - 導數與切線斜率定義. f ( 在的導數 : f ( h 對實函數 f ( 若極限存在 h h 則稱 f ( 在點可微分而此極限值稱為 f ( 在的導數以 f ( 表示 f ( f ( 函數 f ( 在的導數也可以表成 f ( 註 : 為了