2/26/2016

[數值方法] Forward Euler's Method 與 其改良型 求解 ODE

給定 初始值問題(Initial Value Problem, IVP)
\[
y'(t) = f(t,y(t)),\;\;\;\; y(0) = y_0
\] 一般而言,上述 一般的 IVP 問題並沒有辦法寫下解析解,但我們可以退而求其次詢問是否有合適的數值方法求解上述 初始值問題,最為基本的想法是利用所謂的 Forward Euler method:

Forward Euler Method:
\[
y_{n+1} = y_n + h \cdot f(t_n, y_n)
\]其中 $t_n = n h$。令終止時間 $T$ 在 $n$ 步之後到達,則 $T = n h$。

Comment:
1. 上式中 $h$ 稱為 迭代步長 (step size)
2. $t_{n+1} = t_n + h$
3. 儘管 Forward Euler Method 非常簡便,但其數值誤差相當大,以下我們看個實際例子:


Example: 
考慮 IVP $ y'=t-y $ 且 $ y(0) = 0 $ ,
(a) 試求解上述 IVP
(b) 現在令 $h := 0.1$ ,試求用 Forward Euler Method 求解 $y_2$
(Note: $y_2 = y(t_2) = y(2h)$)
(c) 比較 (a) 與 (b)

Solution (a):
首先改寫原式:
\[y'\left( t \right) =  - y\left( t \right) + t
\]由上式可知此為線性常係數 ODE 其解可立即求得
\[\begin{array}{l}
y\left( t \right) = {e^{ - 1}}y\left( 0 \right) + \int_0^t {{e^{ - 1\left( {t - s} \right)}}sds} \\
 \Rightarrow y\left( t \right) = 0 + \underbrace {{e^{ - t}}\int_0^t {{e^s}sds} }_{ = {e^{ - t}}\left( {{e^t}(t - 1) + 1} \right)}\\
 \Rightarrow y\left( t \right) = (t - 1) + {e^{ - t}}
\end{array}\]
Solution (b):
由 $y_0 := y(0)=0$ 出發,則我們可以計算
\[\begin{array}{*{20}{l}}
{{y_1} = {y_0} + h \cdot f({t_0},{y_0})}\\
{ \Rightarrow {y_1} = 0 + \left( {0.1} \right) \cdot \left( {0 \cdot h - {y_0}} \right)}\\
{ \Rightarrow {y_1} = \left( {0.1} \right) \cdot \left( {0 - 0} \right) = 0}
\end{array}\]同理,我們接著計算 $y_2$
\[\begin{array}{l}
{y_2} = {y_1} + h \cdot f({t_1},{y_1})\\
 \Rightarrow {y_2} = {y_1} + h \cdot \left( {1 \cdot h - {y_1}} \right)\\
 \Rightarrow {y_2} = 0 + \left( {0.1} \right) \cdot \left( {0.1 - 0} \right)\\
 \Rightarrow {y_2} = 0.01
\end{array}
\]
Solution (c)
由 (a) 可知
\[\begin{array}{l}
y\left( {2h} \right) = (2h - 1) + {e^{ - 2h}}\\
 \Rightarrow y\left( {0.2} \right) =(0.2 - 1) + {e^{ - 0.2}} \approx {\rm{0.01873}}
\end{array}
\]但由 (b) 可知我們得到的 $y_2 = y(2h) =0.01$ 亦即具有誤差
\[
|y(2h) - y_{2}| = |0.01873 - 0.01| = 0.0873
\]

以下為簡單的 MATLAB code 實現 Forward Euler Method:

t0 = 0; %Initial time
y0 = 0; %Initial value
h = 0.1; %Step size
tfinal = 1; %Final time

F = @(t,y) t - y;

y = y0;
yout = y;
for t = 0: h: tfinal-h
    y = y + h * F(t,y); %Euler Method
    yout =[yout; y];

end

t_interval = 0: h: tfinal
plot(t_interval,  yout, 'o-')




Improved Euler Method
為了改良上述誤差過大的問題,我們可採用以下 Predictor-Corrector 架構來改善 Forward Euler Method,亦即我們首先透過 Forward Euler Method 建構

Predictor Equation
\[
\tilde{y}_{n+1} := y_n + h \cdot f(t_n,y_n)
\]接著再修正上述的預測結果:
Corrector  Equation:
\[
y_{n+1} = y_n + \frac{h}{2}(f(t_n,y_n) + f(t_{n+1}, \tilde{y}_{n+1}))
\]
我們以下再次使用前述例子來看看改良型 Euler 是否有得到更加精準的結果:

Example (Revisit): 
考慮 IVP $y'=t-y$ 且 $y(0)=0$ ,
(a) 令 $h := 0.1$ ,試求用 Improved Euler Method 求解 $y_2$

