顯示具有 Linear System 標籤的文章。 顯示所有文章
顯示具有 Linear System 標籤的文章。 顯示所有文章

4/25/2018

[訊號與系統] LTI系統輸入輸出關係由 Convolution 決定 - 從離散時間觀點

以下我們討論 為何 線性非時變 (Linear Time-Invariant, LTI) 系統 輸入與輸出關係 由 所謂的 convolution 表示。為了避免過多繁雜的數學,以下僅討論離散時間的情況。首先我們需要一些定義的幫助:

========================
Definition: Unit Impulse Function in Discrete-Time (Kronecker Delta)
我們說函數 $\delta: \mathbb{N} \to \{0,1\}$ 為 unit impulse function in discrete-time time 若 $\delta$ 滿足
\[\delta \left[ n \right] = \left\{ \begin{gathered}
  1,\;\;\;\;n = 0 \hfill \\
  0,\;\;\;\; o.w. \hfill \\
\end{gathered}  \right.
\]========================

========================
Definition: Impulse Response
給定任意系統配備輸入 $x[n]$ 與輸出 $y[n]$ 關係為 $y[n] = T\{x[n]\}$其中 $T$ 視為 operator (定義在某函數空間),若輸入為 $x[n]=\delta[n]$ 則 輸出
\[
h[n] := y[n] = T\{\delta[n]\}
\]稱為系統 $T$ 的 脈衝響應 (impulse response)
========================


========================
FACT: 任意離散訊號 $x[n]$ 可由 $\delta$ 做組合疊加,亦即
\[
x[n] = \sum_{k=-\infty}^{\infty}x[i] \delta[n-k]
\]========================
Proof: 證明顯然,在此不做贅述。$\square$



========================
Definition: Linear System
給定系統 $T$ 滿足以下輸入與輸出關係: $y_1[n]=T\{x_1[n]\}$ 且 $y_2[n] = T\{x_2[n]\}$。現在定義 $x[n] =ax_1[n] + bx_2[n]$ 其中 $a,b \in \mathbb{R}$則我們說系統 $T$ 為 linear 若下列條件成立
\[
y[n] = T\{x[n]\}
\]且 $y[n] = a y_1[n] + b y_2[n]$
========================

Remarks:
上述 系統 $T$  其實可看按作 泛函分析中的 operator。故線性系統 就是 泛函分析中的 線性算子 (linear operator)。



========================
Definition: Time-Invariant System
給定系統 $T$ 滿足以下輸入與輸出關係: $y[n]=T\{x[n]\}$。現在定義 $x[n] :=x[n-d]$  其中 $d \in \mathbb{N}$ 則我們說 $T$ 為 time-invariant 若下列條件成立
\[
y[n-d] = T\{x[n-d]\}
\]========================


有了以上定義的幫助,我們現在有辦法給出以下為何 LTI 系統的輸入輸出關係確實為 convolution。


========================
Theorem: LTI input-output Relationship is Governed by Convolution
給定 線性非時變 (LTI)系統 $L$ 滿足
\[
y[n] = L\{x[n]\}
\]其中 $L$ 為 linear operator,則
\[
y[n] = \sum_k h[k] x[n-k]
\]========================
Proof:
由於 任意 輸入 $x[n]$ 皆可被表為 脈衝函數的 組合:
\[
x[n] = \sum_{k = -\infty}^\infty x[k] \delta[n-k]
\]現在將此訊號作為LTI系統 $L$ 的輸入,則由於 $L$ 為 linear operator 我們有
\begin{align*}
  y[n] &= L\left\{ {x\left[ n \right]} \right\} \\
&= L\left\{ {\sum\limits_{k =  - \infty }^\infty  x [k]\delta [n - k]} \right\} \hfill \\
   &= \sum\limits_{k =  - \infty }^\infty  x[k] L\left\{ \delta [n - k] \right\}  \;\;\;\; (*)
\end{align*} 接著因為 $L$ 為 time-invariant 我們可進一步改寫
\[L\left\{ {\delta [n - k]} \right\} = h\left[ {n - k} \right] \;\;\;\; (**)
\]故由 $(*)$ 與 $(**)$ 可得
\[y[n] = \sum\limits_{k =  - \infty }^\infty  {x[k]h\left[ {n - k} \right]} \]即為所求。$\square$

Comments:
上述 \[
y[n] = \sum\limits_{k =  - \infty }^\infty  {x[k]h\left[ {n - k} \right]}
\]一般又稱作 $x[n]$ 與 $h[n]$ convolution sum。記作
\[
y[n] = x[n] * h[n]
\]

4/01/2016

[微分方程] 二階常係數 線性 齊次微分方程

首先回憶標準二階微分方程 (2nd ODE)可表為
\[
y'' = f(t,y',y'')
\]現在我們考慮上述方程的一類重要的子集:二階線性常係數 齊次微分方程 (2nd Order Linear Constant Coefficient Homogeneous Ordinary Differential Equation) ,該子集的方程一般可寫作
\[
y'' + py' + q y = 0,\;\;\;\; -\infty < t < \infty \;\;\;\;\; (\star)
\]其中 $p, q \in \mathbb{R}$。讀者可注意到 我們仍處在線性 ODE 的世界,亦即對 $y,y'$ 與 $y''$ 皆線性。我們的目標是要對上述 ODE 就進行求解,亦即要找到某函數 $\phi(t)$ 對 $t \in (-\infty, \infty)$ 都滿足 $(\star)$ 。

Comments: 讀者也許會對於上述 $(\star)$ 為何稱之為 $y'' = f(t,y'=,y')$ 的一類子集感到疑惑,在此我們邀請讀者仔細檢查 $(\star)$,應可發現若令
\[
f(t,y,y') := -(py' + q y )
\] 則確實可看出 $(\star)$ 式 具有 標準二階 ODE 之形式。

以下定理給了我們一個重要的結果來界定 解的形式 與條件:

====================
Theorem: 令 $\phi_1(t)$ 與 $\phi_2(t)$ 在 $t \in (-\infty, \infty)$ 為對下列 二階微分方程的 IVP 問題:
\[
y'' + p y' + qy = 0, \;\;\;\; y(t_0) = y_0;\;\;\;\; y'(t_0) = v_0
\]的兩組解。則
(1) 對任意常數 $c_1,c_2 \in \mathbb{R}$ 其線性組合 $c_1 \phi_1(t) + c_2 \phi_2(t)$ 亦為上述ODE之解。
(2) 除此之外,若下列行列式函數
\[\Delta \left[ {{\phi _1},{\phi _2}} \right]\left( {{t_0}} \right) := \left| {\begin{array}{*{20}{c}}
{{\phi _1}\left( {{t_0}} \right)}&{{\phi _2}\left( {{t_0}} \right)}\\
{{\phi _1}'\left( {{t_0}} \right)}&{{\phi _2}'\left( {{t_0}} \right)}
\end{array}} \right| \ne 0,\begin{array}{*{20}{c}}
{}&{}
\end{array}{t_0} \in \left( { - \infty ,\infty } \right)
\]則對 上述 ODE 的任意解 $\phi(t)$ 滿足以下條件:存在 $c_1, c_2 \in \mathbb{R}$ 使得
$$
\phi(t) = c_1 \phi_1(t) + c_2 \phi_2(t)
$$====================

Comments:
上述行列式條件一般稱為 Wronskian

Proof:
令 $\phi_1(t)$ 與 $\phi_2(t)$ 在 $t \in (-\infty, \infty)$ 為對下列二階 ODE \[ y''+ p y' + qy = 0 \]的解。首先證明則對任意常數 $c_1,c_2 \in \mathbb{R}$ 其線性組合 $c_1 \phi_1(t) + c_2 \phi_2(t)$ 亦為上述ODE之解。故給定任意常數 $c_1,c_2$,令\[ \phi(t) = c_1 \phi_1(t) + c_2 \phi_2(t) \] 則我們僅需證明上式滿足 ODE,現在觀察 \[\begin{array}{l} \phi '' + p\phi ' + q\phi \\ \begin{array}{*{20}{c}} {}&{} \end{array} = {\left( {{c_1}{\phi _1}(t) + {c_2}{\phi _2}(t)} \right)^{\prime \prime }} + p{\left( {{c_1}{\phi _1}(t) + {c_2}{\phi _2}(t)} \right)^\prime } + q\left( {{c_1}{\phi _1}(t) + {c_2}{\phi _2}(t)} \right)\\ \begin{array}{*{20}{c}} {}&{} \end{array} = {c_1}\left( {{\phi _1}^{\prime \prime }(t) + p{\phi _1}^\prime (t) + q{\phi _1}(t)} \right) + {c_2}\left( {{\phi _2}^{\prime \prime }(t) + p{\phi _2}^\prime (t) + q{\phi _2}(t)} \right)\\ \begin{array}{*{20}{c}} {}&{} \end{array}{c_1}\underbrace {\left( {{\phi _1}''(t) + p{\phi _1}'(t) + q{\phi _1}(t)} \right)}_{ = 0} + {c_2}\underbrace {\left( {{\phi _2}''(t) + p{\phi _2}'(t) + q{\phi _2}(t)} \right)}_{ = 0} = 0 \end{array}
\]上述最後一行等式成立因為我們利用了已知假設  $\phi_1(t)$ 與 $\phi_2(t)$ 在 $t \in (-\infty, \infty)$為對下列二階微分方程 $y'' + p y' + qy = 0$ 的解,亦即
\[\left\{ \begin{array}{l} {\phi _1}''(t) + p{\phi _1}'(t) + q{\phi _1}(t) = 0\\ {\phi _2}''(t) + p{\phi _2}'(t) + q{\phi _2}(t) = 0 \end{array} \right.\]

接著我們證 $(2)$,令 $\phi(t)$ 為上述 ODE 的任意解且固定時刻 $t_0 \in (-\infty, \infty)$,則我們可定義
\[
\phi(t_0) := y_0;\;\;\; \phi'(t_0)"=v_0
\]現在我們要證明存在 $c_1, c_2 \in \mathbb{R}$ 使得 $\phi(t) = c_1 \phi_1(t) + c_2 \phi_2(t)$。觀察若 $\phi \left( t \right) = {c_1}{\phi _1}\left( t \right) + {c_2}{\phi _2}\left( t \right)$ 則 $\phi '\left( t \right) = {c_1}{\phi _1}'\left( t \right) + {c_2}{\phi _2}'\left( t \right)$,故此我們可針對 任意時刻 $t_0$寫下
\[\left\{ \begin{array}{l}
\phi \left( {{t_0}} \right) = {c_1}{\phi _1}\left( {{t_0}} \right) + {c_2}{\phi _2}\left( {{t_0}} \right) = {y_0}\\
\phi '\left( {{t_0}} \right) = {c_1}{\phi _1}'\left( {{t_0}} \right) + {c_2}{\phi _2}'\left( {{t_0}} \right) = {v_0 }
\end{array} \right.\]亦即我們得到兩條線性方程且未知數為 $c_1,c_2$ 回憶高中數學,利用 Cramer Rule可知上述線性方程有 唯一解 若且唯若
\[ {c_1} = \frac{{\left| {\begin{array}{*{20}{c}}
{{y_0}}&{{\phi _2}\left( {{t_0}} \right)}\\
{{v_0}}&{{\phi _2}'\left( {{t_0}} \right)}
\end{array}} \right|}}{{\underbrace {\left| {\begin{array}{*{20}{c}}
{{\phi _1}\left( {{t_0}} \right)}&{{\phi _2}\left( {{t_0}} \right)}\\
{{\phi _1}'\left( {{t_0}} \right)}&{{\phi _2}'\left( {{t_0}} \right)}
\end{array}} \right|}_{ = \Delta [{\phi _1},{\phi _2}]\left( {{t_0}} \right)}}};{c_1} = \frac{{\left| {\begin{array}{*{20}{c}}
{{\phi _1}\left( {{t_0}} \right)}&{{y_0}}\\
{{\phi _1}'\left( {{t_0}} \right)}&{{v_0}}
\end{array}} \right|}}{{\underbrace {\left| {\begin{array}{*{20}{c}}
{{\phi _1}\left( {{t_0}} \right)}&{{\phi _2}\left( {{t_0}} \right)}\\
{{\phi _1}'\left( {{t_0}} \right)}&{{\phi _2}'\left( {{t_0}} \right)}
\end{array}} \right|}_{ =  \Delta [{\phi _1},{\phi _2}]\left( {{t_0}} \right)}}}
\]注意到由假設前提可知,對任意 $t_0 \in (-\infty, \infty)$,我們有 $\Delta[\phi_1,\phi_2](t_0) \neq 0$。故我們取上述 $c_1, c_2$ 則任意解 $\phi$ 皆可表成\[
\phi(t) = c_1 \phi_1(t) + c_2 \phi_2(t)
\]且上述等式對任意 $t \in (-\infty, \infty)$ 成立且滿足初始條件。 $\square$


Comments:
1. 微分方程與線性代數彼此之間存在巨大的連結,比如說上述定理中我們先假設有兩解 $\phi_1, \phi_2$ 接著我們才說其他的解亦可透過此兩解做線性組合。所有滿足 ODE 的解所構成之集合稱之為解空間 (Solution Space) 。現在我們提供一個更簡潔(抽象)的觀點:若我們將微分看作算子(operator) ,亦即
\[\begin{array}{l}
y'' + py' + qy = 0\\
 \Rightarrow \left( {\frac{{{d^2}}}{{d{t^2}}} + p\frac{d}{{dt}} + q} \right)\left( y \right) = 0\\
 \Rightarrow L\left( y \right) = 0
\end{array}
\]則讀者可驗證微分算子 $L(\cdot)$ 為 線性算子,且其解空間
\[
\{\phi: L (\phi) = 0\}
\]空間等價為 線代中的零空間(Null Space)

2. 我們稱 $\phi_1, \phi_2$ 為在 $t \in (-\infty, \infty)$ 的一組解空間的 基底 ( basis) 若且唯若 $\Delta[\phi_1, \phi_2](t_0) \neq 0$ 對任意 $t_0 \in (-\infty, \infty)$,那麼該如何找到此兩解?我們要求此兩解滿足 Wronskian 行列式條件
\[
{\Delta [{\phi _1},{\phi _2}]\left( {{t_0}} \right)} \neq 0
\]事實上此條件等價說明 $\phi_1, \phi_2$ 彼此線性獨立。

3. 由於此系統為二階系統,我們僅需找到兩個彼此線性獨立的解 $\phi_1, \phi_2$,亦即 對任意 $p,q$ 且任意 $t_0 \in (-\infty, \infty)$,下列 Wronkian
\[\Delta \left( {{t_0}} \right): = \left| {\begin{array}{*{20}{c}}
{{\phi _1}\left( {{t_0}} \right)}&{{\phi _2}\left( {{t_0}} \right)}\\
{{\phi _1}'\left( {{t_0}} \right)}&{{\phi _2}'\left( {{t_0}} \right)}
\end{array}} \right| \ne 0
\]

4. 在前述證明中,我們在最後需回頭驗證上述結果是否滿足初始條件 $y(t_0) = y_0, \; y'(t_0) = v_0$,且且我們稱 $\phi(t)$
$$
\phi(t) := c_1 \phi_1(t) + c_2 \phi_2(t),\;\;\; -\infty< t < \infty, \;\;c_1,c_2\in \mathbb{R}
$$ 為 general solution。


儘管我們有了前述定理的幫助,知道對任意二階ODE 我們要先找出兩組線性獨立解,則其餘所有可能的解可透過此兩組線性獨立解基底做線性組合而得,但我們還有一個問題尚待解決:也就是該如何事先找到此兩組線性獨立解 $\phi_1, \phi_2$?一個常用的辦法為 先猜出解的形式,再來決定 待定係數,我們將此法稱作 Guess Method:


透過 Guess Method 來找出 $\phi_1(t), \phi_2(t)$:
由於我們要求解的方程為 ODE ,故我們可猜測 $e^{zt}, \;\;\; t \in (-\infty, \infty)$ 為解候選 (solution candidate),其中 $z$ 為待定常數 (注意到在此 $z$ 不一定要是實數,可以是複數)。則
\[\begin{array}{l}
y'' + py' + qy = 0\\
 \Rightarrow {z^2}\left( {{e^{zt}}} \right) + pz\left( {{e^{zt}}} \right) + q\left( {{e^{zt}}} \right) = 0\\
 \Rightarrow \left( {{z^2} + pz + q} \right){e^{zt}} = 0
\end{array}\]因此 $e^{zt}$ 為 $y''+p y' + qy=0$ 的解 若且唯若
\[
(z^2 + p z + q) e^{zt} = 0
\]因為 $e^{zt} \neq 0$ ,故充分必要條件可簡化成
\[
z^2 + p z + q =0
\]上述 等式 稱為 特性方程式 (Characteristic equation) 或者 輔助方程 (Auxiliary equation),由於此特性方程為二次函數,回憶中學數學裡,我們可對其求根 $z$ 如下
\[
{z_{1,2}} = \frac{{ - p \pm \sqrt {{p^2} - 4q} }}{2}
\]注意到上述 $z$ 現在我們討論以下三種可能

Case 1: $p^2 > 4 q$,則 $z_{1,2}$ 為 實數 且 相異

Case 2: $p^2 < 4q$,則 $z_{1,2}$ 為一對共軛複數

Case 3: $p^2 = 4q$ ,則 $z_1 = z_2$ (重根)

以下我們分別針對此三種情況來求出 我們對二街常係數齊次 ODE 所需要的線性獨立解 $\phi_1$ 與 $\phi_2$:

對於 Case 1,我們有 $\phi_1(t) = e^{z_1 t}$ 且 $\phi_2(t) = e^{z_2 t}$ 且 $z_1 \neq z_2$,故我們檢驗 Wronskian
\[\begin{array}{l}
\Delta \left( {{t_0}} \right) = \left| {\begin{array}{*{20}{c}}
{{\phi _1}\left( {{t_0}} \right)}&{{\phi _2}\left( {{t_0}} \right)}\\
{{\phi _1}'\left( {{t_0}} \right)}&{{\phi _2}'\left( {{t_0}} \right)}
\end{array}} \right|\\
 \Rightarrow \Delta \left( {{t_0}} \right) = \left| {\begin{array}{*{20}{c}}
{{e^{{z_1}{t_0}}}}&{{e^{{z_2}{t_0}}}}\\
{{z_1}{e^{{z_1}{t_0}}}}&{{z_2}{e^{{z_2}{t_0}}}}
\end{array}} \right|\\
 \Rightarrow \Delta \left( {{t_0}} \right) = \left( {{z_2} - {z_1}} \right){e^{{z_1}{t_0}}}{e^{{z_2}{t_0}}}
\end{array}
\] 由於 $z_1 \neq z_2$ 且 $e^{(z_1 + z_2)t_0}$ 永遠不為零,則 $\Delta(t_0) \neq 0$ 對任意 $t_0 \in (-\infty, \infty)$ ,故對於 Case 1,建構線性組合
\[
\phi(t) := c_1 \phi_1(t) + c_2 \phi_2(t), \;\;\;\; t\in (-\infty, \infty)
\]且此即為我們的 general solution。

同理,對於 Case 2 ,此時
\[{z_{1,2}} = \frac{{ - p \pm i\sqrt {4q - {p^2}} }}{2}\]故我們檢驗 Wronskian
\[\begin{array}{*{20}{l}}
{\Delta \left( {{t_0}} \right) = \left| {\begin{array}{*{20}{c}}
{{\phi _1}\left( {{t_0}} \right)}&{{\phi _2}\left( {{t_0}} \right)}\\
{{\phi _1}^\prime \left( {{t_0}} \right)}&{{\phi _2}^\prime \left( {{t_0}} \right)}
\end{array}} \right|}\\
{ \Rightarrow \Delta \left( {{t_0}} \right) = \left| {\begin{array}{*{20}{c}}
{{e^{{z_1}{t_0}}}}&{{e^{{z_2}{t_0}}}}\\
{{z_1}{e^{{z_1}{t_0}}}}&{{z_2}{e^{{z_2}{t_0}}}}
\end{array}} \right|}\\
\begin{array}{l}
\begin{array}{*{20}{c}}
{}&{}
\end{array} = \left( {{z_2} - {z_1}} \right){e^{\left( {{z_1} + {z_2}} \right){t_0}}}\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} = \left( {\frac{{ - p + i\sqrt {4q - {p^2}} }}{2} - \frac{{ - p - i\sqrt {4q - {p^2}} }}{2}} \right){e^{\left( {\frac{{ - p + i\sqrt {4q - {p^2}} }}{2} + \frac{{ - p - i\sqrt {4q - {p^2}} }}{2}} \right){t_0}}}\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} = \left( {i\sqrt {4q - {p^2}} } \right){e^{\left( { - 2p} \right){t_0}}} \ne 0
\end{array}
\end{array}\]故對於 Case 2,我們亦可建構線性組合
\[
\begin{array}{l}
\phi (t): = {c_1}{\phi _1}(t) + {c_2}{\phi _2}(t)\\
 \Rightarrow \phi (t) = {c_1}{e^{{z_1}{t_0}}} + {c_2}{e^{{z_2}{t_0}}}
\end{array}
\]其中 \[{z_{1,2}} = \frac{{ - p \pm i\sqrt {4q - {p^2}} }}{2}
\]此即為我們的 general solution。

對於 Case 3 ,我們陷入重根的情況,此時 $z_1 = z_2 = -p/2$,也就是說我們僅有一解
\[
\phi_1(t)  = e^{z_1t} = e^{-(p/2)t}
\] 但我們知道對於二階 ODE 要構成完整的解空間,需要兩組線性獨立解,故我們需要再找出一組解與其 $\phi_1(t)$線性獨立。故我們猜測有一待定函數 $w(t)$ 使得
\[
\phi_2(t) := e^{z_2 t} w(t) = e^{}
\]亦為對 二皆 ODE 之一解,則 $\phi_2$ 滿足
\[
\phi_2'' + p \phi_2' + q \phi_2 = 0 ,\;\;\; t \in (-\infty, \infty)
\]此表示
\[
\begin{array}{l}
{\phi _2}^{\prime \prime } + p{\phi _2}^\prime  + q{\phi _2} = 0\\
 \Rightarrow {\left( {{e^{{z_2}t}}w(t)} \right)^{\prime \prime }} + p{\left( {{e^{{z_2}t}}w(t)} \right)^\prime } + q\left( {{e^{{z_2}t}}w(t)} \right) = 0\\
 \Rightarrow \left( {{z_2}^{}w'(t) + w''\left( t \right)} \right){e^{{z_2}t}} + {z_2}{e^{{z_2}t}}\left( {{z_2}^{}w(t) + w'\left( t \right)} \right)\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + p\left( {{z_2}{e^{{z_2}t}}w(t) + w'\left( t \right){e^{{z_2}t}}} \right) + q\left( {{e^{{z_2}t}}w(t)} \right) = 0\\
 \Rightarrow \left[ {w''\left( t \right) + \left( {q - \frac{{{p^2}}}{4}} \right)w(t)} \right]{e^{ - \frac{p}{2}t}} = 0
\end{array}\]由於在 Case 3 我們知道 $p^2 = 4q$ 故上述條件退化成
\[\begin{array}{l}
\left[ {w''\left( t \right) + \underbrace {\left( {q - \frac{{{p^2}}}{4}} \right)}_{ = 0}w(t)} \right]{e^{ - \frac{p}{2}t}} = 0\\
 \Rightarrow w''\left( t \right){e^{ - \frac{p}{2}t}} = 0
\end{array}
\]但注意到 $e^{-pt/2} \neq 0$ 故我們要求 $w''(t) =0$對任意 $t \in (-\infty,\infty)$ 故我們僅需選擇 $w(t)$ 使得 $w''(t)=0$ 成立。那麼不難發現可令 $w(t)$ 為多項式形式,故我們猜測以下兩組 $w(t)$ 作為 候選:
\[
w(t) = 1 \;\;\; w(t) = t
\]但注意到若 $w(t) =1$ 則我們的解 $\phi_2(t) = e^{z_2 t} w(t) = e^{z_2 t} =\phi_1(t)$ 不符合線性獨立要求,故我們選 $w(t) = t$ 則
$$
\phi_2(t) := e^{z_2 t} w(t) =  e^{z_2 t} t = t \phi_1(t)
$$ 現在我們檢驗 Wronskian:
\[\begin{array}{l}
\Delta [{\phi _1},{\phi _2}]\left( {{t_0}} \right) = \left| {\begin{array}{*{20}{c}}
{{\phi _1}\left( {{t_0}} \right)}&{{\phi _2}\left( {{t_0}} \right)}\\
{{\phi _1}^\prime \left( {{t_0}} \right)}&{{\phi _2}^\prime \left( {{t_0}} \right)}
\end{array}} \right|\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = \left| {\begin{array}{*{20}{c}}
{{e^{ - pt/2}}}&{t{e^{ - pt/2}}}\\
{\left( {\frac{{ - p}}{2}} \right){e^{ - pt/2}}}&{{e^{ - pt/2}} + \left( {\frac{{ - p}}{2}} \right)t{e^{ - pt/2}}}
\end{array}} \right|\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = {e^{ - pt/2}}\left( {{e^{ - pt/2}} + t\left( {\frac{{ - p}}{2}} \right){e^{ - pt/2}}} \right) - t{e^{ - pt/2}}\left( {\frac{{ - p}}{2}} \right){e^{ - pt/2}}\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = {e^{ - pt}} \ne 0
\end{array}
\]故可知\[
\phi_2(t) = t e^{-pt/2}
\]確實為另一解,且此解與先前的 $\phi_1$ 彼此線性獨立。

3/28/2016

[控制理論] 具有負實部特徵值之 LTV系統 並不保證系統穩定

考慮線性非時變 (Linear Time-Invariant, LTI)系統利用狀態空間表示:
\[
{\bf \dot x} = A {\bf x} + B {\bf u}
\] 回憶在大學部自動控制課程中,我們知道 LTI 系統穩定 的 充分必要條件 為系統矩陣 $A$ 之特徵值具有負實部 (或者等價論述為 極點 pole 落在 複數平面的左半面)。現在我們想問若 系統為 線性時變 (Linear Time-Varying, LTV)系統是否此條件依然成立?

答案是否定的,以下為一個極為出色的反例:考慮線性時變系統 ${\bf \dot x} = A(t) {\bf x} $ 其中
\[A\left( t \right): = \left[ {\begin{array}{*{20}{c}}
{ - 1}&{{e^{2t}}}\\
0&{ - 1}
\end{array}} \right]
\] 且給定初始狀態為 $x_1(0) = x_2(0)=1$ 則由於此系統 $A(t)$ 矩陣為三角矩陣,其特徵值為對角線元素,亦即 $\lambda_{1,2} = -1$,具有負實部。然而,若我們求解此 LTV 系統,亦即觀察
\[{\bf{\dot x}}\left( t \right) = \left[ {\begin{array}{*{20}{c}}
{ - 1}&{{e^{2t}}}\\
0&{ - 1}
\end{array}} \right]\left[ \begin{array}{l}
{x_1}\left( t \right)\\
{x_2}\left( t \right)
\end{array} \right] = \left[ \begin{array}{l}
 - {x_1}\left( t \right) + {e^{2t}}{x_2}\left( t \right)\\
 - {x_2}\left( t \right)
\end{array} \right]\]故我們可首先解得
\[\begin{array}{*{20}{l}}
{{{\dot x}_2}\left( t \right) =  - {x_2}\left( t \right)}\\
\begin{array}{l}
 \Rightarrow {x_2}\left( t \right) = {e^{ - t}}{x_2}\left( 0 \right)\\
 \Rightarrow {x_2}\left( t \right) = {e^{ - t}}
\end{array}
\end{array}
\]再將此 $x_2(t)$ 帶回 $\dot x_1(t)$ 式中,可求解 $x_1$ 如下
\[\begin{array}{*{20}{l}}
{{{\dot x}_1}\left( t \right) =  - {x_1}\left( t \right) + {e^{2t}}{x_2}\left( t \right)}\\
{ \Rightarrow {{\dot x}_1}\left( t \right) =  - {x_1}\left( t \right) + {e^{2t}}{e^{ - t}}}\\
{ \Rightarrow {x_1}\left( t \right) = {e^{ - t}}{x_1}\left( 0 \right) + \int_0^t {{e^{ - \left( {t - \tau } \right)}}{e^\tau }d\tau } }\\
{ \Rightarrow {x_1}\left( t \right) = {e^{ - t}} + {e^{ - \left( t \right)}}\int_0^t {{e^{2\tau }}d\tau } }\\
{ \Rightarrow {x_1}\left( t \right) = \frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}}
\end{array}\]故系統之解為
\[{{\bf{x}}\left( t \right) = \left[ \begin{array}{l}
{e^{ - t}}\\
\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}
\end{array} \right]}\]但注意到若我們計算上述之狀態的 2-norm 且取極限 $t \to \infty$ 會發現
\[\begin{array}{l}
\mathop {\lim }\limits_{t \to \infty } \left\| {{\bf{x}}\left( t \right)} \right\| = \mathop {\lim }\limits_{t \to \infty } \left\| {\left[ \begin{array}{l}
{e^{ - t}}\\
\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}
\end{array} \right]} \right\|\\
 = \mathop {\lim }\limits_{t \to \infty } {\left( {\left[ {\begin{array}{*{20}{c}}
{{e^{ - t}}}&{\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}}
\end{array}} \right]\left[ \begin{array}{l}
{e^{ - t}}\\
\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}
\end{array} \right]} \right)^{1/2}}\\
 = \mathop {\lim }\limits_{t \to \infty } {\left( {{e^{ - 2t}} + {{\left( {\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}} \right)}^2}} \right)^{1/2}} = \infty
\end{array}\]亦即系統狀態發散。

