制御論第二のバックアップ(No.18)

制御論第二 †

担当：津村幸治准教授
1.5単位
- 数理：限定選択C
- システム：限定選択B
10:15-11:45 工学部六号館 63講義室
教科書
- システム制御理論入門実教出版

↑

内容 †

序論
- 古典制御理論から現代制御理論へ
  - 古典制御理論:(19世紀後半～1930年代)
    マクスウェルによるガバナーの解析。蒸気機関の制御。
    ラウス・フルビッツの安定判別法。
    伝達関数によるシステムの表現。
    周波数領域。ナイキストの安定判別。
    時間領域。PID制御。
    ウィナーフィルタ。
  - 古典制御の利点
    実験的に得やすい入出力特性に基づいて解析・設計する。
    計算量は少ない。
    直感的
  - 古典制御の欠点
    複雑なシステムの設計には向かない。
    システムの内部状態を無視している。
    方法論としての工学。理論よりも実践が優先。
  - シャノンの情報理論(1948年)
    「通信」の一般化
    情報量と通信容量の関係。
    情報量についての普遍的法則を発見。
    モデルや表現法に依存しない一般法則。
  - カルマンの最終目標
    制御に関する純粋理論
    望ましい制御を実現するためにはどのような、そしてどれくらいの情報が必要であるか？
    与えられたプラントを、制御という観点から完全に特徴付ける本質的な特徴は何か？
  - カルマン(1959年～1960年代)
    制御問題の本質を明らかにすること。
    内部状態変数。
    可制御性・可観測性。
    状態推定とフィードバック制御

動的システムと状態方程式
- 動的システム
  - 粘性抵抗とダンパのシステム。~~c dx/dt = u(t)~~dx/dt = u(t)/c = f(u(t))~~ある時点でのシステムの状態は、初期状態と入力の履歴で定まる。~~x = (1/c)∫[0,t]u(τ)dτ + x(0)~~力 u(t) を入力、位置 x(t) を出力とするシステムと見なせる。
  - 動的システム(Dynamical System)~~入力と出力のあるシステムで、初期状態と入力の履歴により、任意の時刻の状態が定まるもの。
  - 静的システム(Static System)~~入力と出力のあるシステムで、その時刻での入力だけからその時刻の出力が定まるもの。~~内部状態のないシステム。
  - 因果性(Causality)~~ある時刻の出力が、それより過去の入力のと初期状態のみによって定まること。~~出力が未来の入力に依存しないこと。~~ある時点から先の入力をカットしても、それ以前の出力はカットしないときと同じになること。~~物理的なシステムは全て因果的。音声フィルタなど、因果的でないシステムも作り得る。

状態方程式
- マス・バネ・ダンパー系。~~運動方程式：{M(d/dt)^2 + c(d/dt) + k}x(t) = u(t)~~出力方程式：y(t) = x(t)~~内部状態の次元を増やして、運動方程式を一階微分にする。~~ x = (x[1],x[2]) = (x,dx/dt)~~ (d/dt)x[1] = x[2]~~ (d/dt)x[2] = (u-cx[2]-kx[1])/M~~y = x[1]~~と変形すると。~~x = f(x,u)~~ y=g(x,u) の形で表せる。
- システムの状態方程式表示~~入力：u(t)~~出力：y(t)~~状態変数：x(t)を使って、~~状態方程式：dx/dt = f(x,u,t)~~出力方程式：y = g(x,u,t)の形でシステムを記述する。~~入力・出力・状態変数は多次元になりうる。~~入出力が共に一次元のシステムを、一入出力システム(SISO System)~~入出力が多次元のシステムを、多入出力システム(MIMO System)と呼ぶ。