Solution:
由 $y_0 = 0$ 開始,為了計算 $y_1$ 我們先透過 predictor 計算 $\tilde{y}_{1} $:
\[
\tilde{y}_{n+1} := y_n + h \cdot f(t_n,y_n) = 0
\]接著用Corrector  Equation:
\[\begin{array}{l}
{y_1} = {y_0} + \frac{h}{2}\left( {f({t_0},{y_0}) + f({t_1},{{\tilde y}_1})} \right)\\
 \Rightarrow {y_1} = 0 + \frac{h}{2}\left( {\left( {{t_0} - {y_0}} \right) + \left( {{t_1} - {{\tilde y}_1}} \right)} \right)\\
 \Rightarrow {y_1} = \frac{h}{2}\left( {\left( {0 \cdot h - 0} \right) + \left( {1 \cdot h - 0} \right)} \right)\\
 \Rightarrow {y_1} = \frac{h}{2}h = \frac{1}{2}{\left( {0.1} \right)^2} = 0.005
\end{array}
\]有了 $y_1$ 我們再透過 Predictor 計算 $\tilde{y}_{2} $:
\[\begin{array}{l}
{{\tilde y}_2} = {y_1} + h \cdot f({t_1},{y_1})\\
 \Rightarrow {{\tilde y}_2} = 0.005 + h \cdot \left( {{t_1} - {y_1}} \right)\\
 \Rightarrow {{\tilde y}_2} = 0.005 + h \cdot \left( {1 \cdot h - 0.005} \right)\\
 \Rightarrow {{\tilde y}_2} \approx {\rm{0.01450}}
\end{array}\]再用 Corrector 計算 $y_2$
\[\begin{array}{l}
{y_2} = {y_1} + \frac{h}{2}(f({t_1},{y_1}) + f({t_2},{{\tilde y}_2}))\\
 \Rightarrow {y_2} = 0.005 + \frac{h}{2}(\left( {{t_1} - {y_1}} \right) + \left( {{t_2} - {{\tilde y}_2}} \right))\\
 \Rightarrow {y_2} = 0.005 + \frac{{0.1}}{2}(\left( {1 \cdot h - 0.005} \right) + \left( {2 \cdot h - {\rm{0.01450}}} \right))\\
 \Rightarrow {{y}_2} \approx {\rm{0.0190}}
\end{array}
\]現在回憶前述例子中,我們知道 $0.01873$ 此表明此例中,improved Euler 確實大幅降低計算誤差。

Comments:
1. 實際數值應用中,大多數數值方法採用所謂 Runge-Kutta Method 來獲取更加精確或者穩定的數值解,但未免離題,在此不做贅述。

2/18/2016

[基礎數學] 區別引理,定理,系理以及命題

在一般數學教科書或者文獻之中,我們很常看到作者使用各式各樣的引理(Lemma),命題(Proposition),定理(Theorem),初次接觸時並不容易區別彼此之間的差別,但身為讀者該如何區分其輕重?

以下給出一般的觀點:

最重要的結果,我們一般將其稱作 定理 (Theorem) ,接續該定理所衍生的立即應用,我們稱做 系理(Corollary)。另外用以證明上述定理的小結果稱之為 引理 (Lemma),讀者可以將引理視為 "引出" 定理的 前置結果。

另外如果不是太過重要的結果,則作者一般會將該結果稱作 命題 (Proposition)。最後稍微提一下有時候作者會採用 宣稱 (Claim) 這是最弱的結果,一般不太重要。

Comments:
1. 上述只是一般觀點,並不代表所有引理都不重要,比如說 機率論與 測度論中 最重要的三大收斂性結果之一的 Fatou's Lemma 本身儘管是 引理,但卻是極為重要的結果。(另外兩個重要的收斂性結果為 Dominated Convergence Theorem 與 Monotone Convergence Theorem)

2. 一般而言, 定理是作者認為最為重要的結果,所以讀者在閱讀文獻時,若遇到定理可多留心這個部分。


2/16/2016

[投資理論] 數學能擊敗金融市場嗎?

以下為個人在 University of Wisconsin-Madison 臺灣學生會 2016年 第一場學術沙龍中 分享的簡報:

數學能擊敗金融市場嗎?-從控制理論觀點 (2016, 02. 16)

講者:謝宗翰
講題:數學是否能擊敗金融市場?-從控制理論觀點
簡介:此講題將試圖回答一個基本問題:是否存在一種「必勝法」,使得投資績效具備恆正報酬?我們將從現代投資理論出發,最終止於財務工程與控制理論,過程中,我們將逐步揭示何時可以透過數學幫助我們建構一組可行的「最佳」交易策略。


關於 UW-Madison 臺灣學生會 連結
https://sites.google.com/site/satuwmadison/

[最佳化] C^2 函數一階逼近的餘項積分表示

令 $f: \mathbb{R}^m \to \mathbb{R}$ 為 $C^2$-函數。對 $f$ 在 $y$ 附近使用一階泰勒展開: \[ T_y(x) := f(y) + \nabla f(y)^\top (x - y) \] 則其餘項 $R(x,y)$ 訂為 $$R(...