上述結果闡釋了對於 LTV 系統而言,負實部特徵值 (左半面極點) 不保證系統穩定。

2/21/2015

[最佳估計] 狀態估測器的收斂性

考慮離散時間 LTI 系統
\[\begin{array}{l}
x(k + 1) = Ax(k) + w(k)\\
y(k) = Cx\left( k \right) + v(k)
\end{array}\]其中 $x \in \mathbb{R}^n, u \in \mathbb{R}^m, y \in \mathbb{R}^p$ 且 $x(0) = x_0$;

基本狀態估計問題:給定初始估計誤差 ( i.e., $\bar x(0) \neq x_0$) 且 不考慮 雜訊 無外部干擾的情況 ($w(k) = v(k)=0$),我們想問 $\hat x(k) \to x(k)$ as $k \to \infty$ ?

=================
Theorem: (Convergence of Estimator Cost)
給定無 noise 量測輸出 ${\bf y}( T) = \{Cx(0), CAx(0),..., CA^T x(0)\}$ 則最佳估測器的 cost $V_T^*({\bf y}(T)) $ 在 $T \to \infty$ 時收斂 。
=================

Proof:
由於
\[V_T^{} = \frac{1}{2}\left( \begin{array}{l}
\left| {\hat x\left( 0 \right) - \bar x\left( 0 \right)} \right|_{{{\left( {{P^ - }\left( 0 \right)} \right)}^{ - 1}}}^2\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + \sum\limits_{k = 0}^{T - 1} {|\hat x\left( {k + 1} \right) - A\hat x\left( k \right)|_{{Q^{ - 1}}}^2}  + \sum\limits_{k = 0}^T {\left| {y\left( k \right) - C\hat x\left( k \right)} \right|_{{R^{ - 1}}}^2}
\end{array} \right)\]故我們分三個步驟證明 $V_T^*({\bf y}(T)) $ 在 $T \to \infty$ 時收斂 :
首先證明 sequence $\{V_T^*\}$ 有上界(bounded above);接著我們證明 $\{V_T^*\}$ 為 nondecreasing。則由步驟一與步驟二可推論 $\{V_T^*\}$ 必定收斂。

