顯示具有 動態規劃 標籤的文章。 顯示所有文章
顯示具有 動態規劃 標籤的文章。 顯示所有文章

11/06/2024

[隨筆] A+焦慮的世代

接住A+世代學生

當了老師之後發現要"接住"學生確實不容易,撇開老師自身可能也有需要被接住的問題不談。我這幾年常常感受到這世代的學生們有著很大的徬徨,不太清楚未來的方向,但是卻有著非得要拿到A/A+不可的糾結,於是課優先選甜涼課,實習競賽投好投滿。好像看著同學做了些什麼,自己也要趕快跟上,深怕落後。想到很多年前讀過的太傻十日談,裡面提到許多跟現在當前學生們的焦慮幾乎如出一徹。

對我而言,拿A+大概意思是對於該科目"學得還行",對於內容掌握了一定的基礎,但並非完全精通該課程。學海無崖,正如聖經哥林多前書8:2 説到:
若有人以為自己知道什麼,按他所當知道的,他仍是不知道---哥林多前書8:2
這也是我常常提醒我自己的一句話。身為老師總是希望學生們別讓分數成為唯一的目標,能更深入地去理解和應用所學知識。這樣的觀點也許不符合學生的"期待",但這確實是我一直努力的方向。

職業與學術之間:適合比起優異更加重要

為什麼要拿A+? 或許是為了有漂亮成績單找工作/留學/申請研究所。但是這些知名公司/頂尖名校想要什麼樣的人選?難道是實習超多,獎項拿滿,GPA滿分的人才有機會嗎?我認為對於公司/名校而言真正重視的是想要找到"適合"的人,找到符合該職位或者學術需求的人。那麼適合那個職位的人到底需要些什麼技能我想才是需要被關注的地方。對於公司來說,他們往往重視應聘者能否在特定情況下解決實際問題,或者能否快速融入環境,快速學習。對於留學/申請研究所而言,展示自己具備研究能力,擁有學術熱誠的人,會遠比一張GPA滿分的成績單更具說服力。那麼這些能力該如何具體展示才是重點。



個人選擇 vs 個人夢想的不一致性:

很多人會說要找到自己喜歡的目標/領域,但是多數同學並不清楚自己喜歡的是什麼?有時候礙於家人環境同儕甚或經濟上的壓力,迫不得已折衷選了一條相對妥協的路。如果有機會,我常會建議:
選幾個至少不討厭的目標/領域試試看如何?

這是我常常跟不清楚未來的同學說的一句話,如果不知道自己要做什麼。至少選些"不討厭"的目標或者領域試試看。這可以讓你知道是否願意繼續深入還是提早退出。那麼該怎麼試?如果大學期間,除了本系之外的課程,涉獵自己有興趣的跨領域課程是一種方法,聽聽看各種學/業界的演講,參加有興趣的行業實習,參加社團/比賽也是一種方法,跟老師做做專題也是一種方法。八萬四千法門,選擇一種試試,時間允許的情況下,我想有心人永遠會找到一絲線索。但是永遠記得不是看別人做什麼就跟著做,而是先想清楚自己想要試試看什麼再做也不遲。


以終為始:動態規劃的最優解

如果把你想要的結果給你,你會不會愉快地接受?面對的每個決策,都有代價。確定目標之後,很多人可能會急著從現在開始一步步向目標前進,但是從現在的你看過去,可能的路線何其多,你怎麼知道哪條最正確呢?如果能練習從結果回推,往回看看前一小步所需的條件是什麼或許會更加有效。比如說你設立的一個長期目標:想出國留學,那麼前一步也許是決定領域?未來五年這個領域如何?哪所學校的研究方向更加適合?需要哪些語文能力/研究經驗/特定技能?每一步都從結果回推需求,這樣的動態規劃方式會讓你的前進方向更加具體更符合實際。


原則:尋找你自己

社群媒體的上的多數文章光鮮亮麗,滿手offer的貼文不免讓人稱羨,更甚者不免讓人失去自我,陷入攀比,但是年輕人,"你要保守你的心,勝過一切,因為一生的果效是由心發出的"--箴言4:32。別人的成功之路總是難以複製,那麼問問我們自己能不能勇敢走出屬於自己的路?

最後用詩人 Robert Frost 的 The Road Not Taken 作為結尾:

Two roads diverged in a wood, and I—I took the one less traveled by,
And that has made all the difference.


選擇一條少有人走的路,因為可能他帶給你最美的風景。並且譜寫成只屬於你的故事

後記:我大概有點資格談談這件事, 我本身研究橫跨隨機控制,優化理論與財務工程,走在這條路上的同行者寥寥。但是一路走來,我慶幸自己能有機會欣賞這條少有人走的道路以及沿途的特有風景。



5/06/2011

[最佳控制] 離散時間 穩態 LQR 控制問題 (1)

延續前篇,這次要介紹的是 Discrete Time Linear Quadratic Regulator in Infinite Horizon 或稱 Steady State LQR。

================
LQR Problem (Infinite Horizon LQR):
考慮離散狀態方程:
\[
x(k+1) = A x(k) + B u(k)
\]其中 $x(k) \in \mathbb{R}^n, A\in \mathbb{R}^{n \times n}, B \in \mathbb{R}^{n \times m}, u(k) \in \mathbb{R}^{m \times 1}$且 $(A,B)$ controllable。
定義 Performance index:
\[
J(u) = \displaystyle \sum_{k=0}^{\infty} x^T(k+1) Q x(k+1) + u^T(k) R u(k)
\] 其中 $Q, R$ 必須滿足 $Q^T = Q, Q \succ 0$, $R^T = R, R \succ 0$。 (亦即 $Q, R$ 必須為 對稱 + 正定 矩陣)

試求出一組最佳控制力序列 $u^*$ 使得成本函數 $J(u)$ 最小。
================

Comment:
讀者須注意到 Infinite Horizon 的 LQR問題要求計算 Performance index 為無窮級數和,此解必須保證收斂。以下定理告訴我們何時 此 Performance index 收斂