線形時不変システム
- 十分長い時間、入出力が共に0となっているとする。~~動的システムP：y(t)=P[u(t)]について、~~任意の入出力関係 y_1 = P[u_1], y_2 = P[u_2] の間に、~~線形性：αy_1 + βy_2 = P[αu_1+βu_2]が成り立つとき、~~Pは線形であるという。~~状態方程式・出力方程式が線形。
- システムが任意の時間シフトについて~~y(t+τ) = P[u(t+τ)]となるとき、~~Pは時不変であるという。~~状態方程式・出力方程式に陽にtを含まない。
- 例）dx/dt = ax + u, y = x~~y = ∫[-∞,t]exp(a(t-τ))u(τ)dτ は線形時不変システム。~~t→t+t'と時間をずらすと、~~y(t+t') = ∫[-∞,t+t']exp(a(t+t'-τ))u(τ)dτ = ∫[-∞,t]exp(a(t-τ))u(τ+t')dτ
- システムが y(t) = ∫[-∞,t]p(t-τ)u(τ)dτ と、たたみこみの形で書かれていれば線形時不変システム。

伝達関数
- t<0でx(t)=0、∫[0,∞]|x(t)|e^(-at)dt＜∞となる滑らかな関数はラプラス変換でき、~~x(t)→X(s)=∫[0,∞]x(t)exp(-st)dt と表す。~~逆変換も存在し、x(t)=(1/(2πi))∫[c-i∞,c+i∞]X(s)exp(st)ds
- 二つの関数のたたみこみを変換すると、変換後にはただの乗算になっている。~~線形時不変システムはたたみこみで表せるので、ラプラス変換すればシステム全体を関数一つで表現できる。~~これが伝達関数で、二年冬や三年夏に扱ったもの。

線形時不変システム←[ラプラス変換]→伝達関数線形時変システム←[？？]→？非線形システム←[？？]→？線形時不変システム以外は伝達関数では扱えない。状態空間表現や元々の微分方程式を直接扱う必要がある。

状態空間法
- システムに内部状態を表す変数を導入して、~~内部状態についての連立一階微分方程式と、内部状態から出力を導く方程式で表現する。
- 伝達関数での表現と比べて広い範囲のシステムを表現できる。~~システムが線形ならば状態方程式・観測方程式は線形にでき、時不変ならば方程式に陽に時間変数が現れないようにできる。
- x : 状態ベクトル~~dx/dt = f(x,t,u) ：状態方程式~~y = g(x,t,u)：観測方程式。出力方程式。~~状態方程式と観測方程式を合わせてシステム方程式と呼ぶ。~~講義では主に線形時不変システムを扱う。

線形システムの解
- dx/dt = Ax + Bu~~y = Cx + Du~~Aはtに依存してもよいとする。
- まず同次方程式を解く。~~dx/dt = Ax
- 標準基底e[i]を使って、初期値e[i]のときの解をφ[i]、それを横に並べたものをΦとすると、~~dφ[i]/dt = Aφ[i]~~dΦ/dt = AΦ~~Φ(0) = Iが成立。~~一般解はx(t) = Φ(t)x(0) で表せる。
- 元の方程式の解は、x(t) = Φ(t)x(0) + ∫[0,t]Φ(t-τ)B(τ)u(τ)dτ

線形時不変システムの解
- dx/dt = Ax + Bu~~y = Cx + Du~~Aはtに依存しない定行列
- (d/dt - A)x = exp(At)(d/dt)exp(-At)x(t) = Bu(t)~~x(t) = exp(At)x(0) + ∫[0,t]exp(A(t-τ))Bu(τ)dτ
- 伝達関数表現では、X(s) =(sI-A)^(-1){x(0)+BU(s)}~~(sI-A)^(-1)にあたるのがexp(At)、(たたみこみが積に変わっている)

exp(At)？
- (d/dt)exp(At) = Aexp(At)~~Aの固有ベクトル Au = λu への作用させると、exp(At)u = exp(λt)u~~固有ベクトルの方向にexp(λt)だけ伸ばす。
- detP≠0となる行列について、~~P exp(At) P^(-1) = exp(PAP^(-i)t) ：合同変換は中に素通り。
- AとBが可換ならば、exp(At)exp(Bt)＝exp((A+B)t)が成立。非可換ならばそうなるとは限らない。

行列 $A$ をジョルダン標準型に $TAT^{-1} = J$ に標準化すると、

$T \exp(At)T^{-1} = \exp(Jt)$
（中略）

$L[\exp(Jt)] = (sI-J)^{-1}$

$L[\exp(At)] = L[T^{-1}\exp(Jt)T] = T^{-1}(sI-J)^{-1} T = (sI-A)^{-1}$
状態空間表現での $\exp(At)$ をラプラス変換すると $(sI-A)^{-1}$
線形時不変システムの状態空間表現と伝達関数表現の対応。

安定性
- 平衡点
  状態空間表現でのシステム $(d/dt)x = f(x)$ について、
  $f(x[e]) = 0$ となる点を平衡点と呼ぶ。
- 平衡点まわりの安定性
  任意の正数 $\varepsilon$ に対して、 $x[e]$ からの距離が $\delta$ 未満である初期値から時間を進めたときに、
  それより先の任意の時点で $x[e]$ からの距離が $\varepsilon$ 未満に留まるような正数 $\delta$ が存在するなら、その系は安定。
  さらに $t\to\infty$ で $x[e]$ に収束するなら、その系は漸近安定であるという。
- 線形時不変システム
  $x = \exp(At)x[0]$ について、
  ジョルダン標準系に直して、 $Tx = T\exp(At)T^{-1}Tx[0] = \exp(Jt)Tx[0]$
  $x$ の安定性と $Tx$ の安定性は一致する。
  $J$ の固有値の実部が全て負ならば、 $x,Tx$ は原点に漸近安定。
  実部が正の固有値があれば、不安定。
  固有値の実部が全て正ではなく、実部が0の固有値がサイズ２以上のジョルダンブロックを成すなら不安定。
  それ以外なら安定。
- リアプノフの定理。
  $dx/dt = Ax$ について、
  $XA+A^T X+Q = O$ となる正定値行列 $Q$ と対称行列 $X$ について、
  $A$ が漸近安定であることと $X$ が正定値であることは同値。
  $X$ が一意であることは講義では証明略。
  $X$ として $X = \int_0^{\infty}\exp(A^T t)Q\exp(At)dt$ を選ぶと、正値対称。
  $[\exp(A^T t)Q\exp(At)]_0^{\infty} = \int_0^{\infty}(d/dt)\exp(A^T t)Q\exp(At)dt = A^T X+XA$
  $\exp(At)\to 0 (t\to \infty)$ ならば、左辺は $-Q$ になるから、 $X$ は $XA+A^T X+Q = O$ の正値対称な解になっている。
  もし $A$ が、実部が非負の固有値を持つなら、その固有値・固有ベクトルを $\lambda ,x$ として、
  $XA+A^T X+Q = 0$ から、
  $Re \lambda*(x,Xx) + (x,Qx) = 0$ だが、各項が正なので、
  この式は不成立。よって $A$ には実部が非負の固有値は存在しない。
  - $A^T$ は $A$ の転置

動的システムの構造
- 相似変換
  
  $dx/dt = Ax + Bu$
  
  $y = Cx + Du$
  について、
  $z = Tx$ と変数変換すると、各係数は
  
  $dz/dt = TAT^{-1}z + TBu = Ez + Fu$
  
  $y = CT^{-1}z + Du = Gz + Hu$
  と変換される。
  どちらも同じ入力に対して同じ出力を返す、同じシステムの別な表現になっている。
  入出力関係を変えない変換を、相似変換と呼ぶ。
  状態空間表現には、相似変換の分だけ任意性が残っている。
- モード分解
  状態空間表現の任意性を使って、方程式を理解しやすい形に書き直す。
  $A$ をジョルダン標準形に直すと、状態方程式をジョルダンブロックごとに
  独立した小さな微分方程式に分解できる。
  分解された各方程式は、各固有値に対応するモードと呼ばれる。

可制御性(Controllability)
- 線形時不変システム $\dot{x} = Ax + Bu$ について、
  任意の初期値に対して、適切な入力をとれば有限時間内に、
  $x$ を原点に到達させることができるとき、
  このシステムは可制御であるという。
- この微分方程式の解は、 $x(t) = e^{At}x_0 + \int_0^t e^{A(t-\tau)}Bu(\tau)d\tau$
  時刻 $t=T$ で原点に到達したとすると、
  $-x_0 = \int_0^T e^{-A\tau}Bu(\tau)d\tau$ が成立する。
  この $u$ に対して、初期値 $x_2 = e^{-AT}x_1+x_0$ で始めることを考えれば、
  可制御性とは、任意の入力に対して任意の位置に持っていけるかどうかを表すとも考えられる。
- 定理
  線形時不変システムが可制御であるための必要十分条件は、
  可制御性行列 $U_c = [B,AB,A^2 B,...,A^{n-1}B]$ のランクがnであることである。
  nは状態ベクトルの次元。
- 証明
  - 可制御であるためには $\mathrm{rank} U_c = n$ が必要。
    ケーリー・ハミルトンの定理を使えば、任意のn次正方行列の多項式は、(n-1)次以下の多項式で表せるので、
    $e^{At}$ の各項を(n-1)次以下の多項式で書けば、
    $e^{At} = \sum_{k=1}^n q_{k}(\tau)A^{k-1}$ と表せる。
    これを使って解を変形すると、
    $e^{-At}x(t) - x_0 = \int_0^t e^{(-A\tau)}Bu(\tau)d\tau$
    $= \int_0^t \sum_{k=1}^n q_{k}(-\tau)A^{k-1}Bu(\tau)d\tau$
    $= \sum_{k=1}^n A^{k-1}B \int_0^t q_{k}(-\tau)u(\tau)d\tau$
    $= \sum_{k=1}^n A^{k-1}B h_{k}(t)$
    $= U_c h$
    左辺は任意に動かしうる。任意の左辺に対して、等式を成り立たせる $U_c$ が存在するためには
    $\mathrm{rank} U_c = n$ が必要。
  - $\mathrm{rank} U_c = n$ ならば可制御。
    $W_s = \int_0^s e^{-A\tau}BB^{T}e^{-A^{T}\tau}d\tau$
    入力 $u$ を $u=B^{T}e^{-A^{T}\tau}W_s^{-1}\{-x(0)+e^{-As}x(s)\}$
    と取ってやれば、
    目標の $x_s$ を実現する入力になっている。
    途中で使っている $W_s^{-1}$ が存在すれば可制御。
    もし、 $\det W_s = 0$ となる $s\gt 0$ が存在するなら、
    $W_s a = 0$ となる非零ベクトルaが存在する。
    $a^{T}W_s a = \int_0^s ||B^{T}e^{-A^{T}\tau}a||d\tau = 0$
    ノルムの値は非負なので、等式が成り立つためには、被積分関数は常に零である必要がある。
    $B^{T}e^{-A^{T}\tau}a = 0$
    この式を $\tau$ で $k$ 階微分して $\tau = 0$ とすれば、
    $B^{T}(A^{T})^{k}a = 0$
    $U_{c}^{T}a = 0$
    $\mathrm{rank} U_c = n$ ならば、このようなaは存在しない、
    ゆえに、 $\det W_s$ は常に非零で、このシステムは可制御。

演習問題１
- レポート

可観測性
- ある有限の時間 $[0,s]$ の入出力の測定から、
  内部状態の初期値 $x(0)$ を唯一に決定できるとき、
  そのシステムは可観測(observable)という。
- 定理
  線形時不変システムが可観測であるための必要十分条件は、
  可観測性行列 $U_o=[C,CA,CA^2,...,CA^{n-1}]^{T}$ ( $CA^{k}$ を縦に並べた行列)の
  ランクがnであること。
  証明はプリントで。

正準分解と最小実現
- 線形時不変システム,主システム
  $x \in \mathbb{R}^{n}, u \in \mathbb{R}^{m}, y \in \mathbb{R}^{l}$
  $\dot{x} = Ax + Bu$
  $y = Cx$ に対して、
  双対システム
  $\dot{x} = A^{T}x + C^{T}u$
  $y = B^{T}x$ を考えると、
  主システムの可制御性と双対システムの可観測性は同値。
  主システムの可観測性と双対システムの可制御性は同値。
  片方の可制御性行列は、もう片方の可観測性行列の転置。
- $W = \{ x^d \in \mathbb{R}^n | x^d = U_c h, h \in \mathbb{R}^{mn} \}$ は「可制御な部分空間」
  可制御系については $W=R^n$
  $AW = \{ Ax^d | x^d \in W \}$ は、ケーリー・ハミルトンの定理で $U_c$ の中の $A^n$ の次数下げをすれば、
  $AW \subset W$ が成立。可制御な空間からの時間発展は可制御。
  $M = \{ x^c \in \mathbb{R}^n | x^c = \tilde{U}_c h, h \in \mathbb{R}^{ml}\}$ を定義すると、
  $M$ は可観測な部分空間。
  $A^{T}M \subset M$ も成立。可観測な空間の時間発展は可観測
  $M$ の直交補空間 $M^{\bot}$ を考えると、
  $AM^{\bot} \subset M^{\bot}$
- 以上の集合を使って、以下の集合を定義すれば、
  $R_a = W \cap M^{\bot}$ 可制御で、可観測でない
  $R_b = W \cap M$ 可制御で可観測
  $R_c = W^{\bot} \cap M^{\bot}$ 可制御でなく、可観測でもない
  $R_d = W^{\bot} \cap M$ 可制御でなく、可観測
  $\mathbb{R}^n = R_a \oplus R_b \oplus R_c \oplus R_d$ と直和分解できる。
- それぞれの部分空間の時間発展を考えると、
  $AR_a \subset R_a$
  $AR_b \subset R_a \oplus R_b$
  $AR_c \subset R_a \oplus R_c$
  $AR_d \subset \mathbb{R}^n$ となる。
  それぞれの空間の基底をとって、
  $R_a = \langle x^a_i \rangle_{i=1}^{\nu_a}$
  $R_b = \langle x^b_i \rangle_{i=1}^{\nu_b}$
  $R_c = \langle x^c_i \rangle_{i=1}^{\nu_c}$
  $R_d = \langle x^d_i \rangle_{i=1}^{\nu_d}$
  $\nu_a + \nu_b + \nu_c + \nu_d = n$
  基底を初期値に時間発展を考えれば、
  $Ax^a_i = \sum_{j=1}^{\nu_a}f^{aa}_{ji}x^a_j$
  $Ax^b_i = \sum_{j=1}^{\nu_a}f^{ab}_{ji}x^a_j + \sum_{j=1}^{\nu_b}f^{bb}_{ji}x^b_j$
  $Ax^c_i = \sum_{j=1}^{\nu_a}f^{ac}_{ji}x^a_j + \sum_{j=1}^{\nu_c}f^{cc}_{ji}x^c_j$
  $Ax^d_i = \sum_{j=1}^{\nu_a}f^{ad}_{ji}x^a_j + \sum_{j=1}^{\nu_b}f^{bd}_{ji}x^b_j + \sum_{j=1}^{\nu_c}f^{cd}_{ji}x^c_j + \sum_{j=1}^{\nu_d}f^{dd}_{ji}x^d_j$
- 各部分空間の基底を使って、
  $B = [x^{a}_1,x^{a}_2,...,x^{a}_{\nu_a}]G^a + [x^{b}_1,x^{b}_2,...,x^{b}_{\nu_b}]G^b$
  各部分空間の基底を横に並べた行列を $T$ として変換をすると、
  $AT=TF,B=TG,CT=H$
  $F=\left[\begin{matrix}F^{aa} & F^{ab}& F^{ac} & F^{ad}\\0 & F^{bb} & 0 & F^{bd}\\0 & 0 & F^{cc} & F^{cd}\\0 & 0 & 0 & F^{dd}\end{matrix}\right]$
  (by yambi)
- 以上の行列を使って、
  $\xi = T^{-1}x$ と変換すれば。
  $\dot{\xi} = T^{-1}AT\xi + T^{-1}Bu = F\xi + Gu$
  $y = CT\xi = H\xi$ となる。これを正準分解という。
  これをブロック線図で表すと、入力は可制御な空間である $R_a,R_b$ へ
  出力は可観測な空間である $R_d,R_b$ から繋がり、
  $R_b,R_c$ は独立。 $R_a$ は他全てから影響を受け、
  $R_d$ は他全てに影響を及ぼす。

最小実現
- 正準分解した上で入力から出力への伝達関数を書いてみると、
  伝達関数が表しているのは、可制御かつ可観測な部分空間の成分だけであることがわかる。
- 可制御かつ可観測な部分空間の成分を取り出しても、入出力関係は変わらない。
  これをこのシステムの最小実現と呼ぶ。

MIMOでの零点
- システムを表す行列を[[(sI-A), B],[C,D]]と並べたときに、
  この行列をランク落ちさせるsをこのシステムの不変零点と呼ぶ。
  この定義はSISOでの零点の定義の一般化になっている。

SISOの可制御正準形
- ケーリーハミルトンの定理を使って、Aのかかった連立常微分方程式を、
  $x_n$ についての高階の微分方程式に帰着。
  詳細は書ききれないのでググって見つけてください。http://green.cc.tsukuba.ac.jp/tetsu/shibaura/ac.htmlとか。
- $\phi(A)=O$
  $A x = \dot{x} - Bu$ $A^2 x = \ddot{x} - AB\dot{u} - Bu$ $A^3 x = \dddot{x} - B\ddot{u} - AB\dot{u} - A^2 Bu$

SISOの可観測正準形

極配置問題
- SISOの可制御系 $\dot{x} = Ax + Bu$ について、
  内部状態 $x$ を観測し、 $u=Fx$ とフィードバックすることを考える。
  このようなフィードバックを『状態フィードバック』と呼ぶ。
  $\dot{x} = Ax + BFx = (A+BF)x$
  可制御標準形に変換して、 $\dot{z} = \tilde{A}z + \tilde{B}u = (\tilde{A} + \tilde{B}\tilde{F})z$
  $\tilde{B}$ は最後の行が１で、それ以外がすべて０の縦ベクトルなので、
  $\tilde{F}$ の項は第n行にしか効いてこない。
  $\tilde{A}$ の第n行は、Aの特性多項式の係数。
  $F$ を掛けてフィードバックすることは、特性多項式の係数をずらすことになる。
  特定の場所に伝達関数の極を配置したいときには、その極を根に持つ多項式を作り、
  $\tilde{F}$ の値を調節して、特性多項式の係数をずらして、その多項式を作ってやれば良い。
  こうして、可制御系については閉ループ系の極は任意の位置に配置でき、常に安定化できる。
  $F = \tilde{F}T^{-1}$

レポート問題、12/10まで。
- コンデンサと抵抗を並列に繋げた回路と、コイルと抵抗を並列に繋げた回路を直列に繋いだ回路について、
  回路全体に流れる電流を入力、出力電圧を出力、
  コンデンサにかかる電圧と、コイルに流れる電流を内部状態として、
  可制御性行列と伝達関数を調べ、可制御性との関係を述べよ。
- 長さ・質量の異なる振り子について、振り子の振れ幅とその速度についての可制御性と
  振り子の長さ・重りの質量の関係を調べよ。振り子の振れ幅は微小としてよい。
- SISO可制御正準系の $\tilde{A}$ について、 $|sI-\tilde{A}| = s^n + \sum_{k_1}^{n}a[k]s^{k+1}$ を計算せよ。
- $\dot{x} = [[1,2],[-3,-4]]x + [2,3]u$ について、
  状態フィードバックにより、極を $-2,-3$ に配置したい。
  - 直接的にFを求めよ。
  - 可制御正準系に直してからFを求めよ。

状態推定
- オブザーバ
  $\dot{x} = Ax + Bu$
  $y = Cx$ について、入出力信号を観測して、 $x$ を推定することを考える。
  - $x$ と同じ方程式に従う別な系 $\hat{x}$ を用意すると、
    誤差 $e = \hat{x}-x$ は $\dot{e} = Ae$ に従い、
    $A$ が漸近安定ならば誤差は0へ漸近し、不安定ならば発散してしまう。
  - 代わりに $\hat{x}$ として、 $\dot{\hat{x}} = A \hat{x} Bu + L(c\hat{x}-y)$ を使うことを次次次回考える。

今日は原さんの講義

最適レギュレータ問題
アクチュエータに使うエネルギーと目標地点への制御の重み付き最適化。

可制御系は適当な係数を掛けた状態フィードバックで安定化できる。
$\dot{x}\ =\ Ax\ +\ Bu$
$y\ =\ Cx$
$u\ =\ Fx$ 状態フィードバック
- このシステムの内部状態や出力の変化を明示すれば
  $x\ =\ e^{(A+BF)t}x_0$
  $y\ =\ Ce^{(A+BF)t}x_0$
  $F$ をうまく選んで、任意の初期値に対して
  $\mathrm{lim}_{t\to \infty}x(t)\ =\ 0$ となるようにする。
  システムが可制御であれば、 $A+BF$ の固有値は任意に設定可能。
  では具体的に $F$ はどう決めるのがよいか。
- 評価関数(Performance Index)
  $J\ =\ \int_{0}^{\infty}\left{x^{T}Qx + u^{T}Ru \right}dt$ という評価関数を作って、これを最小化する。
  $Q,R$ は状態，制御入力に対する重み行列。
  $Q$ は半正値対称、 $R$ は正値対称に選ぶ。
  これは $Q = C^{T}C$ と選ぶと $y$ の値で評価したことになる。
  $J\ =\ \int_{0}^{\infty}x_{0}^{T}e^{A^{T}+F^{T}B^{T}}\left{Q + F^{T}RF^{T}\right}e^{(A+BF)t}x_{0}dt$ を最小とする $F_{opt}$ は、
  $F_{\mathrm{opt}} = - R^{-1}BP$ の形で与えられる。
  $P$ はRiccati方程式 $A^{T}P+PA+Q-PBR^{-1}B^{T}P = O$ の解で、正値対称行列。
  可制御で可観測なシステムについては、解は常に存在し、 $A-BF_{\mathrm{opt}}$ は安定になる。
  - $x$ の安定性を仮定すると、任意の正値対称行列 $P$ について、
    $\int_{0}^{\infty}\frac{d}{dt}\left{x^{T}Px\right}dt = -x_{0}Px_{0} = \int_{0}^{\infty}x^{T}\left{(A^{T}+F^{T}B^{T})P + P(A+BF)\right}x$ から、
    $J = x_{0}^{T}Px_{0} + \int_{0}^{\infty}x^{T}\left{(F^{T}+PB^{T}R^{-1})R(F+R^{-1}B^{T}P) + PB^{T}R^{-1}B^{T}P+A^{T}P + PA + Q \right}x$ が成立。
    $F$ を変化させて $J$ を最小化するには、 $F=-R^{-1}B^{T}P$ を選べばよい。

Riccati方程式の解法
- ハミルトン行列 $H=[[A,-BR^{-1}B^{T}],[-Q,-A^T]]$ を定義し、固有値を求める。
  固有値は実軸・虚軸に対して対称に現れる。
  $(A,B)$ が可制御で、 $(A,\sqrt{Q})$ が可観測ならば虚軸上には固有値は無い。
  $H$ の安定な固有値・固有ベクトルn個を並べたものを $\Lambda,[U,V]^{T}$ とすると、
  $[UV^{-1},-I]H[U,V]^{T} = [UV^{-1},-I][U,V]^{T}\Lambda = 0$ が成立。
  $P$ として $P=UV^{-1}$ をとると、PはRiccati方程式の解になっている。
  これでRiccati方程式を2n次正方行列であるハミルトン行列Hの固有値問題に帰着できた。

Riccati方程式の解について、
$F=-R^{-1}B^{T}P$
$A^{T}P+PA+Q-PBR^{-1}B^{T}P=0$
- $Q=P(j\omega I-A)+(-j\omega I-A^{T})P+F^{T}RF$ と変形し、
  $|I-F(j\omega-A)^{-1}B|_R$ を考えると、
  $||I-F(j\omega-A)^{-1}B||_R = B^{T}(-j \omega I-A^{T}^{-1}Q(j\omega I-A)^{-1})B+R \geq R$ が任意の $\omega$ について成立。
  $||\cdot||_R$ は $R$ を計量とする内積。
  $I-F(j\omega-A)^{-1}B$ は、最適レギュレータをかけたときの一巡伝達関数。
  この不等式が意味するところは、任意の周波数の入力に対し、
  外乱に対する強度を示す感度関数の値を１より小さく保てること。

制御入力に対する重みを無限大にして、最適制御をかけたときの閉ループ伝達関数は、
開ループ伝達関数の不安定極を虚軸対称に反転させたもの。
制御対象に対する重みを無限大にして、最適制御をかけたときの閉ループ伝達関数は、
不安定極を虚軸対称に反転させ、他の極は∞に飛ばした形になる。

状態推定
- オブザーバ
  $\dot{x} = Ax + Bu$
  $y = Cx$ の入出力を観測してxを推定する。
- $\dot{\hat{x}} = A \hat{x}+Bu + L(C\hat{x}-y)$ と変数を設定する。
  これを恒等オブザーバと呼ぶ。
  誤差 $e = \dot{\hat{x}}-x$ の挙動を見ると、
  $\dot{e} = Ae + LC\hat{x} - LCx = (A+LC)e$
  系が可制御ならば、 $(A^{T},C^{T})$ が可制御で、
  $(A+LC)^{T}$ の極は任意に配置可能、 $A+LC$ は安定化可能になる。
  $A+LC$ を安定化すれば $\mathrm{lim}_{t\to \infty}e = 0$
  $\hat{x}$ によって入出力だけから $x$ を観測できる。

最小次元オブザーバ
$(A,C)$ が可観測で、 $C\in \mathbb{R}^{m\times n}\ \mathrm{rank}C = m$ だとする。
$C$ に行を追加して、正則な行列 $T$ を作り、
$T=\left[\begin{matrix} C \\ U \end{matrix}\right]$
$\tilde{A}=TAT^{-1}$ とする。 $C=[I\ 0]T$ が成立。 $\tilde{C}=CT^{-1}=[I\ 0]$ とすると、
$(A,C)\to (\tilde{A},\tilde{C})$ は $T$ による相似変換。
$\tilde{A} = \left[\begin{matrix}A_{11}&A_{12}\\A_{21}&A_{22} \end{matrix}\right]$ とすると、元の系が可観測であることから
$\mathrm{rank}\left[\begin{matrix}A_{12}\\A_{12}A_{22}\\ \vdots \\ A_{12}A_{22}^{n-1}\end{matrix}\right] = n-m$ 、
$A_{22},A_{12}$ が可観測であることが得られる。
$A_{22}+LA_{12}$ を安定化するLが存在する。
$\dot{z} = (A_{22}+LA_{12})(z-Ly)+(A_{21}+LA_{11})y + (U+LC)Bu$ とすると、
$e=(U+LC)x-z$ は、
$\dot{e} = \left[\begin{matrix}L&I\end{matrix}\right]\left[\begin{matrix}C\\U\end{matrix}\right](Ax+Bu)-\dot{z}=(A_{22}+LA_{12})e$ に従う。 $A_{22}+LA_{12}$ を安定化すれば、誤差 $e$ は0に収束する。
$\left[\begin{matrix}y\\z\end{matrix}\right]\to \left[\begin{matrix}I&0\\L&I\end{matrix}\right]\left[\begin{matrix}C\\U\end{matrix}\right]x$
$\left[\begin{matrix}C\\U\end{matrix}\right]^{-1}\left[\begin{matrix}I&0\\-L&I\end{matrix}\right]\left[\begin{matrix}y\\z\end{matrix}\right]\to x$

カルマンフィルタ
フィルタ：推定対象に混入する外乱を除去するもの。
$\dot{x} = Ax+Bu+v$
$y=Cx+w$ という系を考える。
$v$ はシステム雑音、 $w$ を観測雑音と呼ぶ。
$v,w$ は白色雑音だと仮定する。
$E[v]=E[w]=0$
$E[v(t)v^{T}(t+\tau)] = Q\delta(\tau)$
$E[w(t)w^{T}(t+\tau)] = R\delta(\tau)$
$E[v(t)w^{T}(\tau)] = 0$
入出力 $u,y$ から状態 $x$ を推定する。
$\tilde{x}(t) = E[x(t) | y(\tau),u(\tau),\tau \leq t]$
- 条件付き期待値
  $f(x,y)$ を $x,y$ の同時確率密度だとして、
  $f(x|y) = \frac{f(x,y)}{\int f(x,y)dx}$ と、条件付き確率密度を定義して、
  条件付き期待値を $E[x|y]=\int xf(x|y)dx$ と定義する。
  何らかの方法で $y$ から $x$ の推定値 $\tilde{x}(y)$ を計算できたとしたとき、
  その二乗誤差の期待値には下限が存在し、 $E[(x-\tilde{x})^2]\geq E[(x-E[x|y])^2]$ となる。
  - 補題：任意の関数 $h(x)$ について、 $E[(x-E[x|y])h(y)] = 0$
- 条件付き期待値 $\tilde{x}$ は、
  $\dot{\tilde{x}} = F\tilde{x}+Gy+Hu$ という形の微分方程式に従う。
  係数行列 $F,G,H$ が満たすべき性質は、
  $E[\tilde{x}] = E[x]$
  $E[\dot{x}] = AE[x]+Bu$
  $E[y] = CE[x]$
  $E[\dot{\tilde{x}}] = E[\dot{x}] = (F+GC)E[x]+Hu$ から
  $A=F+GC,\ B=H$
  $\dot{\tilde{x}} = (A-GC)\tilde{x}+Gy+Bu$
  
  残ったGを決めるために、誤差 $e=x-\tilde{x}$ を考える。
  $\dot{e} = (A-GC)e - Gw + v$
  誤差の共分散行列を $\tilde{P} = E[e(t)e^{T}(t)]$ とすると、
  十分時間が経った後では、 $\tilde{P}$ は
  $(A-GC)\tilde{P}+\tilde{P}(A-GC)^{T}+GRG^{T}+Q=0$ の解になっている。
  変形すると $A\tilde{P}+\tilde{P}A^{T}+(G-\tilde{P}C^{T}R^{-1})R(G-\tilde{P}C^{T}R^{-1})^{T}+Q-\tilde{P}C^{T}R^{-1}C\tilde{P}=0$
  $A-GC$ が安定ならば、 $\tilde{P}$ には下限 $P$ が存在し、
  $AP+PA^{T}+Q-PC^{T}R^{-1}CP=0$ に対して
  $\tilde{P}\geq P$ が成立。
  $G=PCR^{-1}$ と選んでやると $P=\tilde{P}$ となり、
  誤差の分散を最小にする $G$ が得られた。

制御論第二 のバックアップ(No.18)

制御論第二 †

内容 †

制御論第二のバックアップ(No.18)