現在我們開始證明 $V_T^*({\bf y}(T)) $ 有界:
由於我們的目標為 $\min_{\hat{ {\bf x}}(k)} V_k$ 但我們並知道 怎樣的 $\hat x$ 可以幫助我們達成此目標,故我們先暫取 $\hat x(0) :=  x_0$ (並不一定為最佳解!) 則我們有
\[\begin{array}{l}
\left\{ \begin{array}{l}
\hat x(1) = A\hat x(0)\\
\hat x(2) = A\hat x(1) = {A^2}\hat x(0)\\
 \vdots
\end{array} \right.\\
 \Rightarrow y\left( k \right) = C{A^k}{x_0}
\end{array}\]故若將此結果帶入我們的 cost 可得
\[\bar{V}_T^{} = \frac{1}{2}\left| {\hat x\left( 0 \right) - \bar x\left( 0 \right)} \right|_{{{\left( {{P^ - }\left( 0 \right)} \right)}^{ - 1}}}^2 < \infty
\]但注意到此並非最佳 cost,若代入最佳解 則我們必有
\[
V_T^* \le \bar V_T
\]故可推論 sequence $\{V_T^*\}$ 必定 bounded above。

接著我們證明 optimal cost sequence $\{V_T^*\}$ 為 nondecreasing。
給定量測輸出 ${\bf y}( T) = \{Cx(0), CAx(0),..., CA^T x(0)\}$,定義 在時間 $T$ 的最佳狀態 sequence 為
\[\left\{ {\hat x\left( {0} \right),\hat x\left( {1} \right),...,\hat x\left( {T} \right)} \right\}
\]現在我們比對 $T$ 時刻的 optimal cost 與 $T-1$ 時刻的 optimal cost 可知
\[V_T^* - \frac{1}{2}\left( {\left| {y\left( {T - 1} \right) - C\hat x\left( {T - 1} \right)} \right|_{{R^{ - 1}}}^2 + |\hat x\left( T \right) - A\hat x\left( k \right)|_{{R^{ - 1}}}^2} \right) \ge V_{T - 1}^*\]此說明了 sequence $\{V_T^*\}$ 為 nondecreasing。

故 由  optimal cost sequence $\{V_T^*\}$  bounded above 與 $\{V_T^*\}$ 為 nondecreasing,我們可推論當 $T \to \infty$  $\{V_T^*\}$ 必定收斂。 $\square$


注意到 optimal estimator cost $V_T^*$ 的收斂性與 系統可觀測性無關,但若我們要求我們的估計狀態 $\hat x \to x$ 則系統的觀測性將扮演重要腳色,我們將此結果記做以下定理:

==========================
Theorem: Estimator Convergence
考慮控制系統 $(A,C)$ 可觀測 且 $Q,R >0$ 為正定矩陣 且 給定一組無雜訊量測輸出
\[
{\bf y}(T) := \{Cx(0), CAx(0), ..., CA^T x(0)\}
\]則 最佳狀態估測 收斂到原本系統狀態;亦即
\[
\hat x(T) \to x(T) \;\; \text{ as $T \to \infty$}
\]==========================

Proof:
利用 在時刻 $T+n-1$ 的最佳解作為在時刻 $T-1$的 decision variables ,則前述 Theorem 告訴我們可寫
\[\small V_{T + n - 1}^* - \frac{1}{2}\left[ {\sum\limits_{k = T - 1}^{T + n - 2} {|\hat x\left( {k + 1} \right) - A\hat x\left( k \right)|_{{Q^{ - 1}}}^2}  + \sum\limits_{k = T}^{T + n - 1} {\left| {y\left( k \right) - C\hat x\left( k \right)} \right|_{{R^{ - 1}}}^2} } \right] \ge V_{T - 1}^*\]做變數變換 令 $j = k - T$ 我們可得
\[ \small V_{T + n - 1}^* - \frac{1}{2}\left[ {\sum\limits_{j =  - 1}^{n - 2} {|\hat x\left( {T + j + 1} \right) - A\hat x\left( {T + j} \right)|_{{Q^{ - 1}}}^2}  + \sum\limits_{j = 0}^{n - 1} {\left| {y\left( {j + T} \right) - C\hat x\left( {j + T} \right)} \right|_{{R^{ - 1}}}^2} } \right] \ge V_{T - 1}^*\]注意到當 $T \to \infty$時, $\{V_T^*\}$ 收斂 且 $Q^{-1}, R^{-1} >0$ 故上式
\[\begin{array}{l}
\underbrace {V_{T + n - 1}^* - V_{T - 1}^*}_{ \to 0} \ge \frac{1}{2}\left[ \begin{array}{l}
\sum\limits_{j =  - 1}^{n - 2} {|\hat x\left( {T + j + 1} \right) - A\hat x\left( {T + j} \right)|_{{Q^{ - 1}}}^2} \\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + \sum\limits_{j = 0}^{n - 1} {\left| {y\left( {j + T} \right) - C\hat x\left( {j + T} \right)} \right|_{{R^{ - 1}}}^2}
\end{array} \right]\\
 \Rightarrow \sum\limits_{j =  - 1}^{n - 2} {|\hat x\left( {T + j + 1} \right) - A\hat x\left( {T + j} \right)|_{{Q^{ - 1}}}^2}  + \sum\limits_{j = 0}^{n - 1} {\left| {y\left( {j + T} \right) - C\hat x\left( {j + T} \right)} \right|_{{R^{ - 1}}}^2}  \to 0\\
 \Rightarrow \left\{ \begin{array}{l}
\hat x\left( {T + j + 1} \right) - A\hat x\left( {T + j} \right) \to 0,\begin{array}{*{20}{c}}
{}
\end{array}\forall j =  - 1,...,n - 2\\
y\left( {j + T} \right) - C\hat x\left( {j + T} \right) \to 0,\begin{array}{*{20}{c}}
{}
\end{array}\forall j = 0,...,n - 1
\end{array} \right. \ \ \ \ \ \ (**)
\end{array}\]令 $\hat w_T(j) := \hat x(T+j+1|T+n-1) - A \hat x(T+j|T+n-1)$ 並且透過系統方程 $x(k+1) = Ax(k) + w(k)$ 我們有
\[\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{\hat x\left( {T|T + n - {\rm{1}}} \right)}\\
{\hat x\left( {T{\rm{ + 1}}|T + n - {\rm{1}}} \right)}\\
 \vdots \\
{\hat x\left( {T + n - 1|T + n - {\rm{1}}} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
I\\
A\\
 \vdots \\
{{A^{n - 1}}}
\end{array}} \right]\hat x\left( {T|T + n - {\rm{1}}} \right)\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} + \left[ {\begin{array}{*{20}{c}}
0&{}&{}&{}\\
I&0&{}&{}\\
 \vdots & \vdots & \ddots &{}\\
{{A^{n - 2}}}&{{A^{n - 3}}}& \cdots &I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{{\hat w}_T}\left( 0 \right)}\\
{{{\hat w}_T}\left( 1 \right)}\\
 \vdots \\
{{{\hat w}_T}\left( {n - 2} \right)}
\end{array}} \right]\ \ \ \ (*)
\end{array} \]且由於我們的量測輸出 滿足
\[\left\{ \begin{array}{l}
y\left( T \right) = Cx\left( T \right)\\
y\left( {T + 1} \right) = CAx\left( T \right)\\
 \vdots \\
y\left( {T + n - 1} \right) = C{A^{n - 1}}x\left( T \right)
\end{array} \right. \Rightarrow \left[ {\begin{array}{*{20}{c}}
{y\left( T \right)}\\
{y\left( {T + 1} \right)}\\
 \vdots \\
{y\left( {T + n - 1} \right)}
\end{array}} \right] = Ox\left( T \right)\]其中 $O$ 為 observability matrix。現在用上式減去 同乘 $C$ 矩陣 後的 $(*)$ 可得
\[\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{y\left( T \right)}\\
{y\left( {T + 1} \right)}\\
 \vdots \\
{y\left( {T + n - 1} \right)}
\end{array}} \right] - C\left[ {\begin{array}{*{20}{c}}
{\hat x\left( {T|T + n - {\rm{1}}} \right)}\\
{\hat x\left( {T{\rm{ + 1}}|T + n - {\rm{1}}} \right)}\\
 \vdots \\
{\hat x\left( {T + n - 1|T + n - {\rm{1}}} \right)}
\end{array}} \right] = O\left[ {x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right)} \right]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} + \left[ {\begin{array}{*{20}{c}}
0&{}&{}&{}\\
C&0&{}&{}\\
 \vdots & \vdots & \ddots &{}\\
{C{A^{n - 2}}}&{C{A^{n - 3}}}& \cdots &C
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{{\hat w}_T}\left( 0 \right)}\\
{{{\hat w}_T}\left( 1 \right)}\\
 \vdots \\
{{{\hat w}_T}\left( {n - 2} \right)}
\end{array}} \right]
\end{array}\]現在用 $(**)$ 可知
\[\begin{array}{l}
\underbrace {\left[ {\begin{array}{*{20}{c}}
{y\left( T \right)}\\
{y\left( {T + 1} \right)}\\
 \vdots \\
{y\left( {T + n - 1} \right)}
\end{array}} \right] - C\left[ {\begin{array}{*{20}{c}}
{\hat x\left( {T|T + n - {\rm{1}}} \right)}\\
{\hat x\left( {T{\rm{ + 1}}|T + n - {\rm{1}}} \right)}\\
 \vdots \\
{\hat x\left( {T + n - 1|T + n - {\rm{1}}} \right)}
\end{array}} \right]}_{ \to 0} = O\left[ {x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right)} \right]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} + \underbrace {\left[ {\begin{array}{*{20}{c}}
0&{}&{}&{}\\
C&0&{}&{}\\
 \vdots & \vdots & \ddots &{}\\
{C{A^{n - 2}}}&{C{A^{n - 3}}}& \cdots &C
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{{\hat w}_T}\left( 0 \right)}\\
{{{\hat w}_T}\left( 1 \right)}\\
 \vdots \\
{{{\hat w}_T}\left( {n - 2} \right)}
\end{array}} \right]}_{ \to 0}\\
 \Rightarrow O\left[ {x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right)} \right] \to 0
\end{array}\]又因為 observability matrix $O$ 有 linear independent columns,故我們可推知
\[x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right) \to 0 \;\; \text{as $T \to \infty$}\]亦即
\[\hat x\left( {T|T + n - {\rm{1}}} \right) \to x\left( T \right)
\]再者如果我們觀察 $(*)$ 可以發現因為 $\hat w_T (j) \to 0$ 當 $T \to \infty$,故
\[\hat x\left( {T + n - 1|T + n - {\rm{1}}} \right) \to {A^{n - 1}}\hat x\left( {T|T + n - {\rm{1}}} \right) \;\;\;\; \text{ as $T \to \infty$}\]又因為 $A^{n-1} x(T) = x(T + n -1)$ ,我們有
\[\begin{array}{l}
x\left( {T + n - 1} \right) - \hat x\left( {T + n - 1|T + n - {\rm{1}}} \right) \to {A^{n - 1}}x\left( T \right) - {A^{n - 1}}\hat x\left( {T|T + n - {\rm{1}}} \right)\\
 \Rightarrow x\left( {T + n - 1} \right) - \hat x\left( {T + n - 1|T + n - {\rm{1}}} \right) \to {A^{n - 1}}\underbrace {\left[ {x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right)} \right]}_{ \to 0}\\
 \Rightarrow x\left( {T + n - 1} \right) \to \hat x\left( {T + n - 1|T + n - {\rm{1}}} \right)\\
 \Rightarrow x\left( {T + n - 1} \right) \to \hat x\left( {T + n - 1} \right)\\
 \Rightarrow x\left( j \right) \to \hat x\left( j \right),\begin{array}{*{20}{c}}
{}
\end{array}\forall j \to \infty
\end{array}\]至此得證

2/19/2015

[控制理論] 離散線性系統的 追蹤 與 調節 問題