Lemma
考慮離散系統 $x(k+1) = A x(k) + B u(k)$,若 $(A,B)$ 可控制,且選 $Q, R >0$ 為正定矩陣,則上述 infinite horizon LQR 問題保證 閉迴路系統 狀態收斂到 $0$ 且 cost 為有界。

Proof: omitted. (see J. B. Rawlings and D. Q. Mayne, "Model Predictive Control: Theory and Design, p. 24", 2009)


現在我們可以開始求解 Infinite Horizon LQR問題:
Solution
回憶 Steady State Bellman Equation,為了符號簡便起見,我們寫成 functional equation 形式,
\[
I(x) = \displaystyle \min_{u \in \Omega} \{J(x,u) + I(f(x,u)) \}
\] 上式中 $J(x,u)$ 為 Branch cost,亦即 $J(x,u) = x^T Q x + u^T R u$ (並非 $\sum_{k=0}^{\infty} (\cdot)...$)

首先我們猜一組解 $I(x) = x^T P x$ 且矩陣 $P$ 為對稱正定矩陣,亦即滿足 $P^T = P, P \succ 0$。我們之後會找到此 $P$ 應該長甚麼樣子。

將猜測的解代入上述的 Steady State Bellman Equation,故現在我們得到
\[
I(x) = \min_{u \in \Omega} \{J(x,u) + I(f(x,u)) \}
\]注意到 $I(f(x,u) = f(x,u)^T P f(x,u) = (Ax+Bu)^TP(Ax+Bu)$,故我們可得
 \[
\begin{array}{l} I(x) = \mathop {\min }\limits_{u \in \Omega } \{ J(x,u) + I(f(x,u))\} \\ \Rightarrow {x^T}Px = \mathop {\min }\limits_u \left\{ {{x^T}Qx + {u^T}Ru + {{\left( {Ax + Bu} \right)}^T}P\left( {Ax + Bu} \right)} \right\}\\ \Rightarrow {x^T}Px = \mathop {\min }\limits_u \left\{ {{x^T}\left( {Q + {A^T}PA} \right)x + 2{x^T}{A^T}PBu + {u^T}Ru + {u^T}{B^T}PBu} \right\} \end{array}
\]透過一階必要條件 FONC: $ \frac{\partial }{{\partial u}} = 0$ 對上式右邊求解
 \[\begin{array}{l} 2{\left( {{x^T}{A^T}PB} \right)^T} + 2Ru + 2{B^T}PBu = 0\\ \Rightarrow {u^*} = - {\left( {R + {B^T}PB} \right)^{ - 1}}{B^T}PAx \end{array}
\]現在將 $u^*$ 代回 $(*)$  可得 \[\begin{array}{l} {x^T}Px = \mathop {\min }\limits_u \left\{ {{x^T}\left( {Q + {A^T}PA} \right)x + 2{x^T}{A^T}PBu + {u^T}Ru + {u^T}{B^T}PBu} \right\}\\ \Rightarrow {x^T}Px = \left\{ {{x^T}\left\{ {Q + {A^T}PA - {A^T}PB{{\left( {R + {B^T}PB} \right)}^{ - 1}}{B^T}PA} \right\}x} \right\} \end{array}
\]比較左右兩邊可得到 $P$ 必須滿足下式: \[P = Q + {A^T}PA - {A^T}PB{\left( {R + {B^T}PB} \right)^{ - 1}}{B^T}PA\] 此式稱為 Discrete Time Algebraic Ricatti Equation (ARE),一般而言,可利用 MATLA 指令 dare(A,B,Q, R) 求解 P。

由於 $u^* = - {\left( {R + {B^T}PB} \right)^{ - 1}}{B^T}PAx$,其中除了 $P$ 未定之外,其餘所需要的參數都已知且皆與跌代時間無關,故此無窮時間LQR問題得到的 最佳控制力為 Time invariant。

現在我們總結如下:求解無窮時間的LQR問題只要做兩個步驟即可
STEP 1: 求解一次 Algebraic Ricatti Equation 得到 $P$ (利用 MATLAB: dare.m 或者徒手計算)
\[
P = Q + {A^T}PA - {A^T}PB{\left( {R + {B^T}PB} \right)^{ - 1}}{B^T}PA
\]STEP2 : 將 $P$ 代入 ${u^*} =  - {\left( {R + {B^T}PB} \right)^{ - 1}}{B^T}PAx$

下面我們看個例子:

Example:
考慮一個離散時間線性系統狀態方程:
\[
x_1(k+1) = x_2(k) \\
x_2(k+1) = x_1(k) + u(k)
\]且考慮 Cost function:
\[
J = \sum_{k=0}^{\infty}2x_1^2(k) + 2x_1(k)x_2(k) + x_2^2(k) + 3u^2(k)
\] 且控制力具有如下形式:
\[
u(k) = K_1 x_1(k) + K_2 x_2(k)

\] 試求 $K_1, K_2$ 使 上述 Cost function 最小:

Solution
首先定義  $x(k) := [x_1(k), x_2(k)]^T$ ,則我們有
\[x\left( {k + 1} \right) = \underbrace {\left[ {\begin{array}{*{20}{c}}
0&1\\
1&0
\end{array}} \right]}_A\left[ {\begin{array}{*{20}{c}}
{{x_1}(k)}\\
{{x_2}(k)}
\end{array}} \right] + \underbrace {\left[ {\begin{array}{*{20}{c}}
0\\
1
\end{array}} \right]}_Bu(k)
\] 與 Cost function
\[\begin{array}{l}
J = \sum\limits_{k = 0}^\infty {(2x_1^2(} k) + 2{x_1}(k){x_2}(k) + x_2^2(k) + 3{u^2}(k))\\
\begin{array}{*{20}{c}}
{}
\end{array} = x{\left( k \right)^T}\underbrace {\left[ {\begin{array}{*{20}{c}}
2&1\\
1&1
\end{array}} \right]}_Qx\left( k \right) + \underbrace 3_R{u^2}(k)
\end{array}
\] 那麼現在此問題變成 Steady-state LQR problem,故由前述討論可知我們有 Optimal feedback control 為
\[
u^*(k) = -(R+B^T P B)^{-1} B^T PA \cdot x(k)

\] 其中 $P$ 滿足 $P=P^T, P \succ 0$ 可由 ARE
\[
P= A^TPA - A^T PB (R+ B^TPB)^{-1}B^TPA+Q

\]利用 MATLAB 指令 dare(A,B,Q,R) 解得 $P = \left[ {\begin{array}{*{20}{c}}
{3.7841}&{1.6815}\\
{1.6815}&{4.4022}
\end{array}} \right]$ 現在將 $P$ 帶回 $u^*$中
\[u\left( k \right) = - \left[ {\begin{array}{*{20}{c}}
{0.5947}&{0.2272}
\end{array}} \right]x\left( k \right)

\]i.e., $K_1 = -05947, K_2 =-0.2272$. $\square$


5/05/2011

[最佳控制] 離散時間 LQR- Finite Time Horizon

這次要介紹 控制理論中一個重要的結果:
 Discrete Time Linear Quadratic Regulator (LQR) in Finite Time Horizon,
中文翻譯為 離散時間線性二次調節器,我們這邊會針對此問題利用 Dynamic Programming 的方法來逐步求解

================
LQR Problem (Finite Horizon):
考慮狀態方程:
\[
x(k+1) = A x(k) + B u(k)
\]其中 $x(k) \in \mathbb{R}^n, A\in \mathbb{R}^{n \times n}, B \in \mathbb{R}^{n \times m}, u(k) \in \mathbb{R}^{m \times 1}$
且 考慮 Performance index:
\[
J(u) = \displaystyle \sum_{k=0}^{N-1} x^T(k+1) Q x(k+1) + u^T(k) R u(k)
\] 其中 $Q, R$ 必須滿足 $Q^T = Q, Q \succ 0$, $R^T = R, R \succ 0$。 (亦即 $Q, R$ 必須為 對稱 + 正定 矩陣)

試求出一組最佳控制力序列 $u(N-1), u(N-2),... u(0)$ 使得成本函數 $J(u)$ 最小。
================

Comment:
1. LQR 顧名思義是其具有系統狀態方程為線性 $x(k+1) = Ax(k)+Bu(k)$ 與  Performance Index 中的項都為二次式。
\[
J(u) = \displaystyle \sum_{k=0}^{N-1} x^T(k+1) Q x(k+1) + u^T(k) R u(k)
\]
2. 上述對於矩陣 $Q, R$ 的對稱與正定假設是必須的 (之後需求在求解最佳控制力序列的時候需要求解反矩陣,故需要這些性質。)

3. 上述 LQR in Finite Horizon的問題所求得的最佳控制力序列為 Time Varying 。(此性質會在下面求解的時候再度強調。)

4. 若 Performance index 考慮 $N \rightarrow \infty$,亦即
\[
J(u) = \displaystyle \sum_{k=0}^{\infty} x^T(k+1) Q x(k+1) + u^T(k) R u(k)
\]則我們說此問題為 LQR in Infinite Horizon 或稱 Steady State LQR。這類問題會在之後再做介紹。(此類問題所得到的最佳控制序列將不再是 Time Varying,且只需求解一次 Algebraic Ricatti Equation 即可獲得最佳控制力序列)



Solution of Discrete Time LQR problem in Finite Horizon 
現在我們開始進行求解。

這邊我們使用 Dynamic Programming 方法來求解上述 LQR問題。回憶 Dynamic Programming Equation 的定義:
\[I\left( {x\left( l \right),N - l} \right) = \mathop {\min }\limits_{u\left( l \right) \in {\Omega _l}} \left\{ {J\left( {x\left( l \right),u\left( l \right)} \right) + I\left( {x\left( {l + 1} \right),N - \left( {l + 1} \right)} \right)} \right\}
\] 考慮 Optimal Cost of one-step-to-go ($l=N-1$):
\[\begin{array}{l}
 \Rightarrow I\left( {x\left( {N - 1} \right),1} \right) = \min \left\{ {J\left( {x\left( {N - 1} \right),u\left( {N - 1} \right)} \right) + I\left( {x\left( N \right),0} \right)} \right\}\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}&{}&{}
\end{array} = \min \left\{ {{x^T}(N)Qx(N) + {u^T}(N - 1)Ru(N - 1)} \right\}
\end{array}
\]將系統狀態方程 $x(k+1) = A x(k) + B u(k) $ 代入:
\[\begin{array}{l}
 \Rightarrow I\left( {x\left( {N - 1} \right),1} \right) = \\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array}\min \left\{ \begin{array}{l}
{\left[ {Ax\left( {N - 1} \right) + Bu\left( {N - 1} \right)} \right]^T}Q\left[ {Ax\left( {N - 1} \right) + Bu\left( {N - 1} \right)} \right]\\
 + {u^T}(N - 1)Ru(N - 1)
\end{array} \right\}\\
 \Rightarrow I\left( {x\left( {N - 1} \right),1} \right) = \min \left\{ \begin{array}{l}
{x^T}\left( {N - 1} \right){A^T}QAx\left( {N - 1} \right)\\
\begin{array}{*{20}{c}}
{}
\end{array} + 2{x^T}\left( {N - 1} \right){A^T}QBu\left( {N - 1} \right)\\
\begin{array}{*{20}{c}}
{}
\end{array} + {u^T}\left( {N - 1} \right){B^T}QBu\left( {N - 1} \right)\\
\begin{array}{*{20}{c}}
{}
\end{array} + {u^T}(N - 1)Ru(N - 1)
\end{array} \right\}
\end{array}
\] 由一階必要條件 FONC: ${\nabla _{u\left( {N - 1} \right)}}I\left( {x\left( {N - 1} \right),1} \right) = 0$ 可知
\[\begin{array}{l}
2{B^T}{Q^T}Ax\left( {N - 1} \right) + 2{B^T}QBu\left( {N - 1} \right) + 2Ru(N - 1) = 0\\
 \Rightarrow \left[ {R + {B^T}QB} \right]u(N - 1) =  - {B^T}{Q^T}Ax\left( {N - 1} \right)\\
 \Rightarrow u^*(N - 1) =  - \underbrace {{{\left[ {R + {B^T}QB} \right]}^{ - 1}}{B^T}QA}_{K\left( {N - 1} \right)}x\left( {N - 1} \right)
\end{array}
\]接著我們把上述求得的 1-step-to-go 的最佳控制力 $u^*(k)$ 代回 Optimal cost of 1-step-to-go,我們得到 (透過一些代數運算)
\[I\left( {x\left( {N - 1} \right),1} \right) = {x^T}\left( {N - 1} \right)\underbrace {{A^T}\left\{ {Q - QB{{\left[ {R + {B^T}QB} \right]}^{ - 1}}^T{B^T}Q} \right\}A}_{: = M\left( {N - 1} \right)}x\left( {N - 1} \right)
\]在做下一步跌代之前我們先給個 comments

