Introduction to Scientific Computing 科學計算導論
Homework 2 (Due: June 4, 2009)-Rev 2 1. 考慮一耦合的簡諧震盪系統 (http://en.wikipedia.org/wiki/Normal mode),
Figure 1: Coupled Harmonic Oscillator from wiki
假設x1(t),x2(t)分別為m1及m2的位移, 三個彈簧的彈性係數皆為k且m1 = m2 = m. 則我們可以用下列方程組 描述整個系統的行為:
d2
dt2x1 =−k
mx1+ k
m(x2− x1), (1)
d2
dt2x2 =−k
mx2+ k
m(x1− x2), x1(0) = x10, x2(0) = x20, x01(0) = v10, x02(0) = v20,
請寫出一個 Matlab 的程式以內建的 ode45 得到 (1) 的數值解 (numerical solution):
(a) 這個系統有兩個 normal modes, 找出這兩個 normal modes 的共振頻率.
(b) 選取適當的初始條件, 分別產生這兩個 normal modes. 檢查數值計算的結果是否與理論推導結果相近.
(c) 討論k/m = 1, k/m = 4兩系統間的關係. 分別選取恰當的初始條件, 將數值解畫在適當的框架中比較. 或 是經過適當的座標變換將兩系統解畫在一起比較.
(d) Question: 對任意初始條件, 本系統的解是否為週期解? Why?
2. 考慮一生物數學模型
d
dtx = x− 0.1x2− xy 1 + x, d
dty = 1
6y(1.2− y x), x(0) = x0,
y(0) = y0.
請寫出一個 Matlab 的程式以內建的 ode45 得到 數值解 (numerical solution):
(a) 本系統有三個平衡點, 其中一個落在左半平面, 兩個落在右半平面. 找出這三個平衡點及其線性穩定性並與數 值計算結果比較.
(b) 這個系統存在一個 Limit cycle(http://en.wikipedia.org/wiki/Limit cycle). 選取數個不同的初始條件, 將數值解的y對x作圖 (Phase portrait), 可以看到這些數值解收斂到 limit cycle 的情形
(c) 選取適當的初始條件, 使系統的解會收斂至 limit cycle, 將數值解 x(t), y(t)對時間t作圖. 是否有看到週期 性?
(d) (bonus)試以 Poincar´e-Bendixson theorem 說明 limit cycle的存在性.
(http://en.wikipedia.org/wiki/Poincar´e-Bendixson theorem)