Setpoint Tracking
考慮離散線性系統
\[\left\{ \begin{array}{l}
x\left( {k + 1} \right) = Ax\left( k \right) + Bu\left( k \right)\\
y\left( k \right) = Cx\left( k \right)
\end{array} \right.\] $x \in \mathbb{R}^n, y \in \mathbb{R}^p, u \in \mathbb{R}^m$;一般而言,控制系統中常見的 追蹤(setpoint tracking)問題 (或稱 servo problem) 如下:
給定 setpoint  $y_{sp}$ (e.g., 步階訊號 or 常數),我們希望在系統達到穩態(steady-state) 之後,系統的輸出 $ y  = y_{sp} $。

若我們考慮給定的 setpoint $y_{sp} = 0$ 則稱此類問題為 regulation problem。

回憶在最佳控制理論中的 LQR 方法給予我們對於線性系統的 regulation 問題提供一組最佳解,故我們的想問是否能將此法應用在 一般的追蹤問題?

答案是肯定的,僅需引入 deviation variable 做基本 座標轉換 即可。

現在我們令 $y_{sp}$ 為 output setpoint,定義系統穩態時候的狀態 與 控制力為 $(x_s, u_s)$ 則線性系統方程可改寫為
\[\left\{ \begin{array}{l}
{x_s} = A{x_s} + B{u_s}\\
{y_s} = C{x_s}
\end{array} \right.\]注意到對於穩態時 我們希望 $y_s = y_{sp}$ 故 我們有 $y_{sp} = C x_s$ 現在我們將上式改寫為矩陣形式
\[\left\{ \begin{array}{l}
{x_s} = A{x_s} + B{u_s}\\
{y_{sp}} = C{x_s}
\end{array} \right. \Rightarrow \left[ {\begin{array}{*{20}{c}}
{I - A}&{ - B}\\
C&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_s}}\\
{{u_s}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
0\\
{{y_{sp}}}
\end{array}} \right] \ \ \ \ (*)
\]若上述矩陣方程有解 (亦即可以解出 $x_s, u_s$),且系統控制力無拘束,則我們可以定義 deviation variables 如下
\[\left\{ \begin{array}{l}
\tilde x\left( k \right): = x\left( k \right) - {x_s}\\
\tilde u\left( k \right): = u\left( k \right) - {u_s}
\end{array} \right.\]現在觀察
\[\begin{array}{l}
\tilde x\left( {k + 1} \right) = x\left( {k + 1} \right) - {x_s}\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = Ax\left( k \right) + Bu\left( k \right) - \left( {A{x_s} + B{u_s}} \right)\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = A\left( {x\left( k \right) - {x_s}} \right) + B\left( {u\left( k \right) - {u_s}} \right)\\
 \Rightarrow \tilde x\left( {k + 1} \right) = A\tilde x\left( k \right) + B\tilde u\left( k \right)
\end{array}\]此表示我們的 deviation variable 仍然滿足原本給定的線性系統。且此時若考慮使用 deviation variable $\tilde{x}, \tilde{u}$ 則 我們的 setpoint tracking problem 被改寫為 regulation problem ;也就是說我們要找到一組 $\tilde{u}(k)$ 使得 $\tilde x(k) \to 0$ (此等價為 $x(k) \to x_s$,故 $C x(k) \to C x_s = y_{sp}$)。解完此 regulation problem 之後,真實的控制力 $u(k) = \tilde u(k) + u_s$。

Comment:
注意到上述論述建立在 $(n + p) \times (n + m)$ 的矩陣方程
\[\left\{ \begin{array}{l}
{x_s} = A{x_s} + B{u_s}\\
{y_{sp}} = C{x_s}
\end{array} \right. \Rightarrow \left[ {\begin{array}{*{20}{c}}
{I - A}&{ - B}\\
C&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_s}}\\
{{u_s}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
0\\
{{y_{sp}}}
\end{array}} \right] \ \ \ \ (*)
\]有解。但何時才有解/與此解是否唯一的問題 我們並未解決 !! 以下我們將對此點進行討論。


對任意 setpoint $y_{sp}$, 上述矩陣方程 $(*)$ 有解 的充分條件為:$(*)$ 要有 linearly independent rows;亦即 $p \le m$ 也就是說我們需要 控制力的數目 $m$ 與 量測輸出的數目 $p$ 至少要相等)。但在實際情況上卻是常常相反,我們可能會獲得非常多的量測輸出數目 (因為裝了很多 sensor),但實際可以調控的變數卻少於 sensor 數目。為了要解決此問題,我們可選定某矩陣 $H$ 並引入新的 控制變數 $r \in \mathbb{R}^{n_c}$ 為 量測輸出的線性組合;亦即
\[
r := Hy
\]如此一來,若 $p > m$ 情況發生時,我們可選一部份的 輸出 $n_c \le m$ 作為控制變數 並且 對此引入的控制變數 $r$ 也給定所需的 setpoint $r_{sp}$。

另外若 $m > p$ (控制力數目 大於 量測輸出的數目) 時,則對某些 $H$ 與 $r_{sp}$,矩陣方程 $(*)$ 有解 但此時唯一性並不被保證,故我們會希望有唯一解,此時需要對 穩態控制力  $u_s$ 也給定 setpoint $u_{sp}$。


總和以上所述,我們可以建構以下 Steady-State Target Problem
========================
Steady-State Target Problem
考慮最佳化問題
\[
\min_{x_s, u_s} \frac{1}{2} (|u_s - u_{sp}|_{R_s}^2 + |y_s - y_{sp}|_{Q_s}^2)
\]subject to
\[\left\{ \begin{array}{l}
{x_s} = A{x_s} + B{u_s}\\
{y_{sp}} = C{x_s}
\end{array} \right. \Rightarrow \left[ {\begin{array}{*{20}{c}}
{I - A}&{ - B}\\
C&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_s}}\\
{{u_s}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
0\\
{{y_{sp}}}
\end{array}} \right]
\]與 $E u_s \le e$ 與 $F C x_s \le f$。其中 $R_s$ 為 正定矩陣。且對 控制變數的 setpoints $r_{sp}$,上述 target problem 為 feasible 。
========================

有了 steady-state 的解之後我們可以便可以求解 當初我們引入 deviation variable 的 regulation problem:
========================
Dynamic Regulation Problem
考慮以下 cost function
\[\begin{array}{l}
V\left( {\tilde x\left( 0 \right),{\bf{\tilde u}}} \right) = \frac{1}{2}\sum\limits_{k = 0}^{N - 1} {\left| {\tilde x\left( k \right)} \right|_Q^2 + \left| {\tilde u\left( k \right)} \right|_R^2} \\
s.t.\\
\tilde x\left( {k + 1} \right) = A\tilde x\left( k \right) + B\tilde u\left( k \right)
\end{array}\]其中 $\tilde x(0) = \hat x(k) - x_s$ (此初始條件表示  透過 steady-state $x_s$ 平移估計狀態 $\hat x$ 而得。);且所求得的 regulator 將會求解以下的 regulation problem
 \[\begin{array}{l}
\mathop {\min }\limits_{{\bf{\tilde u}}} V\left( {\tilde x\left( 0 \right),{\bf{\tilde u}}} \right)\\
s.t.\\
E\tilde u \le e - E{u_s}\\
FC\tilde x \le f - FC{x_s}
\end{array}\]===================
上述 regulation problem 的 optimal cost 為 $V^*(\tilde x(0))$ 且對應的控制力為 $\tilde {u}^0 (\tilde{x}(0))$



ref: J. B. Rawlings and D. Q. Mayne, "Model Predictive Control: Theory and Design"

2/15/2015

[線性系統] Hautus Lemma 與 控制性/可穩定性質

以下我們介紹 線性系統理論中關於 控制性的 重要結果
==========================
Lemma: Hautus Lemma for Controllability
我們稱一個離散線性系統
\[
x(k+1) = A x(k) + Bu(k)
\] 為 controllable 或稱 (A, B) controllable 若且唯若
\[
rank[\lambda I - A\;\; B] =n, \;\; \forall \lambda \in \mathbb{C}
\]其中 $\mathbb{C}$ 表示 任意 complex number 所形成的集合。
==========================
Proof. omitted.

Comments:
1. 上述 Lemma 等價為 矩陣 $[\lambda I - A\;\; B]$ 有 $n$ 個 independent rows (full row rank) 。

2. 也許讀者會好奇為何需要如此多 控制性 檢驗工具? 為什麼不使用 可控性矩陣(controllability matrix) 檢驗法就好? 事實上此 Hautus Lemma 主要是功用是大幅簡化 理論證明,但一般實際 檢驗 系統的 可控性 仍多仰賴 可控性矩陣(controllability matrix) 檢驗法。

3. 注意到上述條件需要檢驗 "任意" complex number $\lambda$ 此條件明顯過於嚴苛。不過我們可以注意到 若 $\lambda $ 並非為 $A$ 矩陣的特徵值,則 $rank[\lambda I - A\;\; B] =n$ 必然成立,故我們只需檢驗 $\lambda$ 為 $A$ 矩陣的特徵值部分即可。故

==========================
Lemma (Modified): Hautus Lemma for Controllability
一個離散線性系統
\[
x(k+1) = A x(k) + Bu(k)
\]為 controllable 或稱 (A,B) controllable 若且唯若
\[
rank[\lambda I - A\;\; B] =n, \;\; \forall \lambda \in eig(A)
\]==========================
Proof: omitted.


以下我們將討論拓展到 stabilizability。若現在我們將原本系統 $x(k+1) = A x(k) + Bu(k)$ 透過similar transformation (關於相似轉換細節請參閱 [線性系統] 控制性矩陣 與 非奇異轉換 (Controllability matrix & Non-singular transformation)) 轉成以下 partitioned system
\[\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( {k + 1} \right)}\\
{{x_2}\left( {k + 1} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{A_{11}}}&{{A_{12}}}\\
0&{{A_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( k \right)}\\
{{x_2}\left( k \right)}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
{{B_1}}\\
0
\end{array}} \right]u\left( k \right)\]其中 $(A_{11}, B_1)$ 為 controllable。上式又稱為 controllability canonical form。注意到若我們觀察此系統的 contorllability matrix
\[C = \left[ {\begin{array}{*{20}{c}}
{{B_1}}&{{A_{11}}{B_1}}&{{A_{11}}^2{B_1}}& \cdots &{{A_{11}}^n{B_1}}\\
0&0&0& \cdots &0
\end{array}} \right]\]由於 $(A_{11}, B_1)$ 為 controllable 故 controllability matrix $C$ 上排具有 $n_1$ 個 independent row。但是 $C$ 的 下排 $n_2$ rows 均為 $0$ ,故必定不滿足 full row rank test,此系統 uncontrollable。


Exercise: 試證 uncontrollable mode 為 $x_2(k)$。


儘管有 uncontrollable mode ,但我們可以退而求其次,若此 uncontrollable mode 為 stable,則我們仍可控制此系統,此類系統稱作 stabilizable。

===========
Definition: Stabilizability
考慮以下 partitioned system
\[\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( {k + 1} \right)}\\
{{x_2}\left( {k + 1} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{A_{11}}}&{{A_{12}}}\\
0&{{A_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( k \right)}\\
{{x_2}\left( k \right)}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
{{B_1}}\\
0
\end{array}} \right]u\left( k \right)\]其中 $(A_{11}, B_1)$ 為 controllable。若 $A_{22}$ 為 stable 則此 partitioned system 稱作 stabilizable。=========


接著我們可以拓展前述討論的 Hautus Lemma 到 Stabilizability 之中。
================
Lemma: (Hautus Lemma for Stabilizability)
一個離散線性系統
\[
x(k+1) = A x(k) + Bu(k)
\] 為 stabilizable 若且唯若
\[
rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1
\]================

Proof:
$(\Rightarrow)$ 假設系統 stabilizable,我們要證明  $rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1 $成立。

注意到 系統 stabilizable,故由 stabilizability 定義可知
考慮以下 partitioned system
\[\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( {k + 1} \right)}\\
{{x_2}\left( {k + 1} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{A_{11}}}&{{A_{12}}}\\
0&{{A_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( k \right)}\\
{{x_2}\left( k \right)}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
{{B_1}}\\
0
\end{array}} \right]u\left( k \right)\]其中 $(A_{11}, B_1)$ 為 controllable 且 $A_{22}$ 為 stable 。故我們有
 \[rank[\lambda I - A\;\;B] = rank\left[ {\begin{array}{*{20}{c}}
{\lambda I - {A_{11}}}&{ - {A_{12}}}&{{B_1}}\\
0&{\lambda I - {A_{22}}}&0
\end{array}} \right] \ \ \ \ \ (*)
\] 且 若我們觀察 矩陣 $[\lambda I - A_{11}\;\; B_1]$ 與 $[\lambda I - A_{22}]$ 的 row,可發現若 $|\lambda| \ge 1$ 則這些 rows 互為 independent。故由 Hautus Lemma for controllability 可知 ,對  $|\lambda| \ge 1$而言,$(*)$ 的 rows 為 independent 故  $rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1 $


$(\Leftarrow)$ 假設 $rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1 $成立,我們要證明 系統 stabilizable 。亦即 考慮以下 partitioned system
\[\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( {k + 1} \right)}\\
{{x_2}\left( {k + 1} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{A_{11}}}&{{A_{12}}}\\
0&{{A_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( k \right)}\\
{{x_2}\left( k \right)}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
{{B_1}}\\
0
\end{array}} \right]u\left( k \right)\]其中 $(A_{11}, B_1)$ 為 controllable 且 $A_{22}$ 為 stable 注意 partitioned system 為 controllability canonical form,故 $(A_{11}, B_1)$ 已為 controllable ,故我們僅須證明 $A_{22}$ 為 stable 。

由 假設 $rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1 $ 我們可知 對任意 $|\lambda| \ge 1 $,
\[rank[\lambda I - A\;\;B] = rank\left[ {\begin{array}{*{20}{c}}
{\lambda I - {A_{11}}}&{ - {A_{12}}}&{{B_1}}\\
0&{\lambda I - {A_{22}}}&0
\end{array}} \right]\]具有 independent rows 。此暗示了 上述矩陣下排 $[\lambda I - A_{22} ]$ 亦有 full row rank $\forall |\lambda| \ge 1$ 故 $A_{22}$ 的 eigenvalue 必 $<1$ 故 $A_{22} $ 為 stable。

1/15/2015

[線性系統] 離散時間系統的可觀察性質 (Observability)

考慮 離散時間 線性非時變 (Discrete Time Linear Time Invariant, DT-LTI)系統
\[\begin{array}{l}
x(k + 1) = Ax(k)\\
y(k) = Cx(k)
\end{array}\] $x \in \mathbb{R}^n, y \in \mathbb{R}^p$。

我們說上述系統為 可觀察 (observable) 或稱 $(A,C)$ 可觀察 若下列條件成立:
存在 常數 $N < \infty$ ,使得對任意 初始狀態 $x(0)$ 而言,可用 $N$ 組量測輸出 $\{y(0), y(1),...,y(N-1)\}$  uniquely 決定該初始狀態 $x(0)$。

Comment
1. 上述定義可類比 可控制性條件,

2. 事實上若我們無法透過 $n$ 組 量測輸出 來區別 $x(0)$ 則就算給額外再多的量測輸出 e.g., $N>n$ 組 仍無法區別 $x(0)$。(此結果可由 Cayley-Hamilton Theorem 證明。)

以下我們看個 unobservable 的例子
上方的方塊圖 顯示了 子系統 $G(z)$ 的狀態 無法從輸出 $Y(z)$ 觀察到。 (圖中的 $(z)$ 表示對原系統做 Z-transform。)



觀察性基本問題:
透過 sensor 所量測到的輸出 $y$ 是否足夠讓我們找出系統 初始狀態 $x(0)$ uniquely?