Comments :
1.注意到上式中我們求 ${u^*}(N - 1) =  - {\left[ {R + {B^T}QB} \right]^{ - 1}}{B^T}QAx\left( {N - 1} \right)$ 需要計算反矩陣 ${\left[ {R + {B^T}QB} \right]^{ - 1}}$ ,故需檢驗反矩陣是否存在:不過這個問題可以被證明反矩陣確實存在: (因為 $R$ 為對稱正定矩陣,$B^TQB$ 為對稱半正定矩陣,由FACT: 正定矩陣+ 半正定矩陣 = 正定矩陣,且正定矩陣具有 eigenvalue 全為正,故可推知反矩陣存在。)

2. 上式所求得的最佳控制力 ${u^*}(N - 1) =  - {\left[ {R + {B^T}QB} \right]^{ - 1}}{B^T}QAx\left( {N - 1} \right) =  - K\left( {N - 1} \right)x\left( {N - 1} \right)$ 其中的 $K_{N-1}$ 稱作 Optimal Gain Matrix。此控制力具有回授控制形式 Feedback form。 (類似 $u=-Kx$ 的形式)

3. 最佳控制力儘管具有回授控制 (Feedback form)的形式,但在此問題中為時變得 Gain。此性質在下一步跌代會顯現出來。

4. 當我們有了上述第一步跌代的結果,整個 LQR問題就簡單許多,因為之後的跌代只有 $Q$ 矩陣會改變其餘結果均不變。我們直接看下一步跌代便會發現此性質:

Back to computation :
考慮 Optimal Cost of 2-steps-to-go: ($l=N-2$):
\[\begin{array}{l}
I\left( {x\left( {N - 2} \right),2} \right) = \min \left\{ {J\left( {x\left( {N - 2} \right),u\left( {N - 2} \right)} \right) + I\left( {x\left( {N - 1} \right),1} \right)} \right\}\\
 \Rightarrow I\left( {x\left( {N - 2} \right),2} \right)\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} = \min \left\{ \begin{array}{l}
{x^T}(N - 1)Qx(N - 1) + {u^T}(N - 2)Ru(N - 2)\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} + {x^T}\left( {N - 1} \right)M\left( {N - 1} \right)x\left( {N - 1} \right)
\end{array} \right\}\\
 \Rightarrow I\left( {x\left( {N - 2} \right),2} \right)\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} = \min \left\{ {{x^T}(N - 1)\underbrace {\left[ {Q + M\left( {N - 1} \right)} \right]}_{: = \tilde Q}x\left( {N - 1} \right) + {u^T}(N - 2)Ru(N - 2)} \right\}
\end{array}
\]可以發現上式中只有 $Q$ 變成了 $\tilde{Q}$ 其餘參數均固定。故我們可以馬上寫下對應的最佳控制力 $u^*(N-2)$:
\[{u^*}(N - 2) =  - {\left[ {R + {B^T}\tilde QB} \right]^{ - 1}}{B^T}\tilde QAx\left( {N - 2} \right) = K (N-2) x(N-2)
\]其對應的 Optimal cost of 2 steps to go:
\[I\left( {x\left( {N - 2} \right),2} \right) = {x^T}\left( {N - 2} \right)\underbrace {{A^T}\left\{ {\tilde Q - \tilde QB{{\left[ {R + {B^T}\tilde QB} \right]}^{ - 1}}^T{B^T}\tilde Q} \right\}A}_{: = M\left( {N - 2} \right)}x\left( {N - 2} \right)
\]故我們只需要持續重複上述步驟,即可依序解得 $K(N-3), K(N-4),..., K(0)$ 與 $M(N-3), M(N-4),...,M(0)$ 與 $J^*$


Comment:
注意到 $u^*(N-2)$ 的 Gain $K(N-2)$ 與 $u^*(N-1)$ 的 Gain $K(N-1)$ 並不相同,此說明了之前所敘的時變特性(Time Varying Property)

5/05/2010

[動態規劃] 淺談 離散時間動態規劃 (2) -Minimax Dynamic programming

考慮系統狀態方程:
\[
x(k+1) = f(x(k), u(k))
\]現在引入一個 監測函數 (monitor function): $g(x(k))$ 此函數用來擷取我們所關心的系統狀態。

現在我們定義 cost function
\[
J(u) := \displaystyle \min_{k=1,2,...,N} g(x(k))
\]
我們的目標:
藉由選取最佳的 $u(k) \in \Omega_k$ 使得 $\max J(u)$,亦即
\[
\max \displaystyle \min_{k=1,2,...,N} g(x(k))
\]

此時我們的 Bellman equation 需要進行修正:
\[
I(x(l), N-l) = \max_{u(l) \in \Omega_l} \left \{ \min \left \{ g\left( x(l+1), I(x(l+1), N-(l+1)) \right ) \right \} \right\}
\]

Example
考慮如下狀態方程:
\[x\left( {k + 1} \right) = \left[ {\begin{array}{*{20}{c}}
{\frac{1}{2}}&{\frac{1}{4}}\\
{ - \frac{1}{2}}&{\frac{1}{4}}
\end{array}} \right]x\left( k \right) + \left[ {\begin{array}{*{20}{c}}
1\\
1
\end{array}} \right]u\left( k \right)
\] 其中 $x(k) := [x_1(k) \ x_2(k)]^T$;
現在,僅考慮 $N=2$ (只有兩步)的情況,另外控制力 $u(0), u(1) \geq 0$ 且需滿足如下拘束: $u(0) + u(1) \leq 1$,定義監測函數
\[
g(x(k)) = x_1(k)
\]
試求一組 $u^*(0), u^*(1)$ maximize the Floor of $g((k))$。