為何我們關心 初始狀態? 因為一但有初始狀態則其餘任意時刻狀態均可透過狀態方程求解獲得。亦即 給定 $x(0)$ 則
\[\left\{ \begin{array}{l}
x(1) = Ax(0)\\
x(2) = {A^2}x(0)\\
...\\
x\left( N \right) = {A^N}x\left( 0 \right)
\end{array} \right.\]故若給定初始狀態 $x(0)$ 則其餘任意時刻狀態 $x(1), x(2),...,x(N)$均可透過狀態方程 $x(k+1) = Ax(k)$ 獲得。

但現在我們僅給定 $y(0),...,y(N)$ 亦即我們僅知道
\[ \Rightarrow \left\{ \begin{array}{l}
y(0) = Cx(0)\\
y(1) = Cx(1) = CAx(0)\\
y(2) = Cx(2) = C{A^2}x(0)\\
...\\
y\left( N \right) = C{A^N}x\left( 0 \right)
\end{array} \right.\]或者更進一步改寫成矩陣形式
\[\left[ {\begin{array}{*{20}{c}}
{y(0)}\\
{y(1)}\\
 \vdots \\
{y(N)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
C\\
{CA}\\
 \vdots \\
{C{A^N}}
\end{array}} \right]x\left( 0 \right)\] 想問是否可從這些 $\{y(0),...,y(N)\}$ 回推 $x(0)$ (uniquely!)。 如果可以我們稱系統可觀察,若不行我們稱系統不可觀察。

Comment:
回憶在線性代數中,我們說 $Ax = b$ 解存在 若且唯若 $A$ 有 indepenent row (此對應 controllabilility problem);若我們說 $Ax = b$ 有唯一解 (注意 唯一不保證存在!!),若且為若 $A$ 有 independent column (此對應 observability problem)。

故我們要求觀察性矩陣 $O$ ( 其維度 $\dim(O) = Np  \times n$)
\[O = \left[ {\begin{array}{*{20}{c}}
C\\
{CA}\\
 \vdots \\
{C{A^N}}
\end{array}} \right]\] 有 independent column。由 Caley-Hamilton Theorem ,我們可僅考慮 $n$ 個量測輸出,則我們有
\[\left[ {\begin{array}{*{20}{c}}
{y(0)}\\
{y(1)}\\
 \vdots \\
{y(n - 1)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
C\\
{CA}\\
 \vdots \\
{C{A^{n - 1}}}
\end{array}} \right]x\left( 0 \right)\]故觀察性矩陣 $O$ ( 其維度 $\dim(O) = np  \times n$)
\[O = \left[ {\begin{array}{*{20}{c}}
C\\
{CA}\\
 \vdots \\
{C{A^{n - 1}}}
\end{array}} \right]\]需要有 full column rank $n$ 我們總結 可觀察性的測試如下:

Theorem: Observability Rank Test
若系統 $(A,C)$ 為 observable 若且唯若 $rank(O) = n$。

同樣的我們也有 Hautus Lemma for observability

Lemma: Hautus Lemma for observability
一個系統為 observable 若且唯若 對任意 $\lambda \in \mathbb{C}$,
\[rank\left[ {\begin{array}{*{20}{c}}
{\lambda I - A}\\
C
\end{array}} \right] = n\]
注意到上述 Lemma 中,若 $\lambda \notin eig(A)$,則前面 $n$ rows 為 linearly independent ,故我們可已不用檢驗整個複數平面 $\lambda \in \mathbb{C}$,僅僅需要檢驗 $\lambda \in eig(A)$ 的部分即可。故我們得到以下修正引理:

Lemma: Modified Hautus Lemma
一個系統為 observable 若且唯若 對任意 $\lambda \in eig(A)$,
\[rank\left[ {\begin{array}{*{20}{c}}
{\lambda I - A}\\
C
\end{array}} \right] = n\]







7/08/2014

[線性系統] 控制性矩陣 與 非奇異轉換 (Controllability matrix & Non-singular transformation)

延續先前線性系統理論 對於非奇異轉換的討論,由於 轉移函數 用 State space 表示實現的方法並不唯一;e.g., controllable canonical form, observable canonical form, digonal form. 故現在我們再進一步審視此問題

給定轉移函數 $H(s)$,現考慮對此轉移函數的任兩種 狀態空間實現 $\Sigma$ 與 $\tilde \Sigma$
\[\left\{ \begin{array}{l}
\Sigma  = (A,B,C,D)\\
\tilde \Sigma  = (\tilde A,\tilde B,\tilde C,\tilde D)
\end{array} \right.\],亦即
\[
H(s) = H_{\Sigma }(s) = C(sI-A)^{-1}B + D \equiv  \tilde{C} (sI- \tilde A)^{-1} \tilde B + \tilde D = H_{\tilde{\Sigma }}(s)
\]

那麼我們想知道是否存在一個 $n \times n$ 的非奇異轉換矩陣 $T$ 使得 我們有映射 $\Sigma \rightarrow \tilde \Sigma$

由先前文章可知,$\tilde A = T A T^{-1}$,$\tilde B = TB$,$\tilde C = C T^{-1}$,$\tilde D = D$,現在觀察下式
\[\left\{ \begin{array}{l}
\tilde B = TB\\
\tilde A\tilde B = \left( {TA{T^{ - 1}}} \right)TB = TAB\\
{{\tilde A}^2}\tilde B = \left( {TA{T^{ - 1}}} \right)TAB = T{A^2}B\\
 \vdots \\
{{\tilde A}^{n - 1}}\tilde B = \left( {TA{T^{ - 1}}} \right)TAB = T{A^{n - 1}}B
\end{array} \right.
\] 我們可以看出上式中一些運算的規則,現在將其改寫為更簡潔的形式如下
\[\underbrace {\left[ {\begin{array}{*{20}{c}}
{\tilde B}&{\tilde A\tilde B}& \cdots &{{{\tilde A}^{n - 1}}\tilde B}
\end{array}} \right]}_{: = {C_{\tilde \Sigma }}} = T\underbrace {\left[ {\begin{array}{*{20}{c}}
{AB}&{{A^2}B}& \cdots &{{A^{n - 1}}B}
\end{array}} \right]}_{: = {C_\Sigma }}
\] 亦即 $C_{\tilde \Sigma}= T C_{\Sigma}$ . $(\star)$

上式 $C_{\Sigma}$ 與 $C_{\tilde \Sigma}$ 稱為 控制性矩陣 (Controllability matrix)。故 非奇異矩陣 $T$ 可透過上述關係得到。

注意到如果為單輸入單輸出 (SISO) 系統,且假設  $C_{\Sigma}$ 與 $C_{\tilde \Sigma}$ 為方陣。若 $C_{\Sigma}$ 為 non-singular,則我們可以找到非奇異轉換矩陣 $T$
\[
T = C_{\tilde \Sigma}C_{\Sigma}^{-1}
\]

若 多輸入系統,則無法直接求解反矩陣,故我們需先使 $(\star)$ 左右變成方陣:
\[
\underbrace {{C_{\tilde \Sigma }}{C_\Sigma }^T}_{\underbrace {\left( {n \times nm} \right) \times \left( {mn\times n} \right)}_{n \times n}} = T\underbrace {{C_\Sigma }{C_\Sigma }^T}_{\underbrace {\left( {n \times nm} \right) \times \left( {mn \times n} \right)}_{n \times n}}
\]現在 若  ${{C_\Sigma }{C_\Sigma }}$ 為 non-singular,則我們可以找到非奇異轉換矩陣 $T$
\[
T = {C_{\tilde \Sigma }}{C_\Sigma }^T{\left( {{C_\Sigma }{C_\Sigma }^T} \right)^{ - 1}}
\]

故 我們知道如果要有 非奇異矩陣 $T$,則矩陣 $C_{\Sigma} C_{\Sigma }^T$ 必須非奇異,故我們有下列 Controllability Rank conditon:

Controllability Rank Condition
 $C_{\Sigma} C_{\Sigma }^T$ 為非奇異 若且為若 $\text{rank}{C_{\Sigma}} = n$


Comment:
1. 在 MATLAB中,給定動態系統 $A, B$ 矩陣,則我們可以使用  C = ctrb(A,B) 指令來直接幫助我們計算 Controllability Matrix, C,接著再用 rank(C) 指令確認此矩陣是否滿足我們的 Controllability Rank Condition ,如果滿足我們稱此系統為可控制(controllable)。

2. non-singular transform 不改變 Eigenvalues,亦即
\[eig\left( {TA{T^{ - 1}}} \right) = eig\left( A \right)
\]其中 $eig(\cdot)$ 表特徵值。
Proof
令 $T$ 為 nonsingular transformation matrix,且 $\lambda_i$ 為 $TAT^{-1}$ 矩陣對應的 eigenvalue,也就是說 $TAT^{-1}$ 的 eigenvalues 滿足 $\det(\lambda_i I - TAT^{-1}) =0$。故
\[\begin{array}{l}
\det \left( {{\lambda _i}I - TA{T^{ - 1}}} \right) = 0\\
 \Rightarrow \det \left( {{\lambda _i}T{T^{ - 1}} - TA{T^{ - 1}}} \right) = 0\\
 \Rightarrow \det \left( {T\left( {{\lambda _i}I - A} \right){T^{ - 1}}} \right) = 0\\
 \Rightarrow \det \left( T \right)\det \left( {{\lambda _i}I - A} \right)\det \left( {{T^{ - 1}}} \right) = 0
\end{array}\]由於 $T$ 為 nonsingular,故 $T^{-1}$ 存在且 $\det(T) \neq 0$,  $\det(T^{-1}) \neq 0$。故只有
\[
\det(\lambda_i I - A) =0
\]亦即 $\lambda_i$ 亦為 矩陣 $A$ 的 eigenvalue。 $\square$

[線性系統] 實現定理 與 非奇異轉換

這次要介紹 線性系統理論 中的一個重要結果:稱作實現理論 ( Realization Theory )

考慮一個轉移函數 $H(s)$ 可以將其由 狀態空間表示,我們記做 $\Sigma$。
其中
\[\Sigma : = \left\{ \begin{array}{l}
\dot x = Ax + Bu\\
y = Cx + Du
\end{array} \right.
\] 則我們有以下定義:
==============
Definition: Realization
令 $H(s)$ 為給定轉移函數,則我們說 其狀態空間 $\Sigma $  為 $H(s)$ 的實現 (Realization) 若下列條件成立:
\[
C(sI-A)^{-1}B + D = H(s)
\]
=============
Comments:
1. 上述 實現(Realization) 意指可以透過 實體電路 (e.g., OP放大器等) "實現" 狀態方程。
2. 設 $\sum = (A,B,C,D) $ 為 $H(s)$ 的實現,現在定義 $T$ 為任意 $n \times n$ 的非奇異矩陣 (non-singular matrix),則我們可以定義下列 新系統 以狀態空間表示:
\[
\tilde {\sum} = (\tilde A, \tilde B, \tilde C, \tilde D)
\]
其中 $\tilde A = TAT^{-1}$, $\tilde B = TB$, $\tilde C = CT^{-1}$, $\tilde D = D$。

那麼現在我們來看看此新系統的轉移函數為何?

\[\begin{array}{l}
\tilde C{(sI - \tilde A)^{ - 1}}\tilde B + \tilde D = C{T^{ - 1}}{(sI - TA{T^{ - 1}})^{ - 1}}TB + D\\
 \ \ \ \ \ \ \ \ = C{T^{ - 1}}{(sT{T^{ - 1}} - TA{T^{ - 1}})^{ - 1}}TB + D\\
  \ \ \ \ \ \ \ \ = C{T^{ - 1}}{\left( {T\left( {sI - A} \right){T^{ - 1}}} \right)^{ - 1}}TB + D\\
 \ \ \ \ \ \ \ \  = C{T^{ - 1}}\left( {T{{\left( {sI - A} \right)}^{ - 1}}{T^{ - 1}}} \right)TB + D\\
  \ \ \ \ \ \ \ \ = C{\left( {sI - A} \right)^{ - 1}}B + D \\
  \ \ \ \ \ \ \ \ = H(s)
\end{array}\]

上述結果告訴我們

1. 狀態空間表示 若透過 非奇異轉換 (Non-singular transformation),其轉移函數不變 (invariant)

2. 上式non-singular transformation 等價於 將系統以新的狀態變數 $z := T x$ 改寫。
由於 $z = Tx \Rightarrow \dot z = T\dot x \Rightarrow \dot x = {T^{ - 1}}\dot z$,故原系統狀態表示可改寫為
\[\begin{array}{l}
\left\{ \begin{array}{l}
\dot x = Ax + Bu\\
y = Cx + Du
\end{array} \right. \Rightarrow \left\{ \begin{array}{l}
{T^{ - 1}}\dot z = A{T^{ - 1}}z + Bu\\
y = C{T^{ - 1}}z + Du
\end{array} \right.\\
 \Rightarrow \left\{ \begin{array}{l}
\dot z = \underbrace {TA{T^{ - 1}}}_{\tilde A}z + \underbrace {TB}_{\tilde B}u\\
y = \underbrace {C{T^{ - 1}}}_{\tilde C}z + \underbrace D_{\tilde D}u
\end{array} \right.
\end{array}\]

有了上述結果之後,我們知道同一系統的 任意狀態空間實現 都可透過 非奇異轉換 求得相同的轉移函數,那麼現在問題變成怎樣的轉移函數才可以被實現??

以下我們給出一個重要且簡潔的定理來回答這個問題:

=======================
Theorem: Realization Theorem
任意 proper (分母階數大於或者等於分子階數) 轉移函數皆為可實現 (realizable)。 
=======================

那麼問題變成已知 proper 轉移函數可以實現 (有狀態空間表示),那麼該如何實現呢? 我們用下面這個例子來說明:

現在考慮 轉移函數
\[
H(s) = \frac{b_m s^m + b_{m-1}s^{m-1} + ... + b_1 s^1 + b_0}{s^n + a_{n-1}s^{n-1} + a_{n-2}s^{n-2} + ... + a_1 s^1 + a_0} + r
\] 其中 $m < n$ (properness)

在不失一般性的情況,我們設 $m = n-1$,則我們由 Realization Theorem 可知 此轉移函數存在 狀態空間表示 (可以實現),故我們可寫成
\[\begin{array}{l}
A = \left[ {\begin{array}{*{20}{c}}
0&1&0&0& \cdots &0\\
0&0&1&0& \cdots &0\\
0&0&0&1&{0 \cdots }&0\\
 \vdots & \vdots &{}&\begin{array}{l}
0\\
 \vdots
\end{array}& \ddots & \vdots \\
0&0& \cdots & \cdots &0&1\\
{ - {a_0}}&{ - {a_1}}&{ - {a_2}}& \cdots &{ - {a_{n - 2}}}&{ - {a_{n - 1}}}
\end{array}} \right],B = \left[ {\begin{array}{*{20}{c}}
0\\
0\\
0\\
 \vdots \\
0\\
1
\end{array}} \right]\\
C = \left[ {\begin{array}{*{20}{c}}
{{b_0}}&{{b_1}}&{{b_2}}& \cdots &{{b_{m - 1}}}&{{b_m}}
\end{array}} \right]\\
D = r
\end{array}\]
上式實現 稱為 可控典型式 (Controllable Canonical form)。

Comments:
1. 為何上述實現被稱為可控典型式?

觀察上式
\[\begin{array}{l}
\dot x = Ax + Bu = A\\
 \Rightarrow \dot x = \left[ {\begin{array}{*{20}{c}}
0&1&0&0& \cdots &0\\
0&0&1&0& \cdots &0\\
0&0&0&1&{0 \cdots }&0\\
 \vdots & \vdots &{}&{\begin{array}{*{20}{l}}
0\\
 \vdots
\end{array}}& \ddots & \vdots \\
0&0& \cdots & \cdots &0&1\\
{ - {a_0}}&{ - {a_1}}&{ - {a_2}}& \cdots &{ - {a_{n - 2}}}&{ - {a_{n - 1}}}
\end{array}} \right]x + B = \left[ {\begin{array}{*{20}{c}}
0\\
0\\
0\\
 \vdots \\
0\\
1
\end{array}} \right]u
\end{array}
\] 現在計算對應的特徵方程式 $\det( sI - A)$,我們可得
\[
\det(sI-A) = s^n + a_{n-1}s^{n-1} + ... + a_0
\] 現在如果我們讓控制力 $ u = Kx$亦即
\[u = \left[ {\begin{array}{*{20}{c}}
{{k_1}}&{{k_2}}& \cdots &{{k_n}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}}\\
{{x_2}}\\
 \vdots \\
{{x_n}}
\end{array}} \right]
\]則 受控制的動態系統可以改寫為
\[
\dot x = Ax + Bu = Ax + B(Kx ) = (A+BK)x
\]此時
\[\begin{array}{l}
 \Rightarrow A + BK = \left[ {\begin{array}{*{20}{c}}
0&1&0& \cdots &0\\
0&0&1&{0 \cdots }& \vdots \\
0&0&0& \ddots &0\\
 \vdots & \vdots &{}&0&1\\
{ - {a_0}}&{ - {a_1}}& \cdots & \cdots &{ - {a_{n - 1}}}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
0\\
0\\
 \vdots \\
0\\
1
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{k_1}}&{{k_2}}& \cdots &{{k_n}}
\end{array}} \right]\\
\begin{array}{*{20}{c}}
{}
\end{array}\begin{array}{*{20}{c}}
{}
\end{array}\begin{array}{*{20}{c}}
{}
\end{array}\begin{array}{*{20}{c}}
{}
\end{array} = \left[ {\begin{array}{*{20}{c}}
0&1&0& \cdots &0\\
0&0&1&{0 \cdots }& \vdots \\
0&0&0& \ddots &0\\
 \vdots & \vdots &{}&0&1\\
{{k_1} - {a_0}}&{{k_2} - {a_1}}& \cdots & \cdots &{{k_n} - {a_{n - 1}}}
\end{array}} \right]
\end{array}\]上式可以發對每一個參數 $a_i, \forall i =0, ...,n$ 都有一個對應的控制力參數 $k_j, j=1,...,n$來與之調整,故對應的特徵方程 $\det(sI-(A+BK))$ 的特性根根 (亦即 poles)亦會被 $K$ 直接。此poles 的位置將直接影響到系統性能,故如果某動態系統可寫為可控典型式,則我們可透過上述的控制力 $u=Kx $ 直接改變每一個系統的特性根位置。

2.
在 MATLAB 中 由轉移函數轉成狀態空間實現,可以透過指令 tf2ss.m 來達成。在此不贅述

以下我們看個例子:

Example
考慮轉移函數
\[G(s) = \frac{Y(s)}{U(s)}= \frac{{{b_2}{s^2} + {b_1}{s^1} + {b_0}}}{{a_3^{}{s^3} + {a_2}{s^2} + {a_1}{s^1} + {a_0}}} + r
\]其中 $r$ 為常數。試求出 controllable canonical form:
Solution
注意到我們有額外的常數 $r$ 故可知 $D =r$ (此額外的項,表示輸入可直接影響輸出)

故我們只需專心在 strictly proper 的轉移函數部分即可。另外此例由於階數較低,我們可以用推導的方式求得 controllable canonical form。現在我們觀察轉移函數,並將其繪製成方塊圖

其中我們引入中繼函數 $X(s)$,則透過上圖我們可將轉移函數改寫回微分方程如下
\[\left\{ \begin{array}{l}
\frac{{X(s)}}{{U(s)}} = \frac{1}{{{s^3} + {a_2}{s^2} + {a_1}{s^1} + {a_0}}}\\
\frac{{Y(s)}}{{X(s)}} = {b_2}{s^2} + {b_1}{s^1} + {b_0}
\end{array} \right. \Rightarrow \left\{ \begin{array}{l}
{x^{\left( 3 \right)}} + {a_2}\ddot x + {a_1}\dot x + {a_0}x = u\\
{b_2}\ddot x + {b_1}\dot x + {b_0}x = y
\end{array} \right.\]現在我們定義狀態 $x: = {x_1},\dot x: = {x_2},\ddot x: = {x_3}$ 則上式改寫如下
\[\left\{ \begin{array}{l}
a_3^{}{x^{\left( 3 \right)}} + {a_2}\ddot x + {a_1}\dot x + {a_0}x = u\\
{b_2}\ddot x + {b_1}\dot x + {b_0}x = y
\end{array} \right. \Rightarrow \left\{ \begin{array}{l}
{{\dot x}_3} + {a_2}{x_3} + {a_1}{x_2} + {a_0}{x_1} = u\\
{b_2}{x_3} + {b_1}{x_2} + {b_0}{x_1} = y
\end{array} \right.\]且我們有 ${{\dot x}_1} = {x_2},\;{{\dot x}_2} = {x_3}$ 故我們可寫成
\[\begin{array}{l}
\dot x = Ax + Bu\\
y = Cx + Du
\end{array}\]如下
\[\left\{ \begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{{{\dot x}_1}}\\
{{{\dot x}_2}}\\
{{{\dot x}_3}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
0&1&0\\
0&0&1\\
{ - {a_0}}&{ - {a_1}}&{ - {a_2}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}}\\
{{x_2}}\\
{{x_3}}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
0\\
0\\
1
\end{array}} \right]u\\
y = \left[ {\begin{array}{*{20}{c}}
{{b_0}}&{{b_1}}&{{b_2}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}}\\
{{x_2}}\\
{{x_3}}
\end{array}} \right]
\end{array} \right.\]上式即為 controllable canonical form。

現在合併先前我們的 $D=r$ 故可得最終表示為
\[\left\{ \begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{{{\dot x}_1}}\\
{{{\dot x}_2}}\\
{{{\dot x}_3}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
0&1&0\\
0&0&1\\
{ - {a_0}}&{ - {a_1}}&{ - {a_2}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}}\\
{{x_2}}\\
{{x_3}}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
0\\
0\\
1
\end{array}} \right]u\\
y = \left[ {\begin{array}{*{20}{c}}
{{b_0}}&{{b_1}}&{{b_2}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}}\\
{{x_2}}\\
{{x_3}}
\end{array}} \right] + ru \ \ \ \ \ \ \ \square
\end{array} \right.\]

上述可控典型式 與 實現定裡之間關係 我們會留待下一篇文章在做介紹。
[線性系統] Controllability Matrix

另外亦會對非奇異轉換矩陣 $T$ 的求得?? 也就是是否可以找到一個非奇異轉換矩陣 來幫助我們從一個狀態空間的實現 變成 另一個呢?? 做額外補充。

7/06/2014

[線性系統] 線性動態系統的表示法: 轉移函數 與 狀態空間表示

這次要介紹線性系統理論中對於動態系統的表示方法:
一般而言,線性動態系統 可以用 線性微分方程(O.D.E.) 來表達,但在控制理論中亦提供兩種不同的方法來表達動態系統:

一種稱為 轉移函數(transfer function) 表示法 (主要工具為 拉式轉換)
一種稱為 狀態空間(state space) 表示法 (主要工具為 矩陣線性代數)

那麼同一種 動態系統間,不同的表示法 可以互相等價轉換,現在我們先看個例子:

Example: Dynamic System to Transfer function
考慮下列動態系統微分方程
\[
\frac{{{d^3}y}}{{d{t^3}}} + 6\frac{{{d^2}y}}{{d{t^2}}} + 5\frac{{dy}}{{dt}} - 4y = u\left( t \right) + 2\frac{{du\left( t \right)}}{{dt}}
\] 那麼我們可以對其取拉式轉換(Laplace Transform) $\mathcal{L}(\cdot)$ 來求取轉移函數,亦即
\[\begin{array}{l}
{{\cal L}}\left\{ {\frac{{{d^3}y}}{{d{t^3}}} + 6\frac{{{d^2}y}}{{d{t^2}}} + 5\frac{{dy}}{{dt}} - 4y} \right\} = {{\cal L}}\left\{ {u\left( t \right) + 2\frac{{du\left( t \right)}}{{dt}}} \right\}\\
 \Rightarrow \left\{ \begin{array}{l}
{s^3}Y\left( s \right) - {s^2}y\left( 0 \right) - s{y^{\left( 1 \right)}}\left( 0 \right) - {y^{\left( 2 \right)}}\left( 0 \right)\\
 + 6\left( {{s^2}Y\left( s \right) - sy\left( 0 \right) - {y^{\left( 1 \right)}}\left( 0 \right)} \right)\\
 + 5\left( {sY\left( s \right) - y\left( 0 \right)} \right) \\
- 4Y\left( s \right)
\end{array} \right\} = \left\{ \begin{array}{l}
U\left( s \right)\\
 + 2\left( {sU\left( s \right) - u\left( 0 \right)} \right)
\end{array} \right\}
\end{array}
\] 上述中 $s$ 表示 微分器;反之 $s^{-1}$ 稱之為積分器。

現在考慮 $y^{(k)} =0, \forall k =0,1,2,...$ 且 $u(0) =0$ 亦即我們考慮整個動態系統的初始狀態為休止 (initially at rest),且亦無初始控制力,則上述拉式轉換式可得
\[\begin{array}{l}
{s^3}Y\left( s \right) + 6{s^2}Y\left( s \right) + 5sY\left( s \right) - 4Y\left( s \right) = U\left( s \right) + 2sU\left( s \right)\\
 \Rightarrow \left( {{s^3} + 6{s^2} + 5s - 4} \right)Y\left( s \right) = \left( {1 + 2s} \right)U\left( s \right)\\
 \Rightarrow \frac{{Y\left( s \right)}}{{U\left( s \right)}} = \frac{{ 2s + 1}}{{{s^3} + 6{s^2} + 5s - 4}}
\end{array}
\] 我們稱上述 $H(s) := \frac{Y(s)}{U(s)}$ 為轉移函數 transfer function。接著除了 轉移函數表示法之外,我們亦可將其用矩陣的方式表示:這邊僅簡單介紹可控典型式(controllable canonical form):
考慮上述轉移函數
\[H\left( s \right) = \frac{{Y\left( s \right)}}{{U\left( s \right)}} = \frac{{1 + 2s}}{{{s^3} + 6{s^2} + 5s - 4}}\]可將其改寫為以下狀態空間模型
\[\left\{ \begin{array}{l}
\dot x = Ax + Bu = \left[ {\begin{array}{*{20}{c}}
0&1&0\\
0&0&1\\
4&{ - 5}&{ - 6}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}}\\
{{x_2}}\\
{{x_3}}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
0\\
0\\
1
\end{array}} \right]u\\
y = Cx = \left[ {\begin{array}{*{20}{c}}
1&2&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}}\\
{{x_2}}\\
{{x_3}}
\end{array}} \right]
\end{array} \right.\]

Comments:
1. 一般在 MATLAB中,建構轉移函數可以使用 tf(NUM,DEN) 指令,其中 NUM 表示分子係數,DEN表示分母係數:以上例而言,轉移函數透過 MATLAB 建構為

tf( [2 1], [1 6 5 -4])

2. 在MATLAB 中,如果已知轉移函數欲將其轉換到狀態空間模型 (亦即 欲得到 $A,B,C,D$ 矩陣) 有一個非常簡便的指令:[A,B,C,D] = tf2ss(NUM,DEN)



有了上述例子之後我們可以回頭看看 如何從 狀態空間模型 來 求得 轉移函數:

State-Space Method to Transfer Function

現在我們考慮狀態空間模型:
\[\left\{ \begin{array}{l}
\dot x = Ax + Bu\\
y = Cx + Du
\end{array} \right.
\]其中 $x$ 為 $n \times 1$狀態變數向量,$u$ 為 $m \times 1$ 控制力向量,$y$ 為 $r \times 1$ 輸出向量,$A$為 $n \times n$ 矩陣,$B$ 為 $n \times m$ 矩陣, $C$ 為 $r \times n$ 矩陣,$D$ 為 $r \times m $矩陣。

對上式取拉式轉換 $\cal{L}(\cdot)$ 並令初值為零,則可得
\[\left\{ \begin{array}{l}
sX\left( s \right) = AX\left( s \right) + BU\left( s \right)\\
Y\left( s \right) = CX\left( s \right) + DU\left( s \right)
\end{array} \right.
\]現在整理上式可得
\[\begin{array}{l}
\left\{ \begin{array}{l}
\left( {sI - A} \right)X\left( s \right) = BU\left( s \right)\\
Y\left( s \right) = CX\left( s \right) + DU\left( s \right)
\end{array} \right.\\
 \Rightarrow \left\{ \begin{array}{l}
X\left( s \right) = {\left( {sI - A} \right)^{ - 1}}BU\left( s \right)\\
Y\left( s \right) = CX\left( s \right) + DU\left( s \right)
\end{array} \right.\\
 \Rightarrow Y\left( s \right) = C{\left( {sI - A} \right)^{ - 1}}BU\left( s \right) + DU\left( s \right)\\
 \Rightarrow Y\left( s \right) = \left[ {C{{\left( {sI - A} \right)}^{ - 1}}B + D} \right]U\left( s \right)\\
 \Rightarrow \frac{{Y\left( s \right)}}{{U\left( s \right)}} = \underbrace {C{{\left( {sI - A} \right)}^{ - 1}}B + D}_{H\left( s \right)}
\end{array}
\]上述 $H(s)$ 即為轉移函數

且注意到
\[ \Rightarrow \frac{{Y\left( s \right)}}{{U\left( s \right)}} = \underbrace {C{{\left( {sI - A} \right)}^{ - 1}}B + D}_{H\left( s \right)} = C\frac{{adj\left( {sI - A} \right)}}{{\det \left( {sI - A} \right)}}B + D\]故 轉移函數 $H(s)$ 的分母等於 $\det (sI-A)$ 亦即 $\det(sI-A)=0$為系統特徵方程;且 $H(s)$ 的 pole 等於 $A$ 矩陣的 特徵值(eigenvalue)。


Comments:
1. 對線性動態系統而言,狀態空間表示法並非唯一 (故選取的 $A,B,C,D$ 矩陣稱為轉移函數的 實現 realization)。

2. 轉移函數 $H(s)$ 為唯一。亦即轉移函數具備不變性 (invariant).

3. 若 轉移函數 $H(s)$ 分母階數 $\ge$ 分子階數,我們稱此轉移函數為 proper。若 分母階數 $>$ 分子階數,則稱此轉移函數為 strictly proper。若 分母階數 $<$ 分子階數,稱此轉移函數為 improper。同理我們可直接對矩陣形式做判斷
\[ \Rightarrow \frac{{Y\left( s \right)}}{{U\left( s \right)}} = \underbrace {C{{\left( {sI - A} \right)}^{ - 1}}B + D}_{H\left( s \right)}\]若 $D=0$,則 $H(s)$ 為 stirctly proper ;若 $D \neq  0$ 則 $H(s)$ 並非 strictly proper。

對於上述 comments 有興趣的讀者請閱讀
 [線性系統] Realization Theory and Non-singular Transformation




Example 1.: The simplest improper transfer function
\[
H(s) = s
\]亦即,若轉移函數為一 微分器 ,則此時分母為常數 $1$,階數為0階 小於 分子 $s$ 的一階。故為 improper transfer function。

Example 2
考慮系統
\[\left\{ \begin{array}{l}
\dot x = Ax + Bu = \left[ {\begin{array}{*{20}{c}}
{ - 0.1}&1\\
0&{ - 1}
\end{array}} \right]x + \left[ {\begin{array}{*{20}{c}}
1\\
1
\end{array}} \right]u\\
y = Cx = \left[ {\begin{array}{*{20}{c}}
1&0
\end{array}} \right]x
\end{array} \right.
\](a) 試求 系統輸入輸出轉移函數。
(b) 若 $y(t) = 1 - {e^{ - t}}$ 且 $x(0)=0$ 試求對應的 $u(t)$。
Solution
(a):
\[\begin{array}{l}
\frac{{Y\left( s \right)}}{{U\left( s \right)}} = \underbrace {C{{\left( {sI - A} \right)}^{ - 1}}B + D}_{H\left( s \right)}\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \left[ {\begin{array}{*{20}{c}}
1&0
\end{array}} \right]{\left( {\left[ {\begin{array}{*{20}{c}}
s&0\\
0&s
\end{array}} \right] - \left[ {\begin{array}{*{20}{c}}
{ - 0.1}&1\\
0&{ - 1}
\end{array}} \right]} \right)^{ - 1}}\left[ {\begin{array}{*{20}{c}}
1\\
1
\end{array}} \right] + 0\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \left[ {\begin{array}{*{20}{c}}
1&0
\end{array}} \right]{\left[ {\begin{array}{*{20}{c}}
{s + 0.1}&{ - 1}\\
0&{s + 1}
\end{array}} \right]^{ - 1}}\left[ {\begin{array}{*{20}{c}}
1\\
1
\end{array}} \right]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \left[ {\begin{array}{*{20}{c}}
1&0
\end{array}} \right]\frac{1}{{\left( {s + 0.1} \right)\left( {s + 1} \right)}}\left[ {\begin{array}{*{20}{c}}
{s + 1}&1\\
0&{s + 0.1}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
1\\
1
\end{array}} \right]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \frac{{s + 2}}{{\left( {s + 0.1} \right)\left( {s + 1} \right)}}
\end{array}\]
(b):由於 $y(t) = 1 - \frac{10}{9} e^{-t} + \frac{1}{9} e^{-10t}$ 對此取拉式轉換可得
\[y(t) = 1 - {e^{ - t}} \Rightarrow Y\left( s \right) = \frac{1}{s} - \frac{1}{{s + 1}}\]故由輸入與輸出關係
\[\begin{array}{l}
H\left( s \right) = \frac{{Y\left( s \right)}}{{U\left( s \right)}} = \frac{{s + 2}}{{\left( {s + 0.1} \right)\left( {s + 1} \right)}} \Rightarrow Y\left( s \right) = \frac{{s + 2}}{{\left( {s + 0.1} \right)\left( {s + 1} \right)}}U\left( s \right)\\
 \Rightarrow \frac{1}{s} - \frac{1}{{s + 1}} = \frac{{s + 2}}{{\left( {s + 0.1} \right)\left( {s + 1} \right)}}U\left( s \right)\\
 \Rightarrow U\left( s \right) = \frac{{\left( {s + 0.1} \right)\left( {s + 1} \right)}}{{s\left( {s + 2} \right)}} - \frac{{s + 0.1}}{{s + 2}}
\end{array}\]再取反拉式轉換即可求得所需結果。

12/24/2013

[控制理論] 線性化(Linearization)

這次要跟大家介紹的是 線性化 (Linearization) 的概念,讀者建議須先具備基本 Taylor Series 概念,如果不熟悉的讀者可先參閱 [微積分] 泰勒展開式 與 泰勒級數 。

為何要做線性化?
其實線性化的動機很簡單,主要是因為一般在分析動態系統的時候,大部分系統行為都是呈現非線性(EX: 電路系統(二極體 I/V curve),倒單擺、撓性機構、機器人、生物細胞、金融模型...),但這些非線性行為會有一個大的困難,就是難以直接求解其動態行為。且發展成熟的線性系統理論沒有辦法(有效的)應用在上面,但如果能夠透過一些假設/機制,我們可以把原本非線性的系統轉成線性系統,如此一來原本沒辦法使用的線性系統理論便可以派上用場!!

如何做線性化?
至於實際如何做到對任意 非線性函數 (e.g., $\sin, \cos, \exp, x^n$, ...)線性化呢? 簡單來說,就是採用切線 (微分) 的概念,如果我們對關心的某一點對該點取導數,則我們可以得到一條對該點的切線,此切線可以在某種程度上用來近似 該點附近的函數行為。

https://controls.engin.umich.edu/wiki/index.php/LinearizingODEs


----- 以下進入正題 ----

若用數學來描述非線性的系統可以寫成
\[
\dot x(t) = f(x)
\]其中 $x(t) \in \mathbb{R}^n$ 稱作系統狀態(state variable) (這邊考慮 $n$ 維空間,故有 $n$ 個系統狀態變數); $\dot x(t)$ 為系統狀態的一階導數; $f$ 為用以描述動態系統的任意函數

在此我們考慮系統狀態為 $n$ 階。意思就是有 $n$ 個不同的系統狀態,記做  $x \in \mathbb{R}^n$


在介紹線性化之前,我們得先介紹 "平衡點(equilibrium point)":

=====================
Definition: Equilibrium point
若 $f(\bar{x})=0$ ,則系統狀態 $\bar{x} \in \mathbb{R}^n$ 被稱作 平衡點(equilibrium point)。
=====================

Comments
由上述定義可以推知,如果 $x(0)=\bar{x}$ ,則 $x(t)=\bar{x}, \forall t \geq 0$ ;
也就是說一旦 在最一開始( $t=0$ )的時候,系統就處在平衡點的狀態,則對任意時刻 $t \geq 0$,系統狀態會持續處在平衡點的狀態。

====================
Definition: 穩定平衡點 (Stable equilibrium point)
平衡點若被稱為穩定的,或稱 穩定平衡點(stable equilibrium point),則其必須滿足 在任意時刻 之狀態 $x(t)$ 都需收斂到平衡點 $\bar{x}$,亦即
\[
x(t) \rightarrow \bar{x}
\](  $|| x(0)-\bar{x} ||$ 為足夠小 )
反之,若不收斂則稱為 不穩定的平衡點(unstable equilibrium) 其中
 \[
|| x(0)-\bar{x} || := \left ( \displaystyle \sum_{i} (x_i(0)-\bar{x}_i)^2 \right )^{\frac{1}{2}}\] 為 2-norm
==================

Comments:
上述定義指明 所謂的 穩定平衡點是指 考慮 任意時刻的狀態 $x(t), \forall t$,若此狀態都會回到 某個平衡點 $\bar{x}$ 則我們說他是一個 穩定的平衡點。

在介紹完平衡點之後,我們便可介紹所謂的 線性化,誠如先前所說,線性化的基本概念是微分,所以在此我們會假設動態系統充分可微(smooth),故我們可以進行 泰勒展開 (微分近似)。

=================
線性化(Linearization):
注意:我們僅對 穩定平衡點 做線性化。(不穩定的平衡點亦可線性化只是實際用處不大)

現在我們回頭考慮 $n$ 階非線性系統
\[
\dot x(t) = f(x)
\]其系統狀態 $x(t)$ 可表為 平衡點狀態 $\bar{x}$ 加上 狀態(小擾動)增量( $\Delta x(t)$ );注意。在此我們假設擾動 $\Delta x$ 不能太大。
\[
x(t) = \bar{x} + \Delta x(t)
\]則我們可寫下
\[
f_1 (x) = f_1 (x+\Delta x)\]
若此 $f$ 為 平滑函數(smooth) (也就是說可以對其做泰勒展開),則我們可改寫上式如下:

對 $f_1$ 可寫出其泰勒展開式 (對 $0$ 點展開)

$\Rightarrow f_1(x) = f_1(\bar{x}) + \frac{\partial f_1}{\partial x_1}|_{x=\bar{x}} \Delta x_1  + ... + \frac{\partial f_1}{\partial x_n}|_{x=\bar{x}} \Delta x_n + H.O.T$ ....(1)

同樣的,我們也可以對 $f_2...f_n$ 展開。

$\Rightarrow f_2(x) = f_2(\bar{x}) + \frac{\partial f_2}{\partial x_1}|_{x=\bar{x}} \Delta x_1  + ... + \frac{\partial f_2}{\partial x_n}|_{x=\bar{x}} \Delta x_n + H.O.T.$

$\vdots$

$\Rightarrow f_n(x) = f_n(\bar{x}) + \frac{\partial f_n}{\partial x_1}|_{x=\bar{x}} \Delta x_1  + ... + \frac{\partial f_n}{\partial x_n}|_{x=\bar{x}} \Delta x_n + H.O.T$

其中 $H.O.T$ 表示 高階項(Higher Order Terms)

然後因為增量 $\Delta x_1, \Delta x_2...$ 假設為很小的擾動,在高階項的影響可被忽略

現在,回憶我們手邊有的狀態
\[
x(t) = \bar{x} + \Delta x(t)\]
對上式兩邊對時間微分,可得
\[
\dot x(t) = \Delta \dot x(t)\]
再者,因為我們知道 系統為 $\dot x(t) = f(x)$, 由式 (1) 我們可以帶入泰勒展開到 $f(x)$ 之中,最後整理可得線性化之後的 增量(擾動)系統
\[
 \Delta \dot x(t) = A \cdot \Delta x(t)  \ \ \ \  (2) \]
其中 $A$ 為矩陣其第 (i,j) 元素由下式表示
\[
a_{ij} = \frac{\partial f_i}{\partial x_j}|_{x=\bar{x}}\]

上式 $(2)$ 即稱為 線性化後的動態系統。由於此為線性,故所有的線性系統理論 (eigenvalue, controllability, observability) 便可以在其上進行討論




8/25/2013

[系統理論] (弱)平穩隨機過程特性 與 線性非時變系統 (2)- Jointly wide-sense stationary and Frequency domain property

回憶對於 具有脈衝響應為 $h(t)$ 的 LTI 系統 而言,若輸入 $X_t$ 為 WSS 隨機過程,則輸出 $Y_t$ 必定也為一 WSS 隨機過程,且 $Y_t$ 與 $X_t$ 的輸入輸出關係為 convolution integral:
\[{Y_t} = \int_{ - \infty }^\infty  h (t - \tau ){X_\tau }d\tau  = \int_{ - \infty }^\infty  {{X_{t - \tau }}h} (\tau )d\tau
\]且其 輸出隨機過程 $Y_t$ 的 mean function $m_Y(t)$ 與 autocorrelation function $R_Y(\tau)$為
\[\begin{array}{l}
{m_Y}\left( t\right) = E[{Y_t}] = E[\int_{ - \infty }^\infty  {{X_{t - \theta }}h} (\theta )d\theta ]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \int_{ - \infty }^\infty  {\underbrace {E[{X_{t - \theta }}]}_{ = {m_X}\left( {t - \theta } \right)}h} (\theta )d\theta \\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \int_{ - \infty }^\infty  {{m_X}\left( {t - \theta } \right)h} (\theta )d\theta
\end{array}\]與
\[\begin{array}{l}
{R_Y}\left( \tau  \right) = E[{Y_{t + \tau }}{Y_t}]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = E\left[ {\int_{ - \infty }^\infty  {{X_{t + \tau  - \beta }}h} (\beta )d\beta \int_{ - \infty }^\infty  {{X_{t - \theta }}h} (\theta )d\theta } \right]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = E\left[ {\int_{ - \infty }^\infty  {\left( {\int_{ - \infty }^\infty  {{X_{t + \tau  - \beta }}h(\beta )} d\beta } \right){X_{t - \theta }}h(\theta )d\theta } } \right]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \int_{ - \infty }^\infty  {\left( {\int_{ - \infty }^\infty  {\underbrace {E[{X_{t + \tau  - \beta }}{X_{t - \theta }}]}_{ = {R_X}\left( {\left( {t + \tau  - \beta } \right) - \left( {t - \theta } \right)} \right)}} h(\beta )d\beta } \right)h(\theta )d\theta } \\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \int_{ - \infty }^\infty  {\left( {\int_{ - \infty }^\infty  {{R_X}\left( {\tau  - \left( {\beta  - \theta } \right)} \right)} h(\beta )d\beta } \right)h(\theta )d\theta }  \ \ \ \ \ (*)
\end{array}\]

現在我們定義 jointly WSS 隨機過程
============================
Definition (Jointly wide-sense stationary, J-WSS)
考慮 $X_t$ 與 $Y_t$ 為兩 WSS 隨機過程,現給定任意兩時刻 $t_1, t_2$, 若其 cross-crorelation $E[X_{t_1}Y_{t_2}]$ 只與 $t_1 - t_2$ 有關,則我們稱 $X_t$ 與 $Y_t$ 為 Jointly wide-sense stationary。
============================

現在回憶 WSS隨機過程之 auto-correlation function 可以用單變數改寫,同樣的,J-WSS 的 $X_t$ 與 $Y_t$ 亦可改寫成單變數 univariate cross-correlation function ;亦即,令 $t_1 := t + \tau$ 且 $t_2 := t$ 則我們有
\[
R_{XY}(t_1 - t_2) = E[X_{t_1} Y_{t_2}] = E[X_{t+ \tau} Y_{t}] = R_{XY}(t+ \tau - t) = R_XY(\tau)
\]亦即 J-WSS 有 univariate cross-correlation function 定義如下:
\[
R_XY(\tau) := E[X_{t+\tau}Y_t]
\]

上述討論我們可知,若一個 WSS 隨機過程 輸入到 LTI 系統,則輸入與輸出為 J-WSS 且其 cross-correlation function 可寫為
\[\begin{array}{*{20}{l}}
{R_{XY}}\left( \tau  \right) = E\left[ {{X_{t + \tau }}{Y_t}} \right] = E\left[ {{X_{t + \tau }}\int_{ - \infty }^\infty  {{X_{t - \gamma }}h} (\gamma )d\gamma } \right]\\
{\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \int_{ - \infty }^\infty  {\underbrace {E[{X_{t + \tau }}{X_{t - \gamma }}]}_{ = {R_X}\left( {\left( {t + \tau } \right) - \left( {t - \gamma } \right)} \right)}h} (\gamma )d\gamma }\\
{\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \int_{ - \infty }^\infty  {{R_X}\left( {\tau  + \gamma } \right)h} (\gamma )d\gamma  \ \ \ \ \ \ (\star)}
\end{array}\] 上式已經非常接近 convolution 形式,但並非 convolution 形式;因為上式積分中的 integrand 為 $R_X(t+\gamma)h(\gamma)$ 但我們要的形式是 $R_X(t - \text{something})h(\text{something})$ ,故我們現在用變數變換:令 $\gamma := - \alpha$ 可得
\[ \Rightarrow {R_{XY}}\left( \tau  \right) = \int_{ - \infty }^\infty  {{R_X}\left( {\tau  - \alpha } \right)h} ( - \alpha )d\alpha \]此即具備 convolution 形式

另外觀察上式 $(\star)$ 與之前計算的 $Y_t$ 之auto-correlation function $(*)$ :
\[{R_Y}\left( \tau  \right) = \int_{ - \infty }^\infty  {\underbrace {\left( {\int_{ - \infty }^\infty  {{R_X}\left( {\tau  - \left( {\beta  - \gamma } \right)} \right)} h(\beta )d\beta } \right)}_{ = {R_{XY}}\left( {\tau  - \beta } \right)}h(\gamma )d\gamma } \]亦即
\[{R_Y}\left( \tau  \right) = \int_{ - \infty }^\infty  {{R_{XY}}\left( {\tau  - \beta } \right)h(\gamma )d\gamma }  \]也就是說 $R_Y$ 為 $h$ 與 $R_{XY}$ 的 convolution。

現在我們將手邊結果總結如下 $(\star \star)$:
\[\left\{ \begin{array}{l}
{R_Y}\left( \tau  \right) = \int_{ - \infty }^\infty  {h(\gamma ){R_{XY}}\left( {\tau  - \beta } \right)d\gamma } \\
{R_{XY}}\left( \tau  \right) = \int_{ - \infty }^\infty  {h( - \alpha ){R_X}\left( {\tau  - \alpha } \right)} d\alpha
\end{array} \right.\]回憶之前我們曾定義對 correlation function 定義 Fourier transform 如下:
\[\begin{array}{l}
{S_X}\left( f \right): = \int_{ - \infty }^\infty  {{R_X}\left( \tau  \right){e^{ - j2\pi f\tau }}d\tau } \\
H\left( f \right): = \int_{ - \infty }^\infty  {h\left( \tau  \right){e^{ - j2\pi f\tau }}d\tau }
\end{array}\]且由 Fourier transform 的 multiplication property:時域convolution = 頻域相乘。故對 $(\star \star ) $ 取 Fourier transform 可得
\[\left\{ \begin{array}{l}
{R_Y}\left( \tau  \right) = \int_{ - \infty }^\infty  {h(\gamma ){R_{XY}}\left( {\tau  - \beta } \right)d\gamma } \\
{R_{XY}}\left( \tau  \right) = \int_{ - \infty }^\infty  {h( - \alpha ){R_X}\left( {\tau  - \alpha } \right)} d\alpha
\end{array} \right.\mathop  \Rightarrow \limits^{F\left\{ . \right\}} \left\{ \begin{array}{l}
{S_Y}\left( f \right) = H(f){S_{XY}}\left( f \right)\\
{S_{XY}}\left( \tau  \right) = {H^*}\left( f \right){S_X}\left( f \right)
\end{array} \right.\] (上述轉換用了一個FACT: $\cal{F}\{ h(-t) \}$ $= H^*(f)$ 我們會在後面給出證明)。

故輸入與輸出的頻域關係為
\[{S_Y}\left( f \right) = H(f){S_{XY}}\left( f \right) = H(f){H^*}\left( f \right){S_X}\left( f \right) = \left| {H\left( f \right)} \right|^2{S_X}\left( f \right)\]亦即
\[
S_Y(f) = |H(f)|^2 S_X(f)
\]

Comments:
1. $S_X(f)$ 稱為 power spectral density。 (nonnegative, real, and even function)
2. $S_{XY}(f)$ 稱為 cross power spectral density。

=================
FACT:
若 LTI 系統脈衝響應 $h(t)$ 之 Fourier transform 為 $H(f)$,則 $h(-t)$之 Fourier transform 為 $H^*(f)$:
=================
Proof
由Fourier transform 定義
\[H\left( f \right) = \int_{ - \infty }^\infty  {h(\tau ){e^{ - j2\pi f\tau }}d\tau } \]現在對兩邊同取 complex conjugate可得
\[{H^*}\left( f \right) = {\left( {\int_{ - \infty }^\infty  {h(\tau ){e^{ - j2\pi f\tau }}d\tau } } \right)^*} = \int_{ - \infty }^\infty  {{h^*}(\tau ){e^{j2\pi f\tau }}d\tau } \]又因為 LTI 系統脈衝響應 $h(t)$ 為 real function 故 $h^*(t) = h(t)$ 我們可改寫上式
\[{H^*}\left( f \right) = \int_{ - \infty }^\infty  {h(\tau ){e^{j2\pi f\tau }}d\tau } \]做變數變換令 $\theta := - \tau$ 可得
\[{H^*}\left( f \right) = \int_{ - \infty }^\infty  {h(\tau ){e^{j2\pi f\tau }}d\tau }  = \int_{ - \infty }^\infty  {h( - \theta ){e^{ - j2\pi f\theta }}d\theta }  = H\left( { - f} \right)\]亦即 $h(-t)$之 Fourier transform 為 $H(-f) = H^*(f)$。 $\square$。

8/22/2013

[系統理論] (弱)平穩隨機過程特性 與 線性非時變系統 (1) -Time domain property

以下將專注在幾類比較特別的隨機過程,首先是嚴格平穩過程(也就是這一類的 隨機過程 可被稱為 "嚴格平穩")

============================
Definition: (Stationary Processes)
嚴格平穩過程 (Strict-sense stationary (SSS) or Strictly stationary process)
1. 隨機過程 $X_t$ 被稱為 nth-order strictly stationary 若下列條件成立:
給定任意 $n$ 個 時間集合 $\{t_1, ..., t_n \}$ 與時間平移量 $\Delta t$, 其對應的  $X_{t1 + \Delta t}, ..., X_{t_n + \Delta t}$ 之 聯合機率分布(joint probabilities) 與其 時間平移量 $\Delta t$ 無關。亦即 對任意 $n$ 維 集合 $B$,
\[
P((X_{t_1 + \Delta t},...,X_{t_n + \Delta t}) \in B)
\]與 $\Delta t$無關。

2. 考慮一隨機過程 $X_t$,若其對任意 正有限整數 $n$ ,都為 n-th order strictly stationary,則我們稱此隨機過程為 strictly stationary。
==============================
Comment:
1.  "平穩 (stationary)" 一詞 概念上表示 對抗 時間平移的能力也就是 Time invariant 的能力。

2. 若隨機過程為 1st-order strictly stationary,則對任意 $t_1$,其對應於 $P(X_{t_1 + \Delta t})$ 與時間平移量 $\Delta t$ 無關,(亦即其機率分布與 $P(X_{t_1})$ 相等 )

3. Strict stationary 是一個非常強的條件。e.g., 如 comment 1,1st-order strictly stationary 由於要求對任意 $t$ 與 $\Delta t$ 機率分布都相同,故期望值必定也要相同,且取函數之後的期望值亦同,亦即我們有以下結果:
對任意函數 $g(x)$,期望值 $E[g(X_{t_1})] = E[g(X_{t_1 + \Delta t})] $

也就是說 對 1st-order strictly stationary 而言,對任意函數 $g(x)$上式期望值都要成立。可想而知,光要求對任意"函數" 就是非常嚴苛的條件,故我們以下將會介紹一個較弱的 stationary 想法 (只對 $E[X_t]$ 與 $E[X_{t_1}X_{t_2}]$ 有要求 (不需要要求對任意函數"$g(x)$")):

===========================
Definition: 廣義平穩 or 弱平穩過程 (Wide-sense stationary (WSS) or weakly stationary process)
隨機過程 $X_t$ 被稱為 Wide-sense stationary (WSS)若下列條件成立:
1. 其 mean function, $E[X_t]$ 與時間 $t$ 無關 (i.e., $E[X_t]$ 為 constant)
2. 其 auto-correlation function, $ E[X_{t_1}X_{t_2}]$ 只與 $t_1 - t_2$ 有關。亦即
\[
E[X_{t_1}X_{t_2}] = R_X(t_1 - t_2)
\]其中 $R_X(\cdot)$ 表 $X_t$ 的 autocorrelation function
============================
Comment:
WSS 過程 $\nRightarrow$ Strictly stationary。( 特例: 若 WSS 過程 為 Gaussian process,則其必為 stirctly stationary。)


現在我們看幾個 平穩過程的例子

---------------
Example 1
令 $Z$ 為任意隨機變數,且對任意時間 $t$,令 $ X_t := Z$。試判斷此 $X_t$ 是否為 strictly stationary。
---------------

Proof
首先觀察 $Z$ 為與時間無關的隨機變數,故由 strictly stationary 定義:給定 $n$ 為有限任意正整數,且令 $B$ 為 n 維 集合,則其 time shfited joint probability
\[P(({X_{{t_1} + \Delta t}},...,{X_{{t_n} + \Delta t}}) \in B) = P((Z,...,Z) \in B) \ \ \ \ (*)
\]上式成立由於 $X_t$ 定義 : (對任意時間 $t$,令 $ X_t := Z$ $\Rightarrow X_{t + \Delta t} = Z$)。觀察 $(*)$ 式可發現其確實與時間平移的量無關。故 $X_t$ 為 strictly stationary。$\square$

---------------
Example 2 (Correlation function is even)
回憶 WSS 隨機過程的 autocorrelation 條件:
\[
E[X_{t_1}X_{t_2}] := R_X(t_1 - t_2)
\]現在我們讓 $t_1 : = t+\tau $ 且 $t_2 = t$,則我們有單變數的autocorelation ( univariate autocorrelation function)
\[
R_X(\tau) := E[X_{t + \tau} X_t]
\]試證此 univariate autocorrelation function $R_X(\tau)$ 為 even function。
---------------

Proof
回憶 even function 的定義: $\forall \tau,$ $R_x(\tau) = x(-\tau)$
故給定 $\tau$,現在觀察
\[{R_X}(\tau ): = E[{X_{t + \tau }}{X_t}] = E[{X_t}{X_{t + \tau }}]
\]上式透過期望值乘法交換性,接著我們由 atuocorrelation function 的定義 可推知上式為
\[E[{X_t}{X_{t + \tau }}] = {R_X}\left( {t - \left( {t + \tau } \right)} \right) = {R_X}\left( { - \tau } \right)\]亦即 univariate autocorrelation function $R_X(\tau)$ 確實為 even function。$\square$

---------------
Example 3 (Delayed WSS process is WSS)
令 $X_t$ 為 WSS 隨機過程,且 mean function $E[X_t] =0$, autocorrelation function 為 $R_X(\tau)$。現考慮一 延遲的隨機過程 $Y_t := X_{t - t_0}$。試證此 $Y_t$ 為 WSS 隨機過程。
---------------

Proof
欲證 $Y_t $ 為WSS 我們必須證明兩個性質:
1. $E[Y_t] $ 與時間 $t$ 無關
2. $E[Y_{t_1}Y_{t_2}]$ 只與 $t_1 - t_2$ 有關

現在我們先證 1. 我們寫下 $E[Y_t] = E[X_{t-t_0}]$。因為 $X_t$ 為 WSS 隨機過程,且 $E[X_t] = E[X_{t + t_0}] =0$。亦即我們有
\[
E[Y_t] =0
\]接著我們證明 2.
\[\begin{array}{l}
E[{Y_{{t_1}}}{Y_{{t_2}}}] = E[{X_{{t_1} + {t_0}}}{X_{{t_2} + {t_0}}}]\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = {R_X}\left( {\left( {{t_1} + {t_0}} \right) - \left( {{t_2} + {t_0}} \right)} \right) = {R_X}\left( {{t_1} - {t_2}} \right) \ \ \  \ (*)
\end{array}\] 亦即,$Y_t$ 為 WSS 隨機過程。

注意到事實上 $E[Y_{t_1}Y_{t_2}] = R_X(\tau)$。因為  $X_t$ 為 WSS 隨機過程且其 autocorrelation function 為 $ R_X(\tau)$,故不失一般性的情況下我們可令 $t_1 := t + \tau$, $t_2 = t$ 則 $(*)$變為
\[\begin{array}{l}
E[{Y_{{t_1}}}{Y_{{t_2}}}] = E[{X_{{t_1} + {t_0}}}{X_{{t_2} + {t_0}}}]\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = E[{X_{t + \tau  + {t_0}}}{X_{t + {t_0}}}]\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = {R_X}\left( {\left( {t + \tau  + {t_0}} \right) - \left( {t + {t_0}} \right)} \right)\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = {R_X}\left( \tau  \right). \ \ \ \  \square
\end{array}\]



WSS processes through LTI systems
現在我們考慮一個 WSS的隨機過程 $X_t$ ,且讓此 $X_t$ 作為 LTI system 的輸入時,我們想知道 輸出 $Y_t $ 會有甚麼特性。

由系統理論我們知道 當 系統為 LTI ,則輸出 $Y_t$ 為 輸入 $X_t$ 之間的關係為 Convolution integral:
\[
{Y_t} = \int_{ - \infty }^\infty  h (t - \tau ){X_\tau }d\tau  = \int_{ - \infty }^\infty  h (\tau ){X_{t - \tau }}d\tau
\]除了 convolution 之外,我們還有下面重要結果:

========================
Theorem: (WSS signal is preserved by LTI system)
考慮輸入 $X_t$ 為 WSS 隨機過程,且考慮穩定的 LTI 系統 具有脈衝響應 $h(t)$ 則輸出 $Y_t$ 為 與 $h(t)$ 做 convolution:
\[
{Y_t} = \int_{ - \infty }^\infty  h (t - \tau ){X_\tau }d\tau  = \int_{ - \infty }^\infty  h (\tau ){X_{t - \tau }}d\tau
\]且輸出 $Y_t$ 仍為 WSS 隨機過程。
========================
Proof
我們現在證明 $Y_t$ 仍為 WSS 隨機過程,亦即要證明
1. mean function, $E[Y_t]$ 與時間無關
2. correlation function, $E[Y_{t_1}Y_{t_2}]$ 只與 $t_1 - t_2$有關

我們先證 mean function 性質
\[E\left[ {{Y_t}} \right] = E\left[ {\int_{ - \infty }^\infty  h (t - \tau ){X_\tau }d\tau } \right] = E\left[ {\int_{ - \infty }^\infty  h (\tau ){X_{t - \tau }}d\tau } \right]\]我們可改寫
\[E\left[ {{Y_t}} \right] = \int_{ - \infty }^\infty  {E\left[ {h(\tau ){X_{t - \tau }}} \right]} d\tau
\]注意到上式我們使用了積分與期望值對調:此法可由下面將積分寫成 Rieman sum 近似而得 (或者讀者可使用 Fuibini's theorem )故
\[E\left[ {{Y_t}} \right] = \int_{ - \infty }^\infty  {E\left[ {h(\tau ){X_{t - \tau }}} \right]} d\tau  = \int_{ - \infty }^\infty  h (\tau )E\left[ {{X_{t - \tau }}} \right]d\tau
\]由於 $X_t$ 為 WSS,故 $E[X_t] $ 與時間無關,亦即存在一常數 $m$ 使得
\[
E[X_t] = m
\]故
\[E\left[ {{Y_t}} \right] = \int_{ - \infty }^\infty  h (\tau )E\left[ {{X_{t - \tau }}} \right]d\tau  = m\int_{ - \infty }^\infty  h (\tau )d\tau \]亦即與時間 $t$ 無關。

現在我們計算 auto-correlation function $E[Y_{t_1}Y_{t_2}]$,用相同近似方法將積分與 期望值互換順序。
\[\begin{array}{l}
E[{Y_{{t_1}}}{Y_{{t_2}}}] = E\left[ {\int_{ - \infty }^\infty  h (\tau ){X_{{t_1} - \tau }}d\tau \int_{ - \infty }^\infty  h (\theta ){X_{{t_2} - \theta }}d\theta } \right]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}&{}
\end{array} = E\left[ {\int_{ - \infty }^\infty  h (\tau )\left( {\int_{ - \infty }^\infty  h (\theta ){X_{{t_1} - \tau }}{X_{{t_2} - \theta }}d\theta } \right)d\tau } \right]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}&{}
\end{array} = \int_{ - \infty }^\infty  h (\tau )\left( {\int_{ - \infty }^\infty  h (\theta )E\left[ {{X_{{t_1} - \tau }}{X_{{t_2} - \theta }}} \right]d\theta } \right)d\tau \\
\begin{array}{*{20}{c}}
{}&{}&{}&{}&{}
\end{array} = \int_{ - \infty }^\infty  h (\tau )\left( {\int_{ - \infty }^\infty  h (\theta ){R_X}\left( {\left( {{t_1} - \tau } \right) - \left( {{t_2} - \theta } \right)} \right)d\theta } \right)d\tau \\
\begin{array}{*{20}{c}}
{}&{}&{}&{}&{}
\end{array} = \int_{ - \infty }^\infty  h (\tau )\left( {\int_{ - \infty }^\infty  h (\theta ){R_X}\left( {\left( {{t_1} - {t_2}} \right) - \left( {\tau  - \theta } \right)} \right)d\theta } \right)d\tau
\end{array}\]上式顯示 $E[Y_{t_1}Y_{t_2}]$ 只與 $t_1 - t_2$ 有關。
故 若 $X_t$ 為 WSS訊號 輸入 LTI 系統,則其輸出 $Y_t$ 亦為 WSS。$\square$

Comment:
回憶對於 確定訊號作為輸入的 LTI 系統,我們有 Fourier transform 作為工具幫助我們分析。現在由於輸入訊號變為隨機過程,但並無對隨機訊號做出定義的 Fourier tranform,不過所幸,注意到若輸入為 WSS 隨機過程,則我們可針對其 correlation function 做 Fourier transform。所得的結果稱為 spectrum analysis。

10/14/2012

[系統理論] 閉迴路系統的暫態響應 與 eigenvalue/eigenvector 關係

首先考慮  $n$階 線性微分方程如下
\[
\dot {\bf x}(t) = G{\bf x}(t); \;\; { \bf x}(0) ={\bf x}_0
\]上式一般表示為 無外力輸入 $u$ 的系統。且我們可對其求解得到
\[{\bf x}(t) = e^{Gt} {\bf x}_0 \ \ \ \ \ \ (*)
\] 我們可進一步將上式的解 $(*)$ 重新用 矩陣 $G$ 的 eigenvalues/eigenvectors 表示; i.e., 若 $G$ 為 $n \times n$ 則 下式 eigenvalue-eigenvector 關係需被滿足
\[G{{\bf{v}}_i} = {\lambda _i}{{\bf{v}}_i}, \text{ for $i=1,2,...,n$}
\]其中 $\lambda_i$ 為 $G$ 的 eigenvalues 且 $v_i$ 為對應的 eigenvectors。 注意到在此我們假設 $G$ 的 eigenvalues 均相異。現在使用這些 eigenvectors, $v_i$,建構  非奇異轉換矩陣 或稱 modal matrix  $M$ 如下
\[M: = \left[ {\begin{array}{*{20}{c}}
{{{\bf{v}}_1}}&{{{\bf{v}}_2}}& \cdots &{{{\bf{v}}_n}}
\end{array}} \right]
\]注意到 $M$ 為 nonsingular 因為 $\{v_i\}$ 彼此線性獨立 (因為相異 eigenvalue 對應 線性獨立的 eigenvector) ;現在使用 $M$,定義下列新狀態 $\bf z$ 轉換
\[
{\bf x} = M {\bf z}
\]將上述新的狀態關系代入原系統 $\dot {\bf x}(t) = G{\bf x}(t)$ 我們可以改寫如下
\[\begin{array}{l}
{\bf{\dot x}}(t) = G{\bf{x}}(t)\\
 \Rightarrow M{\bf{\dot z}}\left( t \right) = GM{\bf{z}}(t)\\
 \Rightarrow {\bf{\dot z}}\left( t \right) = {M^{ - 1}}GM{\bf{z}}(t)
\end{array}
\]其中初始值 (Initial Condition, I.C.) 為 $\begin{array}{l}
{\bf{z}}(0) = {M^{ - 1}}{\bf{x}}(0)\\
\end{array}$
上述轉換又稱對角化,故我們得到 $G$ 如下
\[{M^{ - 1}}GM = \left[ {\begin{array}{*{20}{c}}
{{\lambda _1}}&{}&{}&{}\\
{}&{{\lambda _2}}&{}&{}\\
{}&{}& \ddots &{}\\
{}&{}&{}&{{\lambda _n}}
\end{array}} \right]: = \Lambda \]因此,若我們求解 ${{\bf{\dot z}}\left( t \right) = {M^{ - 1}}GM{\bf{z}}(t)}$ 可得
\[{\bf{z}}(t) = {e^{{M^{ - 1}}GMt}}{\bf{z}}(0) = {e^{\Lambda t}}{\bf{z}}(0)
\]現在將上述結果轉回原本的狀態 $\bf x$
\[{\bf{x}}(t) = M\underbrace {{\bf{z}}(t)}_{ = {e^{\Lambda t}}{\bf{z}}(0)} = M{e^{\Lambda t}}\underbrace {{\bf{z}}(0)}_{ = {M^{ - 1}}{\bf{x}}(0)} = M{e^{\Lambda t}}{M^{ - 1}}{\bf{x}}(0) \ \ \ \ \ (**)
\] 其中
\[\begin{array}{*{20}{l}}
{M{e^{\Lambda t}} = \left[ {\begin{array}{*{20}{c}}
{{{\bf{v}}_1}}&{{{\bf{v}}_2}}&{...}&{{{\bf{v}}_n}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{e^{{\lambda _1}t}}}&0& \cdots &0\\
0&{{e^{{\lambda _2}t}}}&{}& \vdots \\
 \vdots &{}& \ddots &0\\
0& \cdots &0&{{e^{{\lambda _n}t}}}
\end{array}} \right]}\\
{\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = \left[ {\begin{array}{*{20}{c}}
{{e^{{\lambda _1}t}}{{\bf{v}}_1}}&{{e^{{\lambda _2}t}}{{\bf{v}}_2}}&{...}&{{e^{{\lambda _n}t}}{{\bf{v}}_n}}
\end{array}} \right]}
\end{array}
\]另外我們觀察 $M^{-1}$ 定義 $L:= M^{-1}$ 具有 rows 向量為 $l_i$; i.e.,
\[L = {M^{ - 1}} = \left[ {\begin{array}{*{20}{c}}
{{{{{l}}}_1}}\\
{{{{l}}_2}}\\
 \vdots \\
{{{{l}}_n}}
\end{array}} \right]\]
因此我們可以更進一步改寫 $(**)$ 如下
\[\begin{array}{l}
{\bf{x}}(t) = M{e^{\Lambda t}}{M^{ - 1}}{\bf{x}}(0)\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = \left[ {\begin{array}{*{20}{c}}
{{e^{{\lambda _1}t}}{{\bf{v}}_1}}&{{e^{{\lambda _2}t}}{{\bf{v}}_2}}&{...}&{{e^{{\lambda _n}t}}{{\bf{v}}_n}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{l_1}}\\
{{l_1}}\\
 \vdots \\
{{l_n}}
\end{array}} \right]{\bf{x}}(0)
\end{array}
\]或者更簡潔的表示成矩陣的形式
\[{\bf{x}}(t) = \sum\limits_{i = 1}^n {{e^{{\lambda _i}t}}{{\bf{v}}_i}{l_i}} {\bf{x}}(0) \ \ \ \ (\star)
\] 也就是說
\[{\bf x}(t) = e^{Gt} {\bf x}_0 = \sum\limits_{i = 1}^n {{e^{{\lambda _1}t}}{{\bf{v}}_i}{\alpha _i}}\]
現在我們觀察 $(\star)$ 整理結果如下:上述的(自由響應 (free response) )解 與下列三者有關:

  1. ${\bf eigenvalues}$:用以決定 系統自由響應 的 衰減率/增長率 
  2. ${\bf eigenvectors}$: 用以決定 系統自由響應的 " shape "
  3. ${\bf Initial Condition.}$: 用以決定哪一個系統的 mode 參與自由響應的程度

[數學分析] Mean Value Theorem 與 Lipschitz 連續性:有界導數蘊含 Lipschitz

數學分析上的 Lipschitz condition 與 mean value theorem 有一些有趣的關係:在區間上具有均勻有界導數(uniformly bounded derivative)的可導函數 $f$ 必為 Lipschitz 連續。為此首先回憶基本一維的...