Solution
首先考慮 Optimal cost of one-step-to-go:  ($l=N-1 =1$)
\[
\begin{array}{l}
I\left( {x\left( 1 \right),1} \right) = \mathop {\max }\limits_{u\left( 1 \right) \in {\Omega _1}} \left\{ {g\left( {x\left( 2 \right)} \right)} \right\} = \mathop {\max }\limits_{u(1) \in \left[ {0,1 - u(0)} \right]} \left\{ {{x_1}\left( 2 \right)} \right\}\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}&{}&{}
\end{array} = \mathop {\max }\limits_{u(1) \in \left[ {0,1 - u(0)} \right]} \left\{ {\frac{1}{2}{x_1}\left( 1 \right) + \frac{1}{4}{x_2}\left( 1 \right) + u\left( 1 \right)} \right\}
\end{array}
\]觀察上式可知欲得到 Maximum,我們需選
\[{u^*}(1) = 1 - u(0)
\]故將上述 $u^*(1)$ 代回,我們可得對應的 Optimal cost to go :
\[I\left( {x\left( 1 \right),1} \right) = \frac{1}{2}{x_1}\left( 1 \right) + \frac{1}{4}{x_2}\left( 1 \right) + 1 - u(0)
\]現在我們計算 Two-steps-to-go: ($l=N-2 =0$)
\[\begin{array}{l}
I\left( {x\left( l \right),N - l} \right) = \mathop {\max }\limits_{u\left( l \right) \in {\Omega _l}} \left\{ {\min \left\{ {g\left( {x\left( {l + 1} \right)} \right),I\left( {x\left( {l + 1} \right),N - l - 1} \right)} \right\}} \right\}\\
 \Rightarrow I\left( {x\left( 0 \right),2} \right) = \mathop {\max }\limits_{u\left( 0 \right) \in {\Omega _0}} \left\{ {\min \left\{ {{x_1}\left( 0 \right),I\left( {x\left( 1 \right),1} \right)} \right\}} \right\}\\
 \Rightarrow I\left( {x\left( 0 \right),2} \right) = \mathop {\max }\limits_{u\left( 0 \right) \in \left[ {0,1} \right]} \left\{ {\min \left\{ {{x_1}\left( 0 \right),\frac{1}{2}{x_1}\left( 1 \right) + \frac{1}{4}{x_2}\left( 1 \right) + 1 - u(0)} \right\}} \right\} \ \ (*)
\end{array}
\] 由系統狀態方程可知
\[\left\{ \begin{array}{l}
{x_1}\left( 1 \right) = \frac{1}{2}{x_1}\left( 0 \right) + \frac{1}{4}{x_2}\left( 0 \right) + u\left( 0 \right)\\
{x_2}\left( 1 \right) =  - \frac{1}{2}{x_1}\left( 0 \right) + \frac{1}{4}{x_2}\left( 0 \right) + u\left( 0 \right)
\end{array} \right.
\]故將上述結果代回 $(*)$ 可得
\[ \Rightarrow I\left( {x\left( 0 \right),2} \right) = \mathop {\max }\limits_{u\left( 0 \right) \in \left[ {0,1} \right]} \left\{ {\min \left\{ {A + u\left( 0 \right),B - \frac{1}{4}u\left( 0 \right)} \right\}} \right\}
\] 其中
\[\left\{ \begin{array}{l}
A = \frac{1}{2}{x_1}\left( 0 \right) + \frac{1}{4}{x_2}\left( 0 \right)\\
B = 1 + \frac{1}{8}{x_1}\left( 0 \right) + \frac{3}{{16}}{x_2}\left( 0 \right)
\end{array} \right.
\] 為了找出 Maximum,我們定義 $F(u(0)) :={\min \left\{ {A + u\left( 0 \right),B - \frac{1}{4}u\left( 0 \right)} \right\}}$ 那麼考慮以下兩種情況:
CASE 1: $A \geq B$:則 $u^*(0)=0, u^*(1) = 1$
CASE 2: $A <B$ :則
\[{u^*}\left( 0 \right) = \left\{ \begin{array}{l}
\frac{4}{5}\left( {B - A} \right),\begin{array}{*{20}{c}}
{}
\end{array}if\begin{array}{*{20}{c}}
{}
\end{array}\frac{4}{5}\left( {B - A} \right) \in \left[ {0,1} \right]\\
1,\begin{array}{*{20}{c}}
{}
\end{array}otherwise
\end{array} \right.\] 且 $u^*(1) = 1 - u^*(0)$

5/07/2009

[動態規劃] 淺談 離散時間動態規劃 (1) - Bellman equation in Infinite Horizon

這次要介紹 無窮時間的 Bellman Equation,亦即我們的 Cost function 為 branch cost 加到無窮大的情況:
\[
J(u) := \displaystyle \sum_{k=0}^{\infty} J(x(k), u(k)) + \Phi(x(N))
\]
那麼現在問題變成 儘管我們手上有 有限時間的 Bellman Equation (請參考前篇文章),但對於此類無窮時間的問題該如何處理??

亦即如果今天我們 cost function 的 $N \rightarrow \infty$ 該怎麼處理? 我們稱這一類問題叫做 Steady State Dynamic Programming 或稱 Dynamic Programming in Infinite Horizon。

無窮時間的動態規劃問題 (Dynamic Programming Problem in Infinite Horizon):

考慮 Performance index (cost function)
\[
J(u) := \displaystyle \sum_{k=0}^{\infty} J(x(k), u(k)) + \Phi(x(N))
\]狀態方程( state equation)
\[
x(k+1) = f(x(k), u(k)), \ x(0) \ \text{is given} \\
\]其中 $x(k)$ 為系統在第 $k$ 時刻的 狀態,$f:\mathbb{R}^n \times \mathbb{R}^m \rightarrow \mathbb{R}^n$,
與控制力拘束條件
\[
u(k) \in \Omega
\]其中$\Omega$ 為拘束條件

此時對應的 Bellman Equation (or Dynamic Programming Equation) 寫為如下的 functional form
\[
I(x) = \min_{u \in \Omega} \{ J(x,u) + I(f(x,u)) \}
\]亦即與 $N$ 無關
上式稱為 Steady State Bellman Equation 或者 Bellman Equation in Infinite Horizon。


現在我們看個例子:

Example
考慮系統狀態方程表示如下:
\[x(k+1) = x(k) - u(k)
\]且 cost function 為
\[
J(u) = \displaystyle \sum_{k=0}^{\infty} (x(k+1)-u(k))^2 + x^2(k+1)
\]試求 Optimal $u^*$ 與其對應的 optimal cost to go

Solution
首先我們寫下 Steady State Bellman Equation:
\[
I(x) = \min_{u \in \Omega} \{J(x,u) + I(f(x,u)) \}
\]現在我們要求解上式,故我們首先猜一個解 (事實上此問題在之前文章中我們解過有限時間的 Bellman equation,當時我們解出 Optimal cost to go 具有 $I(x) = \alpha x^2$ 的形式,有興趣的讀者可前往前篇文章檢驗),故我們現在很合理的可以猜這 $I(x) = \alpha x^2$ 其中 $\alpha$ 為代定係數 $\alpha \geq 0$。
故我們將此解代入 Steady State Bellman Equation,可得
\[\begin{array}{l}
I(x) = {\min _{u \in \Omega }}\{ J(x,u) + I(f(x,u))\} \\
 \Rightarrow \alpha {x^2} = {\min _{u \in \Omega }}\{ {(x - u)^2} + {x^2} + \alpha {\left( {x - u} \right)^2}\} \\
 \Rightarrow \alpha {x^2} = {\min _{u \in \Omega }}\{ {x^2} + \left( {1 + \alpha } \right){\left( {x - u} \right)^2}\}
\end{array}
\]由 FONC: $\frac{\partial }{{\partial u}} = 0$,我們可求解上式右方最佳控制力 $u^*$ :
\[
u^* = x
\]故代回上式我們得到
\[\begin{array}{l}
\alpha {x^2} = {x^2}\\
 \Rightarrow \alpha  = 1
\end{array}
\] 故如果我們選 $\alpha =1$即可解得 Steady-State Bellman Equation。 $\square$

上述問題可以進一步拓展到 $\mathbb{R}^n$ 空間,這一類問題在控制理論中稱為無窮時間的線性二次調節器 (Linear Quadratic Regulator, LQR) 問題,有興趣的讀者請參考:
[最佳控制] 離散時間 穩態線性二次調節器 Discrete Time Linear Quadratic Regulator in Infinite Horizon via Dynamic Programming


4/28/2009

[動態規劃] 淺談 離散時間動態規劃 (0) - Principle of Optimality & Bellman equation in Finite Horzion

這次要跟大家介紹離散時間的動態規劃 (Dynamic Programming, D.P.) 問題:
此為最佳化理論中的一種方法,由 Professor Bellman 的 Principle of Optimality 所奠基。

基本想法如下: 從最佳解往回推 (backward),逐步求解最佳值。再從初始位置沿著以求得的最佳路徑前進到最佳解。

我們現在把動態規劃問題寫下:

(有限時間) 動態規劃問題 (Dynamic Programming Problem in Finite Horizon):

考慮 Performance index (cost function)
\[
J(u) := \displaystyle \sum_{k=0}^{N-1} J(x(k), u(k)) + \Phi(x(N))
\]與離散時間的狀態方程( state equation)
\[
x(k+1) = f(x(k), u(k)), \ x(0) \ \text{is given} \\
\]其中 $x(k) \in \mathbb{R}^n$ 為系統在第 $k$ 時刻的 狀態,$u(k) \in \mathbb{R}^m$ 為控制力,函數 $f:\mathbb{R}^n \times \mathbb{R}^m \rightarrow \mathbb{R}^n$,

另外我們仍需考慮控制力拘束條件
\[
u(k) \in \Omega_k(x(k))
\]其中 $\Omega_k$ 為第 $k$ 時刻的拘束條件

Comment:
1. 關於狀態方程:上式中的 $f(x(k),u(k))$ 稱作 autonomous state equation
但事實上DP問題亦可直接考慮 Non-autonomous state equation,(亦即函數中考慮時間 $k$)
\[
f((x(k), u(k), k))
\] Example: Autonomous and Non-autonomous system
Autonomous state equation:  $x(k+1) = x(k) + u(k)$
Non-autonomous state equation: $x(k+1) = k \cdot x(k)$

2. 關於 Performance Index
\[
J(u) := \displaystyle \sum_{k=0}^{N-1} J(x(k), u(k)) + \Psi (x(N))
\] 其中 $J(x(k), u(k))$ 稱為 Branch cost。 $\Psi (x(N))$ 稱為 Terminal cost。


現在我們介紹 Bellman's Principle of Optimality

==========================
Principle of Optimality
假設我們有一組 optimal sequence
\[
u^*(0), u^*(1), ..., u^*(N-1)
\]令 $x^*(0), x^*(1), ..., x^*(N)$ 為 induced states,亦即對 $k=0,1,...,N-1$
\[
u(k) \in \Omega_k(x^*(k))
\]且 $x^*(k+1) = f(x^*(k),u^*(k))$  (也就是說上述的 optimal sequence 最佳化整個問題。)
那麼現在考慮一個 子問題 (subproblem):
考慮 $k = l$ (現在考慮由 時刻 $l$ 開始 (不是從 $0$ 開始) ),初始狀態現在為 $x(l)=x^*(l)$ (亦即此時的初始狀態 $x(l)$ 落在最佳解上 $x^*(l)$),此時對應的 Cost function 為
\[
J_l(u) := \displaystyle \sum_{k=l}^{N-1} J(x(k), u(k)) + \Psi (x(N))
\]則 Principle of Optimality 告訴我們
\[
u^*(l), u^*(l+1), ..., u^*(N-1)
\]為對此 subproblem 的 optimal sequence.
==========================
Comment:
上述的 Principle of Optimality 其實表達的想法很簡單,就是如果已有一組最佳解序列在手,那麼從任意時刻開始的解也必須是最佳解。比如說我有一組最佳控制力序列
$u^*(0), u^*(1), ..., u^*(10)$ 那麼考慮從 $k=5, 6,...10$ 的子問題,則對應的控制力序列
\[
u^*(5), u^*(6),...u^*(10)
\]亦為此子問題的最佳解。

或者更淺白的說,假設今天要坐飛機到美國,最佳航程為
台灣 => 香港 =>美國
如果今天如果不從台灣開始,而是從香港開始(且香港已經在最佳路徑上),則到美國的這條路徑仍維持最佳解。(香港=>美國)


Proof: Bellman's Principle of Optimality

用矛盾證法 (Proof by contradiciton),
考慮一 子問題(subproblem) 為 從狀態 $x^*(l)$ 開始,且此狀態落在最佳解上,但 $u^*(k), \ k=l,l+1,...,N-1$ 並非 Optimal sequence。

既然上述不是 subproblem 的 Optimal sequence,則必定有其他人比此 $u^*$ 這組 sequence 可以得到更低的 cost,也就是說存在一組
\[
\hat{u}(l), \hat{u}(l+1),... \hat{u}(N-1)
\]使得 $J_l(\hat{u}) < J_l (u^*)$

現在我們定義 Optimal sequence: $\tilde{u}(k)$
\[\tilde u(k): = \left\{ \begin{array}{l}
{u^*}\left( k \right),\begin{array}{*{20}{c}}
{}
\end{array}k < l\\
\hat u\left( k \right),\begin{array}{*{20}{c}}
{}
\end{array}l \le k \le N - 1
\end{array} \right.
\] 那麼我們現在觀察 整個最佳化問題 (非子問題 ) 的 cost;亦即對 $k=0$ 開始 (而非 $k=l$),則我們有
\[\begin{array}{l}
J\left( {\tilde u} \right) = \sum\limits_{k = 0}^{N - 1} {J\left( {x\left( k \right),u\left( k \right)} \right)}  + \Psi \left( {x\left( N \right)} \right)\\
 \Rightarrow J\left( {\tilde u} \right) = \underbrace {\sum\limits_{k = 0}^{l - 1} {J\left( {x\left( k \right),u\left( k \right)} \right)} }_{{J_0}\left( {{u^*}} \right)} + \underbrace {\sum\limits_{k = l}^{N - 1} {J\left( {x\left( k \right),u\left( k \right)} \right)} }_{{J_l}\left( {\hat u} \right)} + \Psi \left( {x\left( N \right)} \right)\\
 \Rightarrow J\left( {\tilde u} \right) < {J_0}\left( {{u^*}} \right) + {J_l}\left( {{u^*}} \right) + \Psi \left( {x\left( N \right)} \right)\\
 \Rightarrow J\left( {\tilde u} \right) < J\left( {{u^*}} \right)
\end{array}
\] 上述結果說明 $\tilde{u}$ 的序列為 Optimal sequence,但此結果與 Principle of Optimality 中原本假設的 最佳序列 $u^*(0), u^*(1),..., u^*(N-1)$ 矛盾。故得證。 $\square$


有了上述的 Principle of Optimality,我們便可以建構 Dynamic Programming Equaiton (DPE) (DPE 又稱為 Bellman equation)


Bellman Equation
首先考慮我們的 Subproblem:
由初始狀態 $x(l)$,定義 Optimal cost-to-go 記作 $I(x(l), N-l)$ (亦即考慮 Subproblem 由 $x(l)$ 開始 到 $x(N)$ 的 optimal cost )
\[I(x(l),N - l): = \mathop {\min }\limits_{\scriptstyle u(k) \in {\Omega _k},\atop
\scriptstyle k \ge l} {J_l}\left( u \right)
\]由 $J_l(u)$ 定義可知,${J_l}(u): = \sum\limits_{k = l}^{N - 1} J (x(k),u(k)) + \Psi (x(N))$,故
\[
\begin{array}{l}
I(x(l),N - l): = \mathop {\min }\limits_{\begin{array}{*{20}{c}}
{u(k) \in {\Omega _k},}\\
{k \ge l}
\end{array}} \left[ {\sum\limits_{k = l}^{N - 1} J (x(k),u(k)) + \Psi (x(N))} \right]\\
 \Rightarrow I(x(l),N - l) = \mathop {\min }\limits_{\begin{array}{*{20}{c}}
{u(k) \in {\Omega _k},}\\
{k \ge l}
\end{array}} \left[ \begin{array}{l}
J(x(l),u(l))\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + \sum\limits_{k = l + 1}^{N - 1} {J(x(k),u(k))}  + \Psi (x(N))
\end{array} \right]
\end{array}
\]由 Principle of Optimality 可知,上式可改寫為
\[
\small{\begin{array}{l}
 \Rightarrow I(x(l),N - l): = \mathop {\min }\limits_{u(k) \in {\Omega _k}} \left[ \begin{array}{l}
J(x(l),u(l))\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + \mathop {\min }\limits_{\begin{array}{*{20}{c}}
{u(k) \in {\Omega _k},}\\
{k > l}
\end{array}} \left\{ {\sum\limits_{k = l + 1}^{N - 1} {J(x(k),u(k))}  + \Psi (x(N))} \right\}
\end{array} \right]\\
 \Rightarrow I(x(l),N - l): = \mathop {\min }\limits_{u(k) \in {\Omega _k}} \left[ {J(x(l),u(l)) + I\left( {x(l + 1),N - \left( {l + 1} \right)} \right)} \right]
\end{array}}
\]
上式稱為 Bellman Equation 或稱 Dynamic Programming Equation
\[
I(x(l), N-l) = \min_{u(l) \in \Omega_l} \{J(x(l), u(l)) + I(x(l+1), N-(l+1)) \}
\]

Example
考慮系統狀態方程表示如下:
\[x(k+1) = x(k) - u(k)
\]且 cost function 為
\[
J(u) = \displaystyle \sum_{k=0}^{N-1} (x(k+1)-u(k))^2 + x^2(k+1)
\]試求 Optimal $u^*(N-1), u^*(N-2),...$ 與其對應的 optimal cost to go

Solution
為簡便起見,這邊只求 $u(N-1)$ 與其對應的 optimal cost to go。

故首先由 $x(N-1)$ 開始,亦即 $l= N-1$,由 DPE 可知
\[\begin{array}{l}
I(x(l),N - l): = \mathop {\min }\limits_{u(k) \in {\Omega _k}} \left[ {J(x(l),u(l)) + I\left( {x(l + 1),N - \left( {l + 1} \right)} \right)} \right]\\
 \Rightarrow I(x(N - 1),1) = \mathop {\min }\limits_{u(k) \in {\Omega _k}} \left[ {J(x(N - 1),u(N - 1)) + I\left( {x(N),0} \right)} \right]\\
 \Rightarrow I(x(N - 1),1) = \mathop {\min }\limits_{u(k) \in {\Omega _k}} \left[ {{{\left[ {x(N) - u(N - 1)} \right]}^2} + {x^2}\left( N \right) + {\rm{0}}} \right]\\
 \Rightarrow I(x(N - 1),1) = \mathop {\min }\limits_{u(k) \in {\Omega _k}} \left[ {{{\left[ {x(N) - u(N - 1)} \right]}^2} + {x^2}\left( N \right)} \right]
\end{array}
\]現在再由系統狀態方程 $x\left( N \right) = x(N - 1) - u(N - 1) $ 代入上式可得
\[I(x(N - 1),1) = \mathop {\min }\limits_{u(N - 1)} \left[ {{{\left[ {x(N - 1) - 2u(N - 1)} \right]}^2} + {{\left[ {x(N - 1) - u(N - 1)} \right]}^2}} \right]
\] 由 FONC 求最佳解:
\[\begin{array}{l}
\frac{{\partial I(x(N - 1),1)}}{{\partial u(N - 1)}} = 0 \\
 \Rightarrow {u^*}(N - 1) = \frac{3}{5}x(N - 1)
\end{array}
\]其對應的Optimal cost-to-go
\[\begin{array}{l}
 \Rightarrow {\left. {I(x(N - 1),1)} \right|_{{u^*}}} = {\left[ {x(N - 1) - 2\frac{3}{5}x(N - 1)} \right]^2}\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}&{}&{}&{}&{}
\end{array} + {\left[ {x(N - 1) - \frac{3}{5}x(N - 1)} \right]^2}\\
 \Rightarrow {\left. {I(x(N - 1),1)} \right|_{{u^*}}} = \underbrace {\frac{1}{5}}_{{K_{N - 1}}}{x^2}(N - 1)
\end{array}
\] 上式子中的 $K_{N-1}$ 被稱作 feedback-gain at $N-1$ stage。

接著我們計算 $x(N-2)$ 亦即 $l=N-2$ 的 Optimal cost-to-go:
\[\begin{array}{*{20}{l}}
{I(x(l),N - l): = \mathop {\min }\limits_{u(k) \in {\Omega _k}} \left[ {J(x(l),u(l)) + I\left( {x(l + 1),N - \left( {l + 1} \right)} \right)} \right]}\\
{ \Rightarrow I(x(N - 2),2) = \mathop {\min }\limits_{u(N - 2) \in {\Omega _{N - 2}}} \left[ {J(x(N - 2),u(N - 2)) + I\left( {x(N - 1),1} \right)} \right]}\\
{ \Rightarrow I(x(N - 2),2) = \mathop {\min }\limits_{u(N - 2)} \left[ \begin{array}{l}
{\left( {x\left( {N - 1} \right) - u\left( {N - 2} \right)} \right)^2}\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + {x^2}\left( {N - 1} \right) + \frac{1}{5}{x^2}\left( {N - 1} \right)
\end{array} \right]}\\
{ \Rightarrow I(x(N - 2),2) = \mathop {\min }\limits_{u(N - 2)} \left[ \begin{array}{l}
{\left( {x\left( {N - 2} \right) - 2u\left( {N - 2} \right)} \right)^2}\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + {\left[ {x\left( {N - 2} \right) - u\left( {N - 2} \right)} \right]^2}\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + \frac{1}{5}{\left[ {x\left( {N - 2} \right) - u\left( {N - 2} \right)} \right]^2}
\end{array} \right]}
\end{array}
\] 由FONC: $\frac{\partial }{{\partial u(N - 2)}} = 0$,可知
\[\begin{array}{l}
\left\{ \begin{array}{l}
 - 4\left( {x\left( {N - 2} \right) - 2u\left( {N - 2} \right)} \right)\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} - 2\left[ {x\left( {N - 2} \right) - u\left( {N - 2} \right)} \right] - \frac{2}{5}\left[ {x\left( {N - 2} \right) - u\left( {N - 2} \right)} \right]
\end{array} \right\} = 0\\
 \Rightarrow {u^*}\left( {N - 2} \right) = \frac{8}{{13}}x\left( {N - 2} \right)
\end{array}
\] 將上式代回 $I(x(N-2),2)$ 可得對應的 Optimal Cost to go 為
\[\begin{array}{l}
I(x(N - 2),2) = \left[ {\frac{9}{{{{13}^2}}}{x^2}\left( {N - 2} \right) + \frac{{25}}{{{{13}^2}}}{x^2}\left( {N - 2} \right) + \frac{5}{{{{13}^2}}}{x^2}\left( {N - 2} \right)} \right]\\
 \Rightarrow I(x(N - 2),2) = \frac{{39}}{{{{13}^2}}}{x^2}\left( {N - 2} \right)
\end{array}
\] $\square$

延伸閱讀:
[動態規劃] 淺談 離散時間動態規劃 (1) - Bellman equation in Infinite Horzion

[Call for Papers] IEEE TCST Special Issue: Advanced Control Technologies for Financial Systems

Submission deadline: February 1, 2027 https://ieeecss.org/publication/transactions-control-systems-technology/special-issues