保護中: 例題

このコンテンツはパスワードで保護されています。閲覧するには以下にパスワードを入力してください。

What’s new? (2024/10/28)

2025.4.7 『RCPSP法による造船工程計画の実践』のサポートページの準備を始めました。しばらくはパスワード入力が必要です。

2024.10.28 『ゼミナール 制御技術入門』の修正版ver.1.1が近々刊行されます。

2024.9.5 『ゼミナール 制御技術入門』のフォローアップサイトを、大まかに仕上げました。また、正誤表に1件追加しました。

2024.5.6 正誤表を作成しました。
正誤表

2024.5.2 『ゼミナール 制御技術入門』のフォローアップサイトの整備が遅れて申し訳ありません。ようやく5章まで要点を中心に仕上げました。また、古いプログラムをアップしていたので、一部修正し更新しました。

2024.1.9 次がXに出ました。


2023.12.30 本サイトの全面改訂を始めます。

問・演習問題の解答(1章)

1. 物理法則から状態方程式を導く


問1.1 十分大きなkに対して,v(1)=0.3679v(2)=0.0498v(3)=0.0067

問1.2
#1:ドアは急に衝突音をともなって閉まる。
#2:ドアは速やかに静かに閉まる。
#3:ドアはなかなか閉まらない。

問1.3 答は図1.12。

問1.4 座標変換
\left[\begin{array}{c} x_a(t) \\ x(t) \end{array}\right] = \left[\begin{array}{cc} 0 & I_2 \\ 1 & 0 \end{array}\right] \left[\begin{array}{c} x(t) \\ x_a(t) \end{array}\right]
を行うと
\left[\begin{array}{c} \dot{x}_a(t) \\ \dot{x}(t) \end{array}\right]= \left[\begin{array}{c|cc} -\frac{1}{T_a} & 0 & 0 \\\hline 0 & 0 & 1\\ K\omega_n^2c_a & -\omega_n^2 & -2\zeta\omega_n \end{array}\right] \left[\begin{array}{c} x_a(t) \\ x(t) \end{array}\right]+ \left[\begin{array}{cc} \frac{K_a}{T_a} \\\hline 0 \\ 0 \end{array}\right] u(t)

y(t)= \left[\begin{array}{c|cc} 0 & 1 & 0 \end{array}\right] \left[\begin{array}{c} x_a(t) \\ x(t) \end{array}\right]

【1】 回路方程式はつぎのように整理できる。

\begin{array}{ll} L\dot{i}(t)=-(\frac{R_1R_2}{R_1+R_2}+\frac{R_3R_4}{R_3+R_4})i(t)-(\frac{R_1}{R_1+R_2}-\frac{R_3}{R_3+R_4})v(t) +u(t) \\[1mm] C\dot{v}(t)=-(\frac{R_2}{R_1+R_2}-\frac{R_4}{R_3+R_4})i(t) -(\frac{1}{R_1+R_2}+\frac{1}{R_3+R_4})v(t)\\[1mm] \end{array}

R_1R_4=R_2R_3が成り立つとき

\begin{array}{ll} L\dot{i}(t)=-(\frac{R_1R_2}{R_1+R_2}+\frac{R_3R_4}{R_3+R_4})i(t)+u(t) \\ C\dot{v}(t)=-(\frac{1}{R_1+R_2}+\frac{1}{R_3+R_4})v(t) \\ \end{array}

ブロック線図は図3.4参照。

【2】 v_1(t)=\dot{x}_1(t)v_2(t)=\dot{x}_2(t)とおくと

\begin{array}{ll} m_1\dot{v}_1(t)=-k(x_1(t)-x_2(t))+f_1(t)+w_1(t) \\ m_2\dot{v}_2(t)=-k(x_2(t)-x_1(t))+w_2(t)\\ \end{array}

【3】 状態方程式??を考えると,平衡状態では
つぎが成り立つ。

0=-\frac{K_TK_E}{JR}\omega^*+\frac{K_T}{JR}v^*-\frac{1}{J}\tau_L^*

式??からこれを辺々引いて,つぎの状態方程式を得る。

\frac{d}{dt}(\omega(t)-\omega^*) = -\frac{K_TK_E}{JR}(\omega(t)-\omega^*) +\frac{K_T}{JR}(v(t)-v^*)-\frac{1}{J}(\tau_L(t)-\tau_L^*)

【4】 シミュレーション結果はつぎのとおり。


8. 伝達関数から状態方程式を導く

線形系を知るということは,インパルス応答を知るということにほかならない。というのは,どのような入力が与えられても,インパルス応答とのたたみこみ積分を行って出力(零状態応答)が得られるからである。インパルス応答をラプラス変換したものを「伝達関数」と呼ぶ。ステップ応答の微分がインパルス応答であり,周波数応答の逆フーリエ変換がやはりインパルス応答であるから,実際には,「システム同定」実験を通してステップ応答や周波数応答を得て(伝達関数を陰に陽に求めて),さまざまな制御方式を設計することがよく行われている。

一方,これまで私たちは,物理法則に基づいて状態方程式と出力方程式からなる状態空間表現を得ていたが,もし同定実験を通して得られた伝達関数から状態空間表現を求めることができれば,これまで学んださまざまな特性解析・コントローラ設計の方法を適用できるであろう。インパルス応答や伝達関数は線形系の入出力関係を表すので「外部記述」と呼ばれ,これに対応する状態空間表現は「内部記述」と呼ばれる。そして,外部記述から内部記述を求めることを「実現」という。

本章では,状態空間表現の構造のより深い理解のためもあって,この実現問題を扱う。実は,与えられた外部記述に対応する状態空間表現は無数に存在する。これは状態空間の次数と座標軸の取り方に自由度が存在するためである。特に,最小次数の状態空間表現を求めることを「最小実現」という。興味深いのは,最小次数が可制御かつ可観測の性質により特徴付けられることである。また,座標軸の適切な取り方の一つに,特に信号処理分野で有用性が認められている「平衡実現」と呼ばれるものがある。

7. 非線形の運動方程式から始める

これまで扱ってきた制御対象は,その運動方程式が線形の微分方程式で表される「線形系」としていた。ところが実際には,例えば,ロボットアームの運動方程式に見られるように,非線形の微分方程式で表される「非線形系」が多い。このとき,これまで学んできた線形系に対する制御理論(LQ制御やLQI制御)はどのように役立つのだろうか。

まず,運動方程式(高階の非線形微分方程式)から「非線形状態方程式」を1階の連立非線形微分方程式として求めておく。つぎに,物理的な釣り合いの状態,すなわち「平衡状態」を考える。その平衡状態付近での運動を,「線形化」された状態方程式で表し,これに対して例えばLQI制御を適用する。

このアプローチは1入力1出力系ばかりでなく多入力多出力系に対しても適用できて,これまで数多くの制御問題の解決が図られてきた。その実績こそが,前章までに述べてきた線形システム制御論の有用性を示唆し,これを学ぶことの強い動機付けを与えている。

そこで,本章では,運動方程式が非線形の微分方程式で表される対象の制御問題として,水タンクの水位制御問題を考え,LQI制御を適用しながら線形システム制御論に基づくアプローチを説明する。皆さんは,計算機を用いて,LQIコントローラの数値計算と閉ループ系の時間応答シミュレーションに挑戦してほしい。そして,できれば適当な制御実験を試みてほしい。自分で考えたアルゴリズムで物をうまく動かすことができたときは,きっと小躍りしたくなるであろう。

6. 定値外乱を抑制する

これまで,安定性の定義,安定化の可能性とその方策を学んできた。その結果,私たちは,閉ループ系の零入力応答を適切に安定化する,すなわち,ある初期状態から出発して速やかに零状態(平衡状態)にもっていくことはできる。この零入力応答は,適当な強さのインパルス入力に対する零状態応答に等しいことを思い出してほしい。

例えば,一定速度で回転しているモータを考え,これに瞬間的な負荷がかかって,その速度が乱されたとする。このときは,状態フィードバックによってもとの速度に復帰させることができる。しかし,モータで何か対象物を回転させるときは,モータ軸への持続的な負荷が生じ,極端な場合はモータは回転することさえできないであろう。同様に,飛行機,船舶,車などのビークルに強い向い風が吹けば定速を維持することは困難となるであろう。

また,モータやビークルでは,速度の設定値を変えることが制御の大切な役割なのに,そのための手段はまったく学んでいない。

したがって,これまでの知識だけでは多くの現実問題を解くことはできない。実際の問題では,必ずといっていいほど,さまざまな「定値外乱」の影響を受ける状況のもとで,「設定値」の保持や,ときには変更を要求されるからである。このような問題を解決するためには,これまでの制御方式に加えて,まず「積分動作」の導入による構造的な改良が必要となる。本章では,定値外乱の影響の除去や,設定値の変更に対応できるように,積分動作を入れて安定化を行う方策について学ぶ。その安定化にLQ制御を適用するとき,その制御方式を簡潔に表すために「LQI制御」と呼ぶ。

6.1 1次系における定値外乱の影響を抑制する
6.1.1 定値外乱の影響

1章で述べた直流モータの状態方程式(1.11)すなわち

\displaystyle{(1)\quad \dot{\omega}(t)=-\frac{K_TK_E}{JR}\omega(t)+\frac{K_T}{JR}v(t)-\frac{1}{J}\tau_L(t) }

を思い出そう。これは,二つの入力変数v(t),\,\tau_L(t)をもち,印加電圧v(t)は操作できるが,負荷トルク\tau_L(t)は操作できない外乱入力であった。ここで,負荷トルクは時間によらずつねに一定とし,その値は測定できないとする。

以下では,この状況を一般化した次式で表される1次系を考える。

\displaystyle{(2)\quad \dot{x}(t)=ax(t)+bu(t)+w }

ここで,u(t)操作変数w定値外乱(constant disturbance)と呼ばれる。そのブロック線図を図6.1に示す。


図6.1 外乱入力をもつ1次系

いま,1次系(2)に対して,状態フィードバック

\displaystyle{(3)\quad u(t)=-fx(t) }

を行うと,閉ループ系

\displaystyle{(4)\quad \dot{x}(t)=(a-bf)x(t)+w }

を得る。ここで,a-bf<0とする。この時間応答は

\displaystyle{(5)\quad x(t)=e^{(a-bf)t}x(0)+\int_0^t e^{(a-bf)(t-\tau)}w\,d\tau }

となる。第2項は

\displaystyle{(6)\quad w\left[\frac{1}{-(a-bf)}e^{(a-bf)(t-\tau)}\right]_0^t =\frac{w}{-(a-bf)}(1-e^{(a-bf)t}) }

と計算できるので,a-bf<0を考慮すると次式を得る。

\displaystyle{(7)\quad x(t)\ \rightarrow\ \frac{w}{-(a-bf)} \quad (t\rightarrow\infty) }

すなわち,定値外乱が加わると

\displaystyle{(8)\quad x(t)\ \rightarrow\ 0 \quad (t\rightarrow\infty) }

とはならない。このことをシミュレーション例で確認しよう。

シミュレーション 6.1
1次系\dot{x}(t)=x(t)+u(t)+wの状態フィードバックu(t)=-2x(t)による閉ループ系において,x(0)=1.5w=1が加わるときの時間応答を,定値外乱がない場合(破線)とともに図6.2に示す。


図6.2 定値外乱の影響

このような定値外乱の影響を取り除くために,式(3)の代わりに

\displaystyle{(9)\quad u(t)=-fx(t)+v(t) }

のような方策を用いる。このときの閉ループ系

\displaystyle{(10)\quad \dot{x}(t)=(a-bf)x(t)+bv(t)+w }

の時間応答は

\displaystyle{(11)\quad x(t)=e^{(a-bf)t}x(0)+\int_0^t e^{(a-bf)(t-\tau)} \underbrace{(bv(t)+w)}_{to be 0} \,d\tau }

となる。これより,外乱を打ち消すように

\displaystyle{(12)\quad v(t)=-\frac{1}{b}w }

と選べば,式(8)を達成できる。しかしながら,wの値は未知としているので,式(12)は実際には使用できない。そこで,どのようにして式(12)の右辺を近似的に作り出すかが問題となる。

図6.2の下に次の問を設定していましたが,紙幅の関係で省略しました。

図6.2で,x(\infty)=-\displaystyle{\frac{w}{a-bf}}であることを確かめなさい。

6.1.2 積分動作を加える

この方法の基礎となるアイデアは

\displaystyle{(13)\quad x_I(t)=\int_0^tx(\tau)\,d\tau \rightarrow constant \quad (t\rightarrow\infty) }

を達成できればx(t)\rightarrow 0 \ (t\rightarrow\infty)が成り立つ,すなわち,積分値が一定値に収束すれば,その被積分項は零に収束するというものである。このように何らかの変数を積分する仕組みのことを,積分動作(integral action)と呼んでいる。

そこで,式(13)を満足するx_I(t)をフィードバックして,定値外乱を打ち消すことが考えられる。すなわち,式(12)の代わりに,v(t)=-f_Ix_I(t)を用いて(f_Iはスケーリングファクタ),結果的に

\displaystyle{(14)\quad v(t) \rightarrow -\frac{1}{b}w \quad (t\rightarrow\infty) \quad \Leftrightarrow \quad x_I(t) \rightarrow \frac{w}{f_Ib} \quad (t\rightarrow\infty) }

を達成することができないか検討しよう。

まず,式(2)と,式(13)のx_I(t)の定義から得られる積分器(integrator)

\displaystyle{(15)\quad \dot{x}_I(t)=x(t) \quad (x_I(0)=0) }

を合わせて

\displaystyle{(16)\quad \underbrace{ \left[\begin{array}{c} \dot{x}(t) \\ \dot{x}_I(t) \end{array}\right] }_{\dot{x}_E(t)} = %\underbrace{ \left[\begin{array}{cc} a & 0 \\ 1 & 0 \end{array}\right] %}_{A_E} \underbrace{ \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] }_{x_E(t)} + %\underbrace{ \left[\begin{array}{c} b \\ 0 \end{array}\right] %}_{B_E} u(t) + \underbrace{ \left[\begin{array}{c} w \\ 0 \end{array}\right] }_{w_E} }

を構成する。式(16)に

\displaystyle{(17)\quad u(t)=-fx(t) %\underbrace{ -f_Ix_I(t) %}_{+v(t)} =- \underbrace{ \left[\begin{array}{cc} f & f_I \end{array}\right] }_{F_E} \underbrace{ \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] }_{x_E(t)} }

を代入すると,閉ループ系は

\displaystyle{(18)\quad \underbrace{ \left[\begin{array}{c} \dot{x}(t) \\ \dot{x}_I(t) \end{array}\right] }_{\dot{x}_E(t)} = \underbrace{ \left[\begin{array}{cc} a-bf & -bf_I \\ 1 & 0 \end{array}\right] }_{A_{EF}} \underbrace{ \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] }_{x_E(t)} + \underbrace{ \left[\begin{array}{c} w \\ 0 \end{array}\right] }_{w_E} }

となる。そのブロック線図を図6.3に示す。


図6.3 積分動作を入れる

さて,A_{EF}は安定行列であるとする。式(18)の解は

\displaystyle{(19)\quad x_E(t)=e^{A_{EF}t}{x_E}(0) +\int_0^t e^{A_{EF}(t-\tau)}{w_E}\,d\tau }

となる。ここで,A_{EF}が正則であることに注意して,第2項は

\displaystyle{(20)\quad \left[-e^{A_{EF}(t-\tau)}A_{EF}^{-1}\right]_0^t{w_E} =-(I-e^{A_{EF}t})A_{EF}^{-1}{w_E} }

と計算できるから

\displaystyle{(21)\quad x_E(t)=e^{A_{EF}t}{x_E}(0) -(I-e^{A_{EF}t})A_{EF}^{-1}{w_E} }

と表される。したがって

\displaystyle{(22)\quad x_E(t)\ \rightarrow\ -A_{EF}^{-1}{w_E}\quad (t\rightarrow\infty) }

を得る。ここで A_{EF}^{-1}= \displaystyle{\frac{1}{bf_I}} \left[\begin{array}{cc} 0 & bf_I \\ -1 & a-bf \end{array}\right]だから

\displaystyle{(23)\quad \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] \ \rightarrow\ \left[\begin{array}{c} 0 \\ \displaystyle{\frac{w}{bf_I}} \end{array}\right] \quad (t\rightarrow\infty) }

となり,外乱の影響が抑制されていることと,式(14)が成り立つことが確認できる。また,式(17)からつぎが成り立つ。

\displaystyle{(24)\quad u(t) \ \rightarrow\ -\frac{w}{b} \quad (t\rightarrow\infty) }

シミュレーション6.2
1次系\dot{x}(t)=x(t)+u(t)+wの状態フィードバックu(t)=-3x(t)-x_I(t)による閉ループ系において,x(0)=1w=1が加わるときの時間応答を図6.4に示す(破線はx_I(t))。


図6.4 定値外乱の抑制

図6.4の下に次の問を設定していましたが,紙幅の関係で省略しました。

図6.4で,x_I(\infty)=\displaystyle{\frac{w}{bf_I}}u(\infty)=-\displaystyle{\frac{w}{b}}であることを確かめなさい。

また,(19),(20),(21)においe^{A_{EF}t}=\exp(A_{EF}t)であることに注意してください。

6.1.3 設定値を変更する

前節では,定値外乱の加わる1次系

\displaystyle{(25)\quad \dot{x}(t)=ax(t)+bu(t)+w\quad (b\ne0) }

に対して,制御目的

\displaystyle{(26)\quad x(t)\ \rightarrow\ 0 \quad (t\rightarrow\infty) }

を達成する方法を学んだ。ここでは,任意の定数rに対して,制御目的

\displaystyle{(27)\quad x(t)\ \rightarrow\ r \quad (t\rightarrow\infty) }

を達成する問題を考える。このr設定値(set value)と呼ぶ(参照(reference)値,目標(target)値,指令(command)値と呼ぶこともある)。


図6.5 外乱を抑制したうえで設定値を変更する

この問題の解決のためには,外乱を抑制するための図6.4を,図6.5に示すように変更すればよい。実際,式(27)をx(t)-r \ \rightarrow\ 0 \quad (t\rightarrow\infty)と書き直してみれば

\displaystyle{(28)\quad \dot{x}_I(t)=x(t)-r \qquad (x_I(0)=0) }

のような積分動作を導入すればよいことに気づくであろう。

式(25)と式(26)を合わせて

\displaystyle{(29)\quad %\underbrace{ \left[\begin{array}{c} \dot{x}(t) \\ \dot{x}_I(t) \end{array}\right] %}_{\dot{x}_E(t)} = %\underbrace{ \left[\begin{array}{cc} a & 0 \\ 1 & 0 \end{array}\right] %}_{A_E} %\underbrace{ \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] %}_{x_E(t)} + %\underbrace{ \left[\begin{array}{c} b \\ 0 \end{array}\right] %}_{B_E} u(t) + %\underbrace{ \left[\begin{array}{c} w \\ -r \end{array}\right] %}_{w_E} }

を構成する。いま,式(29)に

\displaystyle{(30)\quad u(t)=-fx(t)-f_Ix_I(t) =- %\underbrace{ \left[\begin{array}{cc} f & f_I \end{array}\right] %}_{F_E} %\underbrace{ \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] %}_{x_E(t)} }

を代入すると

\displaystyle{(31)\quad %\underbrace{ \left[\begin{array}{c} \dot{x}(t) \\ \dot{x}_I(t) \end{array}\right] %}_{\dot{x}_E(t)} = \underbrace{ \left[\begin{array}{cc} a-bf & -bf_I \\ 1 & 0 \end{array}\right] }_{A_{EF}} %\underbrace{ \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] %}_{x_E(t)} + %\underbrace{ \left[\begin{array}{c} w \\ -r \end{array}\right] %}_{w_E} }

となる。A_{EF}が安定行列のとき,6.1.2項と同様にして

\displaystyle{(32)\quad \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] \ \rightarrow\ \underbrace{ \displaystyle{\frac{1}{bf_I}} \left[\begin{array}{cc} 0 & bf_I \\ -1 & a-bf \end{array}\right] }_{A_{EF}^{-1}} \left[\begin{array}{c} -w \\ r \end{array}\right] \quad (t\rightarrow\infty) }

を得る。この第1行は,x(t)\,\rightarrow\,r\ (t\rightarrow\infty)となっていて,制御目的(27)が成り立つ。また,第2行から

\displaystyle{(33)\quad -f_Ix_I(t)\ \rightarrow\ -\frac{1}{b}w -\frac{a-bf}{b}r \quad (t\rightarrow\infty) }

を得る。ここで,設定値rは既知だから,第2項はあらかじめフィードフォワードすることができる。この場合の制御方式はつぎのように表される。

\displaystyle{(34)\quad u(t)=-fx(t)+f_I\underbrace{\int_0^t(r-x(\tau))\,d\tau}_{-x_I(t)} +\underbrace{\left(-\frac{a-bf}{b}\right)}_{f_r}r }

また,つぎが成り立つ。

\displaystyle{(35)\quad u(t) \ \rightarrow\ u_\infty=-\frac{1}{b}(w+ar) \quad (t\rightarrow\infty) }

それでは,フィードフォワード項による速応性改善の効果をシミュレーションによって確かめよう。

シミュレーション6.3
1次系\dot{x}(t)=x(t)+u(t)+wの状態フィードバックu(t)=-3x(t)-x_I(t)による閉ループ系において,x(0)=0w=0のとき,設定値r=1へ追従する様子を図6.6,図6.7に示す(破線はx_I(t))。


図6.6 設定値への追従(フィードフォワードなし)


図6.7 設定値への追従(フィードフォワードあり)

問6.1
図6.6と図6.7で,x_I(\infty)=\displaystyle{\frac{w}{bf_I}+\frac{a-bf}{bf_I}r}u(\infty)=\displaystyle{-\frac{1}{b}(w+ar)}であることを確かめなさい。

6.2 積分型コントローラ
6.2.1 積分動作を加えた状態フィードバック

次式で表される可制御かつ可観測で,m入力p出力をもつn次系を考える。

\displaystyle{(36)\quad \dot{x}(t)=Ax(t)+Bu(t)+w }

\displaystyle{(37)\quad y(t)=C_Mx(t) }

ここで,状態方程式には,操作入力u(t)のほかに,定値外乱wが加わっていること,C行列をC_Mと書いたことに注意する。いま,出力変数の一部やそれらの組合せからなる新しいm個の被制御変数(controlled variables)を

\displaystyle{(38)\quad z(t)=C_Sy(t)=\underbrace{C_SC_M}_{C}x(t) }

のように選び((A,C)は可観測対とする),定値外乱に無関係に,制御目的

\displaystyle{(39)\quad z(t)\ \rightarrow\ r \quad (t\rightarrow\infty) }

を達成したい。ここで,定数ベクトルrm個の設定値からなる。もし制御目的(39)が物理的に可能とすると,ある状態x_\inftyと入力u_\inftyが確定し

\displaystyle{(40)\quad \left[\begin{array}{c} -w \\ r \end{array}\right] = \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right] \left[\begin{array}{c} x_\infty \\ u_\infty \end{array}\right] }

の関係を満足しているはずである。したがって,どのようなwrに対しても,x_\inftyu_\inftyが定まるように,被制御変数(38)を

\displaystyle{(41)\quad {\rm rank}\, \underbrace{ \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right] }_{S}=n+m }

が成り立つように選ぶものとする。


図6.8 積分動作を加えた状態フィードバックによる閉ループ系

さて,制御目的(39)を達成するために,図6.8に示すような,つぎの積分動作を加えた状態フィードバックを考える。

\displaystyle{(42)\quad u(t)=-Fx(t) -F_I\underbrace{\int_0^t(z(\tau)-r)\,d\tau}_{x_I(t)} }

ここで,第2項は積分動作を表している。このようにx_I(t)を定義すると

\displaystyle{(43)\quad \dot{x}_I(t)=z(t)-r \qquad (x_I(0)=0) }

を得る。式(36)と式(43)を合わせて

\displaystyle{(44)\quad %\underbrace{ \left[\begin{array}{c} \dot{x}(t) \\ \dot{x}_I(t) \end{array}\right] %}_{\dot{x}_E(t)} = %\underbrace{ \left[\begin{array}{cc} A & 0 \\ C & 0 \end{array}\right] %}_{A_E} %\underbrace{ \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] %}_{x_E(t)} + %\underbrace{ \left[\begin{array}{c} B \\ 0 \end{array}\right] %}_{B_E} u(t) + %\underbrace{ \left[\begin{array}{c} w \\ -r \end{array}\right] %}_{w_E} }

を得る。式(44)に,式(42)すなわち

\displaystyle{(45)\quad u(t)=- %\underbrace{ \left[\begin{array}{cc} F & F_I \end{array}\right] %}_{F_E} %\underbrace{ \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] %}_{x_E(t)} }

を代入すると,閉ループ系は,つぎのように表される。

\displaystyle{(46)\quad %\underbrace{ \left[\begin{array}{c} \dot{x}(t) \\ \dot{x}_I(t) \end{array}\right] %}_{\dot{x}_E(t)} = \underbrace{ \left[\begin{array}{cc} A-BF & -BF_I \\ C & 0 \end{array}\right] }_{A_{EF}} %\underbrace{ \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] %}_{x_E(t)} + %\underbrace{ \left[\begin{array}{c} w \\ -r \end{array}\right] %}_{w_E} }

いま,A_{EF}は安定行列であるとする。このとき,A_{EF}は正則であり,つぎのように書ける。

\displaystyle{(47)\quad A_{EF} = \underbrace{ \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right] }_{S} - \left[\begin{array}{cc} B \\ 0 \end{array}\right] \left[\begin{array}{cc} F & F_I+I_m \end{array}\right] }

よって,A_{EF}の逆行列は,公式

(P-XQ^{-1}Y)^{-1}=P^{-1}+P^{-1}X(Q-YP^{-1}X)^{-1}YP^{-1}

を用いて

\displaystyle{(48)\quad \begin{array}{ll} A_{EF}^{-1} =& S^{-1}+ S^{-1} \left[\begin{array}{cc} B \\ 0 \end{array}\right] \left(I_m- \left[\begin{array}{cc} F & F_I+I_m \end{array}\right] S^{-1} \left[\begin{array}{cc} B \\ 0 \end{array}\right] \right)^{-1} \\[5mm] &\times \left[\begin{array}{cc} F & F_I+I_m \end{array}\right] S^{-1} \end{array} }

ここで,\displaystyle{S^{-1} \left[\begin{array}{cc} B \\ 0 \end{array}\right] = \left[\begin{array}{cc} 0 \\ I_m \end{array}\right] } に注意して,整理すると

\displaystyle{(49)\quad A_{EF}^{-1} = \left[\begin{array}{cc} I_n & 0 \\ -F_I^{-1}F & -F_I^{-1} \end{array}\right] S^{-1} }

のように計算される。ところで,式(46)から

\displaystyle{(50)\quad \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] \ \rightarrow\ A_{EF}^{-1} \left[\begin{array}{c} -w \\ r \end{array}\right] \quad (t\rightarrow\infty) }

を得る。この第1ブロック行は,式(40)より

\displaystyle{(51)\quad x(t)\ \rightarrow\ \left[\begin{array}{cc} I_n & 0 \end{array}\right] S^{-1} \left[\begin{array}{c} -w \\ r \end{array}\right] =x_\infty \quad (t\rightarrow\infty) }

となって,式(40)の第2ブロック行から制御目的(39)が成り立つ。また,式(50)の第2ブロック行から

\displaystyle{(52)\quad -F_Ix_{I}(t)\ \rightarrow\ \left[\begin{array}{cc} F & I_m \end{array}\right] S^{-1} \left[\begin{array}{c} -w \\ r \end{array}\right] = Fx_\infty+u_\infty \quad (t\rightarrow\infty) }

を得る。ここで,設定値rは既知だから,rに関係した項\left[\begin{array}{cc} F & I_m \end{array}\right] S^{-1} \left[\begin{array}{c} 0 \\ r \end{array}\right]をフィードフォワードして,速応性を改善できる。すなわち,制御目的(39)を達成する制御方式は

\displaystyle{(53)\quad u(t)=-Fx(t) +F_I\underbrace{\int_0^t(r-z(\tau))\,d\tau}_{-x_I(t)} +F_rr }

のように表され,ここで,F_rはつぎのように決定できる。

\displaystyle{(54)\quad F_r= \left[\begin{array}{cc} F & I_m \end{array}\right] S^{-1} \left[\begin{array}{c} 0 \\ I_m \end{array}\right] }

6.2.2 積分動作を加えたオブザーバベース\,コントローラ

式(42)の代わりに,実際には,状態オブザーバ

\displaystyle{(55)\quad \dot{\hat{x}}(t)=(A-HC_M)\hat{x}(t)+Hy(t)+Bu(t) }

の出力を用いて

\displaystyle{(56)\quad u(t)=-F\hat{x}(t)-F_Ix_I(t) }

を実施することになる。ここでの積分動作を加えたオブザーバベース\,コントローラは,式(56)を式(55)に代入した

\displaystyle{(57)\quad \dot{\hat{x}}(t)=(A-HC_M-BF)\hat{x}(t)-BF_Ix_I(t)+Hy(t) }

と,式(43)を合わせて,つぎのように表される。

\displaystyle{(58)\quad \begin{array}{rl} \left[\begin{array}{c} \dot{\hat{x}}(t) \\ \dot{x}_I(t) \end{array}\right] =& \underbrace{ \left[\begin{array}{cc} A-HC_M-BF & -BF_I \\ 0(m\times n) & 0(m\times m) \end{array}\right] }_{A_K} \left[\begin{array}{c} \hat{x}(t) \\ x_I(t) \end{array}\right]\\[10mm] &+ \underbrace{ \left[\begin{array}{cc} H & 0(n\times m)\\ C_S & -I_m \end{array}\right] }_{B_K} \left[\begin{array}{c} y(t) \\ r \end{array}\right] \end{array} }

\displaystyle{(59)\quad u(t)= \underbrace{- \left[\begin{array}{cc} F & F_I \end{array}\right] }_{C_K} \left[\begin{array}{c} \hat{x}(t) \\ x_I(t) \end{array}\right] }

これによる閉ループ系は,式(56)を式(36)に代入した

\displaystyle{(60)\quad \dot{x}(t)=Ax(t)-BF\hat{x}(t)-BF_Ix_I(t)+w(t) }

と,式(43),(57)を合わせて

\displaystyle{(61)\quad \left[\begin{array}{c} \dot{x}(t) \\ \dot{x}_I(t) \\ \dot{\hat{x}}(t) \end{array}\right] = \left[\begin{array}{ccc} A & -BF_I & -BF \\ C & 0 & 0 \\ HC_M & -BF_I & A-HC_M-BF \end{array}\right] \left[\begin{array}{c} x(t) \\ x_I(t) \\ \hat{x}(t) \end{array}\right] + \left[\begin{array}{c} w \\ -r \\ 0 \end{array}\right] }

のように表される。そのブロック線図を図6.9に示す。


図6.9積分動作を加えたオブザーバベース\,コントローラによる閉ループ系

いま,閉ループ系(61)に,座標変換

\displaystyle{(62)\quad \left[\begin{array}{c} x(t) \\ x_I(t) \\ e(t) \end{array}\right] = \left[\begin{array}{ccc} I_n & 0 & 0 \\ 0 & I_m & 0 \\ -I_n & 0 & I_n \end{array}\right] \left[\begin{array}{c} x(t) \\ x_I(t) \\ \hat{x}(t) \end{array}\right] }

を行えば

\displaystyle{(63)\quad \left[\begin{array}{c} \dot{x}(t) \\ \dot{x}_I(t) \\\hline \dot{e}(t) \end{array}\right] = \underbrace{ \left[\begin{array}{cc|c} A-BF & -BF_I & -BF \\ C & 0 & 0 \\\hline 0 & 0 & A-HC_M \end{array}\right] }_{ A_{EF}'= \left[\begin{array}{c|c} A_{EF} & - \left[\begin{array}{cc} B \\ 0 \end{array}\right] F \\[5mm] \hline 0 & \widehat{A} \end{array}\right] } \left[\begin{array}{c} x(t) \\ x_I(t) \\\hline e(t) \end{array}\right] + \left[\begin{array}{c} w \\ -r \\\hline -w \end{array}\right] }

となり,閉ループ系の固有値は,積分動作を加えた状態フィードバックだけの閉ループ系の固有値と状態オブザーバの固有値からなる

ここで,A_{EF}は安定行列であるとする。このとき,式(63)より

\displaystyle{(64)\quad \left[\begin{array}{c} x(t) \\ x_I(t) \\\hline e(t) \end{array}\right] \ \rightarrow\ \underbrace{ \left[\begin{array}{c|c} A_{EF}^{-1} & -A_{EF}^{-1} \left[\begin{array}{cc} B \\ 0 \end{array}\right] F\widehat{A}^{-1} \\[5mm]\hline 0 & \widehat{A}^{-1} \end{array}\right] }_{A_{EF}'\,^{-1}} \left[\begin{array}{c} -w \\ r \\\hline w \end{array}\right] \quad (t\rightarrow\infty) }

すなわち

\displaystyle{(65)\quad x(t)\ \rightarrow\ x_\infty \quad (t\rightarrow\infty) }

\displaystyle{(66)\quad -F_Ix_{I}(t)\ \rightarrow\ F(x_\infty+e_\infty)+u_\infty \quad (t\rightarrow\infty) }

\displaystyle{(67)\quad e(t)\ \rightarrow\ e_\infty=\widehat{A}^{-1}w \quad (t\rightarrow\infty) }

を得る。したがって,定値外乱が存在するときは状態オブザーバに関して定常偏差が残るにもかかわらず,制御目的(39)が成り立つことがわかる。また,フィードフォーワードゲインF_rは,式(66)の導出過程から,式(54)と同様に決めてよい。


6.3 LQI設計

これまで,閉ループ系(18),(31),(46),(63)において,A_{EF}は安定行列であるとしていた。ここでは,これを満足させるための具体的手段として,前章で学んだLQ制御を使うことを考えたい。LQ制御の議論における閉ループ系は自励系(入力をもたない系)を前提にしていたが,本章における閉ループ系は入力をもつことに注意が必要である。この前提を満足させるために,定常状態との差をとって得られる偏差系(error system)が用いられる。

制御目的(39)が達成されたとき成り立つ式(40)より

\displaystyle{(68)\quad \left[\begin{array}{c} 0 \\ 0 \end{array}\right] = \left[\begin{array}{cc} A & 0 \\ C & 0 \end{array}\right] \left[\begin{array}{c} x_\infty \\ x_{I\infty} \end{array}\right] + \left[\begin{array}{c} B \\ 0 \end{array}\right] u_\infty + \left[\begin{array}{c} w \\ -r \end{array}\right] }

を得る(x_{I\infty}は定数ベクトル)。まず,(44)から式(68)を引いて,つぎの偏差系を得る。

偏差系E1:
\displaystyle{(69)\quad \frac{d}{dt} %\underbrace{ \left[\begin{array}{c} x(t)-x_\infty \\ x_I(t)-x_{I\infty} \end{array}\right] %}_{x_E(t)-x_{E\infty}} = \underbrace{ \left[\begin{array}{cc} A & 0 \\ C & 0 \end{array}\right] }_{A_{E1}} %\underbrace{ \left[\begin{array}{c} x(t)-x_\infty \\ x_I(t)-x_{I\infty} \end{array}\right] %}_{x_E(t)-x_{E\infty}} + \underbrace{ \left[\begin{array}{c} B \\ 0 \end{array}\right] }_{B_{E1}} (u(t)-u_\infty) }

この両辺を微分すれば,状態変数の中の定数ベクトルを除くことができて

偏差系E2:
\displaystyle{(70)\quad \frac{d}{dt} %\underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ z(t)-r \end{array}\right] %}_{\dot{x}_E(t)} = \underbrace{ \left[\begin{array}{cc} A & 0 \\ C & 0 \end{array}\right] }_{A_{E2}} %\underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ z(t)-r \end{array}\right] %}_{\dot{x}_E(t)} + \underbrace{ \left[\begin{array}{c} B \\ 0 \end{array}\right] }_{B_{E2}} {\dot u}(t) }

を得る。さらに,式(36)と(38)をまとめた

\displaystyle{(71)\quad \left[\begin{array}{c} {\dot x}(t)-w \\ z(t) \end{array}\right] = \underbrace{ \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right] }_{S} \left[\begin{array}{c} x(t) \\ u(t) \end{array}\right] }

から式(40)を引いて,つぎの関係式が成り立つ。

\displaystyle{(72)\quad \left[\begin{array}{c} {\dot x}(t) \\ z(t)-r \end{array}\right] = \underbrace{ \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right] }_{S} \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] }

これに基づいて,偏差系E2に座標変換を行えば

偏差系E3:
\displaystyle{(73)\quad \frac{d}{dt} %\underbrace{ \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] %}_{{\tilde x}_E(t)-{\tilde x}_{E\infty}} = \underbrace{ \left[\begin{array}{cc} A & B \\ 0 & 0 \end{array}\right] }_{A_{E3}} %\underbrace{ \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] %}_{{\tilde x}_E(t)-{\tilde x}_{E\infty}} + \underbrace{ \left[\begin{array}{c} 0 \\ I_m \end{array}\right] }_{B_{E3}} {\dot u}(t) }

を得る。ここで,つぎの関係式を用いた。

\displaystyle{(74)\quad \underbrace{ \left[\begin{array}{cc} A & 0 \\ C & 0 \end{array}\right] }_{A_{E2}} \underbrace{ \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right] }_{S} = \underbrace{ \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{cc} A & B \\ 0 & 0 \end{array}\right] }_{A_{E3}} }

\displaystyle{(75)\quad \underbrace{ \left[\begin{array}{c} B \\ 0 \end{array}\right] }_{B_{E2}} = \underbrace{ \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} 0 \\ I_m \end{array}\right] }_{B_{E3}} }

\displaystyle{(76)\quad \underbrace{ \left[\begin{array}{cc} 0 & I_m \end{array}\right] }_{C_{E2}} \underbrace{ \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right] }_{S} = \underbrace{ \left[\begin{array}{cc} C & 0 \end{array}\right] }_{C_{E3}} }

問6.2
(A_{E3},B_{E3})は可制御対であることを示しなさい。

問6.3
(A_{E3},C_{E3})は可観測対であることを示しなさい。

例えば,偏差系E1に対して,状態フィードバック

\displaystyle{(77)\quad u(t)-u_\infty=- \left[\begin{array}{cc} F & F_I \end{array}\right] \left[\begin{array}{c} x(t)-x_\infty \\ x_I(t)-x_{I\infty} \end{array}\right] }

を用いて,閉ループ系を安定化すれば,A_{EF}= \left[\begin{array}{cc} A-BF & -BF_I \\ C & 0 \end{array}\right]を安定行列とすることができる。これは,偏差系E2に対しても同様である。

以下に,偏差系E3をLQ制御により安定化して,積分動作を加えたオブザーバベース\,コントローラを構成する手順を示す。

アルゴリズム6.1 <LQI設計>

入力パラメータ1
A(n\times n),\,B(n\times m),\,C_M(p\times n)
ただし,m\le p(A,B)は可制御対,(A,C_M)は可観測対

出力パラメータ
A_K(n+m\times n+m),\,B_K(n+m\times p+m),\,C_K(m\times n+m)

ステップ1 被制御変数の決定

\left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right]が正則となるように(C=C_SC_M),セレクタ行列C_S(m\times p)を決める(一般に,多入力多出力系の場合,どの操作変数でどの被制御変数を制御するのかについて,物理的に実現可能な1対1対応を考えることが重要である。その際,被制御変数はフィードバックされるので観測量の中から選ばれなけばならない)。

ステップ2 偏差系の安定化

偏差系

\displaystyle{(78)\quad \frac{d}{dt} %\underbrace{ \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] %}_{{\tilde x}_E(t)-{\tilde x}_{E\infty}} = \underbrace{ \left[\begin{array}{cc} A & B \\ 0 & 0 \end{array}\right] }_{A_{E}} %\underbrace{ \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] %}_{{\tilde x}_E(t)-{\tilde x}_{E\infty}} + \underbrace{ \left[\begin{array}{c} 0 \\ I_m \end{array}\right] }_{B_{E}} {\dot u}(t) }

を,状態フィードバック

\displaystyle{(79)\quad {\dot u}(t)=- \left[\begin{array}{cc} K & K_I \end{array}\right] \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] }

によるLQ制御で安定化する。その際,評価関数としては

\displaystyle{(80)\quad \int_0^\infty \left( \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right]^T Q_E \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] +{\dot u}^T(t)R_E{\dot u}(t)\right)\,dt }

を用いる。ただし,(A_E,Q_E^{\frac{1}{2}})は可観測対とする。

さらに,FF_Iを,次式から計算する。

\displaystyle{(81)\quad \left[\begin{array}{cc} F & F_I \end{array}\right] = \left[\begin{array}{cc} K & K_I \end{array}\right] \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right]^{-1} }

ステップ3 オブザーバゲインHの決定

例えば,アルゴリズム5.2のステップ2を用いる。
行列B'を,(A,B')が可制御対となるように選び,重み行列W>0V>0を指定し

\displaystyle{(82)\quad \Gamma A^T+A\Gamma-\Gamma C_M^TV^{-1}C_M\Gamma+B'WB'^T=0 }

を解いて,\Gamma>0を求め,オブザーバゲインHをつぎのように定める。

\displaystyle{(83)\quad H=V^{-1}C_M\Gamma }

ステップ4 LQIコントローラの構成

C_SFF_IHから,つぎを構成する。

\displaystyle{(84)\quad \dot{x}_K(t)=A_Kx_K(t)+B_K \left[\begin{array}{c} y(t) \\ r \end{array}\right] }

\displaystyle{(85)\quad u(t)=C_Kx_K(t) }

ただし

\displaystyle{(86)\quad A_K= \left[\begin{array}{cc} A-HC_M-BF & -BF_I \\ 0\,(m\times n) & 0\,(m\times m) \end{array}\right] }

\displaystyle{(87)\quad B_K= \left[\begin{array}{cc} H & 0\,(n\times m)\\ C_S & -I_m \end{array}\right] }

\displaystyle{(88)\quad C_K=- \left[\begin{array}{cc} F & F_I \end{array}\right] }

この積分動作を加えたオブザーバベース\,コントローラを簡潔にLQIコントローラ}(LQI controller),これによる制御方式をLQI制御(LQ control with integral action)と呼ぶ。

問6.4
式(81)の妥当性について説明しなさい。

前章までに,状態フィードバックと状態オブザーバを組み合わせた安定化コントローラ(オブザーバベース\,コントローラ)を,どのような条件のもとで構成できるかについて検討し,それが可制御性と可観測性(可安定性と可検出性)であることを述べた。そこでつぎに,状態フィードバックと状態オブザーバの構成を合理的に行うことを考えたい。

まず,状態フィードバックについては,代表的な方法として最小2乗規範に基づく方法がある。これは,状態フィードバックによる閉ループ系における時間応答の2乗面積が最小となるように,状態フィードバックゲインを定めるものである。ここで,状態変数の時間的振る舞いの2乗面積ばかりでなく,どのような入力を用いて閉ループ系を構成したかを評価するために,入力変数の時間的振る舞いの2乗面積についても同時に考えることが肝要である。こうして得られる制御方式を,状態フィードバックによる「LQ制御」と呼ぶ。この状態フィードバックゲインは,「リッカチ方程式」と呼ばれる行列に関する2次方程式を解いて求めることができる。

つぎに,状態オブザーバについてであるが,オブザーバベース\,コントローラによる閉ループ系における時間応答の2乗面積の和を最小化する立場から考える。この値は,理想的な状態フィードバックだけの場合と比べて,状態オブザーバを用いるので増加(劣化)する。そこで,この増分が小さくなるようにオブザーバゲインを定めることが行われる。

本章で述べる制御方式は「最適制御」と呼ばれるものの一つで,歴史的には「LQG理論」として確率システムの制御の立場から研究されてきた。しかし,ここでは確定的な議論にとどめて,主要な結果を説明した。

4. センサが足りないが大丈夫か

前章で,状態フィードバックによって安定化可能である条件(可安定性)や閉ループ系の固有値を任意に設定できる条件(可制御性)を学んだ。ただし,すでに指摘しておいたように,状態フィードバックという制御手段は,すべての状態変数にセンサを付けること(y=x)が前提となる。しかしながら,これはいつも満足されているとは限らないので,状態変数を推定することも考えておかなければならない。現代制御は,その手段として「状態オブザーバ」を提供している。これは一つの動的システムで,これを状態空間表現で表すとき,アナログ回路や計算機で実現される。その入力は制御対象の入力と出力であり,また出力は制御対象の状態の推定値(時間の経過とともに状態の真値に近ずく)である。したがって,状態オブザーバの出力に状態フィードバックのゲイン行列をかけて,状態フィードバックを近似的に実現する方法,すなわち「オブザーバベース\,コントローラ」がよく採用される。
\par
それでは,状態オブザーバは,制御対象の状態空間表現が与えられるとき,いつでも構成可能なのであろうか。答えはノーである。状態オブザーバが構成可能な条件は「可検出性」として知られている。可検出性は,前章で学んだ可安定性と密接な関係(双対関係)にある。これについて学ぶのが本章の一つの目標である。可安定性を強めた概念として可制御性があったように,可検出性を強めた概念として「可観測性」があり,これについても本章で学ぶ。最後に,2~4章の重要な結論として,制御対象の安定化可能条件が可安定性と可検出性であり,その手段がオブザーバベース\,コントローラであることを述べる。

3. いつ安定化できるのか

前章で,状態方程式から安定性を調べる方法を学んだ。すなわち,対象が安定であるとはA行列の固有値が複素左半平面内にあるときであった。それでは,仮に対象が不安定なとき,どのような制御方法で,どのような条件のもとで安定化を図ることができるのだろうか。ここで,不安定な対象とは,自然に存在する物としては考えにくいが,人工物としては,倒立振子や磁気吸引浮上,方向不安定な各種ビークル,核融合炉のプラズマ閉じ込めなどに見られる。また,制御手段としては,「状態フィードバック」というものを考える。これは,すべての状態変数の線形結合を各入力にフィードバックするもので,各状態変数の変化を即座に入力に反映させることができるという意味で,最も強力な手段である。ただし,すべての状態変数にセンサを付けること(y=x)が前提となるので,しばしば理想的な手段である。

さて,状態フィードバックによって安定化可能である条件は「可安定性」と呼ばれるが,これを直接明らかにするのは容易ではない。そこで,状態フィードバックによってA行列のすべての固有値を複素平面のどこにでも移動できる条件に等価な「可制御性」に着目する。これにより,線形性が保たれる限りにおいては,どのように速い時間応答でも実現できることになる。これに対し,可安定性は,A行列の不安定な固有値のみを複素左半平面内に移動できればよく,安定な固有値は移動できなくともよい。このような可安定性と可制御性の微妙な違いを,最終的には数学的理論展開を追って理解してほしい。

特異値分解(singular-value decomposition)

——————————————————————————————————————–
定理 Aをサイズm\times nの行列({\rm rank}\,A=k)とする。このとき,サイズm\times mの直交行列Uとサイズn\times nの直交行列Vが存在して

(1)   \begin{equation*} A=U\Sigma V^T \end{equation*}

が成り立つ。ここで,サイズm\times nの行列\Sigmaは,次を満たす(X(m\times n)は,Xがサイズn\times nの行列であることを示す(\timesは行数と列数の区切)。)。

(2)   \begin{equation*} \Sigma= \left[\begin{array}{cc} \Sigma_1(k\times k) & 0\,(k\times n-k)\\ 0\,(m-k\times k) & 0\,(m-k\times n-k) \end{array}\right] \end{equation*}

(3)   \begin{equation*} \Sigma_1={\rm diag}\{\sigma_1,\cdots,\sigma_k\} \quad (\sigma_1\ge\cdots\ge\sigma_k>0) \end{equation*}

を満足するものである。
証明 A^TA\ge 0だから,仮定{\rm rank}\,A=kより,A^TAk個の正固有値とn-k個の零固有値をもち,互いに直交する固有ベクトルをもつ。そこで,A^TAの固有値の正の平方根を,大きい順に,\sigma_1\ge\cdots\ge\sigma_k>\sigma_{k+1}=\cdots=\sigma_{n}=0のように表し,対応する固有ベクトルv_1,\cdots,v_nv_i^Tv_i=1\ (i=j),\ 0\ (i\ne j)を満足するようにとることができる。いま,\Sigma_1を式3のように,また

(4)   \begin{equation*} V= [\begin{array}{cc} V_1 & V_2 \end{array}] = [\begin{array}{ccc|ccc} v_1 & \cdots & v_k & v_{k+1}& \cdots & v_n \end{array}] \end{equation*}

とおくと,Vは直交行列となり,つぎが成り立つ。

(5)   \begin{equation*} A^TAV_1=V_1\Sigma_1^2,\quad \end{equation*}

(6)   \begin{equation*} A^TAV_2=0 \end{equation*}

6の第2式の左から,V_2^Tをかけて

(7)   \begin{equation*} V_2^TA^TAV_2=0\Rightarrow AV_2=0 \end{equation*}

また,U_1=AV_1\Sigma_1^{-1}とおくと,式6の第1式からU_1^TU_1=I_{k}を得る。そこで,U_2U= [\begin{array}{cc} U_1 & U_2 \end{array}]が直交行列となるように選ぶと

(8)   \begin{equation*} U^TAV= \left[\begin{array}{cc} U_1^TAV_1 & U_1^TAV_2\\ U_2^TAV_1 & U_2^TAV_2 \end{array}\right] = \left[\begin{array}{cc} \Sigma_1 & 0\\ U_2^TU_1\Sigma_1 & 0 \end{array}\right] =\Sigma \end{equation*}

が成り立つ。これより,式1を得る。
——————————————————————————————————————–

——————————————————————————————————————–
例題 行列A= \left[\begin{array}{ccc} 1 & 0 & 0 \\ 0 & 1 & 1 \end{array}\right]の特異値分解は

(9)   \begin{equation*} A= \underbrace{ \left[\begin{array}{cc} 0 & 1 \\ 1 & 0 \end{array}\right] }_{U} \underbrace{ \left[\begin{array}{ccc} \sqrt{2} & 0 & 0 \\ 0 & 1 & 0 \end{array}\right] }_{\Sigma} \underbrace{ \left[\begin{array}{ccc} 0 & 1 & 0 \\ \frac{1}{\sqrt{2}} & 0 & -\frac{1}{\sqrt{2}} \\ \frac{1}{\sqrt{2}} & 0 & \frac{1}{\sqrt{2}} \end{array}\right]^T }_{V^T} \end{equation*}

のように与えられることを確かめなさい。
解答 A^TAのサイズは3\times 3であるが,AA^Tのサイズは2\times 2であるので,AA^Tを計算すると\left[\begin{array}{cc} 1 & 0 \\ 0 & 2 \end{array}\right]となる(サイズm\times nの行列Aの特異値を手計算で求めるには,A^TAAA^Tのが同じ非零固有値をもつことから,サイズの小さいほうの固有値計算を行えばよい。)。これから,AA^Tの固有値は1,2で,その正の平方根1,\sqrt{2}が特異値であるが,定理??\Sigma_1の対角成分の特異値は大きい順に並べる約束であるので

(10)   \begin{equation*} \underbrace{ \left[\begin{array}{cc} 1 & 0 \\ 0 & 2 \end{array}\right] }_{AA^T} \underbrace{ \left[\begin{array}{cc} 0 & 1 \\ 1 & 0 \end{array}\right] }_{U} = \underbrace{ \left[\begin{array}{cc} 0 & 1 \\ 1 & 0 \end{array}\right] }_{U} \underbrace{ \left[\begin{array}{cc} 2 & 0 \\ 0 & 1 \end{array}\right] }_{\Sigma_1^2} \end{equation*}

のように,U\Sigma_1が決まる。つぎに,Vについては,視察によって

(11)   \begin{equation*} \underbrace{ \left[\begin{array}{ccc} 1 & 0 & 0 \\ 0 & 1 & 1 \\ 0 & 1 & 1 \end{array}\right] }_{A^TA} \underbrace{ \left[\begin{array}{ccc} 0 & 1 & 0 \\ x & 0 & -y \\ y & 0 & x \end{array}\right] }_{V} = \underbrace{ \left[\begin{array}{ccc} 0 & 1 & 0 \\ x & 0 & -y \\ y & 0 & x \end{array}\right] }_{V} \underbrace{ \left[\begin{array}{ccc} 2 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 0 \end{array}\right] }_{{\rm diag}\{\Sigma_1^2,0\}} \end{equation*}

としてよい。ここで\left[\begin{array}{ccc}x \\y \end{array}\right]\left[\begin{array}{ccc}-y \\x \end{array}\right]が直交しており,x^2+y^2=1の制約がある。式(11)からx=yが出て,x=y=\frac{1}{\sqrt{2}}と定まる。
——————————————————————————————————————–

それでは,階段化アルゴリズムについて述べる。

——————————————————————————————————————–
アルゴリズム<階段化アルゴリズム(staircase algorithm)>
入力パラメータ m,\,n,\,A(n\times n),\,B(n\times m)
出力パラメータ k,\,m_1,\cdots,m_k,\,T(n\times n)
ステップ1 初期化

j=0,\,s=0,\,m_0=m,\,T^{(1)}=I_n,\,B^{(1)}=B,\,A^{(1)}=Aとおく。

ステップ2 B^{(j)}の階数決定

jj+1に更新し,B^{(j)}(n-s\times m_{j-1})に対して,つぎの特異値分解を行い,B^{(j)}の階数m_jU^{(j)}(n-s\times n-s)を求める。

(12)   \begin{equation*} B^{(j)}=U^{(j)}\Sigma^{(j)}V^{(j)}\,^T \end{equation*}

(13)   \begin{equation*} \Sigma^{(j)}= \left[\begin{array}{cc} \Sigma_1^{(j)} & 0 \\ 0 & 0 \end{array}\right] \end{equation*}

もしm_j=n-sまたはm_j=0ならば,k=j,\,T=T^{(j)}と設定して終了する。そうでないときは,ステップ3へ行く。

ステップ3 座標変換

W^{(j)}(n\times n)を,U^{(j)}(n-s\times n-s)を用いて

(14)   \begin{equation*} W^{(j)}= \left[\begin{array}{cc} I_s & 0 \\ 0 & U^{(j)} \end{array}\right] \end{equation*}

と設定し,A^{(j+1)}(n\times n)T^{(j+1)}(n\times n)を,次式から求める。

(15)   \begin{equation*} A^{(j+1)}=W^{(j)}\,^TA^{(j)}W^{(j)},\ T^{(j+1)}=T^{(j)}W^{(j)} \end{equation*}

ステップ4 B^{(j)}の切り出し

ss+m_jに更新し,A^{(j+1)}(n\times n)に対してつぎの分割を行い,B^{(j+1)}(n-s\times m_j)を得て,ステップ2へ戻る。

(16)   \begin{equation*} A^{(j+1)}= \left[\begin{array}{ccc} X(s\times s-m_j) & X(s\times m_j) & X(s\times n-s) \\ 0 & B^{(j+1)} & X \end{array}\right] \end{equation*}

——————————————————————————————————————–

——————————————————————————————————————–
例題
階段化アルゴリズムを用いて,つぎの可安定性を調べなさい。

(17)   \begin{equation*} A= \left[\begin{array}{ccc} 0 &  0 & 0 \\ 0 & -1 & 1 \\ 0 & 0  & 0 \end{array}\right] ,\ B= \left[\begin{array}{cc} 1 & 0 \\ 0 & 1 \\ 0 & 1 \end{array}\right] \end{equation*}

解答 階段化アルゴリズムの各ステップは,つぎのように計算される。

ステップ1j=0,\,s=0,\,m_0=2,\,T^{(1)}=I_3,\,B^{(1)}=B,\,A^{(1)}=Aとおく。
ステップ2j=1とし,B^{(1)}に対する特異値分解

(18)   \begin{equation*} B^{(1)}= \underbrace{ \left[\begin{array}{ccc} 0 & 1 & 1 \\ \frac{1}{\sqrt{2}} & 0 & -\frac{1}{\sqrt{2}} \\ \frac{1}{\sqrt{2}} & 0 & \frac{1}{\sqrt{2}} \end{array}\right] }_{U^{(1)}} \underbrace{ \left[\begin{array}{cc} \sqrt{2} & 0 \\ 0 & 1 \\ 0 & 0 \end{array}\right] }_{\Sigma^{(1)}} \underbrace{ \left[\begin{array}{cc} 0 & 1 \\ 1 & 0 \end{array}\right]^T }_{V^{(1)T}} \end{equation*}

から,B^{(1)}の階数m_1=2U^{(1)}が求まる。
ステップ3m_1\ne3かつm_1\ne0なので,ステップ4へ行く。
ステップ4W^{(1)}=U^{(1)}とし,A^{(2)}=W^{(1)}\,^TA^{(1)}W^{(1)}を求めると

(19)   \begin{equation*} A^{(2)}= \left[\begin{array}{cc|c} 0 & 0 & 1  \\ 0 & 0 & 0  \\ \hline 0 & 0 & -1 \end{array}\right] \end{equation*}

となり,T^{(2)}=W^{(1)}とする。
ステップ5s=2とする。A^{(2)}に対する分割から,B^{(2)}= [\begin{array}{cc} 0 & 0 \end{array}]を得る。
ステップ2′j=2とする。B^{(2)}は零行列なので,その階数はm_2=0である。
ステップ3′m_2=0なので,k=2,\,T=T^{(2)}のように設定する。

以上から,式(24)は,つぎのように得られた。

(20)   \begin{equation*} \left[\begin{array}{cc} T^TB & T^TAT \end{array}\right]= \left[\begin{array}{ccc} B_1 & A_1 & X \\ 0   & B_2 & A_2 \end{array}\right]= \left[\begin{array}{cc|cc|c} 0 & \sqrt{2} &  0 & 0 & 1 \\ 1 & 0        &  0 & 0 & 0 \\ \hline 0 & 0        &  0 & 0 & -1 \end{array}\right] \end{equation*}

ここで,B_2は零行列で,A_2=-1<0より,条件S0が成り立つ。よって,(A,B)は可安定である。
——————————————————————————————————————–

まず,可安定性と可制御性は座標変換により不変であることを示す。
——————————————————————————————————————–
定理 (A,B)が可安定対[可制御対]ならば,(TAT^{-1},TB)も可安定対[可制御対]である。
証明 n次系\dot{x}(t)=Ax(t)+Bu(t)に対するu(t)=-Fx(t)による安定化可能性は,座標変換後の\dot{x}'(t)=TAT^{-1}x'(t)+TBu(t)に対するu(t)=-FT^{-1}x'(t)による安定化可能性と等価である。実際

(21)   \begin{equation*} {\rm det}(\lambda-TAT^{-1}+TB\cdot FT^{-1})= {\rm det}(\lambda-A+BF) \end{equation*}

が成り立ち,閉ループ系の固有値は,座標変換によって不変だからである。

可制御性については,シルベスターの公式

{\rm rank}\,X+{\rm rank}\,Y-n\le{\rm rank}\,XY\le{\rm min}\{{\rm rank}\,X,{\rm rank}\,Y\} (nXの列数)

を用いて

(22)   \begin{equation*} \begin{array}{ll} {\rm rank}\, [\begin{array}{cccc} TB & TAT^{-1}TB & \cdots & (TAT^{-1})^{n-1}TB \end{array}] \\ ={\rm rank}\, [\begin{array}{cccc} B & AB &\cdots & A^{n-1}B \end{array}] \end{array} \end{equation*}

が成り立つことから,定理の主張を得る。
——————————————————————————————————————–

この定理に基づいて,可安定性や可制御性の判別をしやすいように座標変換行列Tを選ぶことを考える。すなわち,n次系\dot{x}(t)=Ax(t)+Bu(t)に対して,つぎを満足するように座標変換できるかどうかを調べる。

(23)   \begin{equation*} \left[\begin{array}{c} \dot{x}_1'(t) \\ \dot{x}_2'(t) \end{array}\right] = \underbrace{ \left[\begin{array}{cc} A_{\#} & X   \\ 0      & A_k \end{array}\right] }_{TAT^{-1}} \left[\begin{array}{c} x_1'(t) \\ x_2'(t) \end{array}\right] + \underbrace{ \left[\begin{array}{c} B_{\#} \\ 0 \end{array}\right] }_{TB} u(t) \end{equation*}

ただし,(A_{\#},B_{\#})は可制御対とする。このブロック線図を描いてみればわかるように,状態x_2'は入力uの影響を受けることはない。したがって,可制御性は成立しないこと,そしてA_kが安定行列かどうかが可安定性と関係していることは想像がつくであろう。以下では,この判定法が導かれること示す。

「式(23)のように座標変換できれば,A_kが安定行列のとき可安定性が成り立ち,そうでないときは可安定性は成り立たない。また,式(23)のように座標変換できなければ,可制御性が成り立っている。」

つぎの結果は,階段化アルゴリズムを用いて得られる。
——————————————————————————————————————–
定理
サイズn\times nの行列Aとサイズn\times mの行列Bに対して,つぎのように変換する直交行列Tが存在する。

(24)   \begin{equation*} \left[\begin{array}{cc} T^TB & T^TAT \end{array}\right]= \left[\begin{array}{c|cccc|c} B_1 & A_1 & X      & \cdots  & X       & X \\ 0   & B_2 & A_2    & \cdots  & X       & X \\ \vdots &     & \ddots & \ddots  & \vdots  & \vdots \\ 0   & 0   & \cdots & B_{k-1} & A_{k-1} & X      \\ \hline 0   & 0   & \cdots & 0       & B_k     & A_k \end{array}\right] \end{equation*}

ここで,kは正整数,またB_jのサイズをm_j\times m_{j-1}とするとき(m_0=m,\ m_1+\cdots+m_k=n),つぎが成り立つ(B_1,\cdots,B_{k-1}は横長の形状となり,行フルランクをもつ。また,B_kは横長で,行フルランクをもつ行列または零行列のどちらかである。)。

(25)   \begin{equation*} {\rm rank}\,B_j=m_j\quad (j=1,\cdots,k-1) \end{equation*}

(26)   \begin{equation*} {\rm rank}\,B_k=m_k\ または\ 0 \end{equation*}

(27)   \begin{equation*} m\ge m_1 \ge \cdots \ge m_{k-1} \ge m_k \end{equation*}

このとき,(A,B)が可制御対であるための必要十分条件は

条件C0: 式(24)において,B_k\ne 0

が成り立つことである。また,(A,B)が可安定対であるための必要十分条件は、つぎが成り立つことである。

条件S0: 式(24)において,B_k=0のとき,A_kは安定行列
——————————————————————————————————————–

B_k=0のとき,式(23)のような形式が得られていることは,式\(ref{eq3.4.3})のTを上のT^{-1}=T^Tと読み替えれば明らかであろう。また,可安定性と可制御性と相違は,前者がB_k\ne0ばかりでなくB_k=0を許すところにあることに注意する。可制御性の条件C0と可安定性の条件S0の妥当性は,あとで示す。

それでは,まず可制御性について,新しい条件を追加した定理を述べる。
——————————————————————————————————————–
定理 可制御性が成り立つための等価な条件

定義: 任意初期状態を,任意有限時間内に,任意状態に移動可能

条件C1: \displaystyle{\int_0^t \exp(A\tau)BB^T\exp(A^T\tau)\,d\tau>0 \quad (\forall t>0)}

条件C2: {\rm rank}\, \left[\begin{array}{cccc} B & AB & \cdots & A^{n-1}B \end{array}\right]=n

条件C3: Fを選んで,A-BFの固有値を任意に設定可能

条件C4: {\rm rank}\, \left[\begin{array}{cc} B & A-\lambda I_n \end{array}\right] =n (\lambdaAのすべての固有値)

条件C5: B^Tw=0,\ A^Tw=\lambda w \Rightarrow w=0 (\lambdaAのすべての固有値)

ここで,条件C4を満足するAの固有値\lambdaは可制御固有値,満足しない固有値は不可制御固有値と呼ばれる。
——————————————————————————————————————–

前節までに,定義\Leftrightarrow条件C1\Leftrightarrow条件C2\Rightarrow条件C3を示したので,以下では,条件C3\Rightarrow条件C0\Rightarrow条件C4\Leftrightarrow条件C5\Rightarrow条件C2を示す。

——————————————————————————————————————–
証明 <条件C3\Rightarrow条件C0>
B_k=0とすると,式(24)を用いて

(28)   \begin{equation*} {\rm det}(\lambda I_n-T^TAT+T^TB\cdot FT) ={\rm det}\left(\lambda I_n- \left[\begin{array}{cc} X & X \\ 0 & A_k \end{array}\right] \right) \end{equation*}

となり,閉ループ系においてA_kの固有値はそのまま残るが,C3に矛盾である。

証明 <条件C0\Rightarrow条件C4>
\lambdaを任意の複素数とする。[\begin{array}{cc} T^TB & T^TAT-\lambda I_n \end{array}]に,式24を適用すると

(29)   \begin{equation*} \left[\begin{array}{c|cccc|c} B_1 & A_1-\lambda I_{m_1} & X      & \cdots  & X       & X \\ 0   & B_2 & A_2-\lambda I_{m_2}    & \cdots  & X       & X \\ \vdots &     & \ddots & \ddots  & \vdots  & \vdots \\ 0   & 0   & \cdots & B_{k-1} & A_{k-1}-\lambda I_{m_{k-1}} & X  \\ \hline 0   & 0   & \cdots & 0       & B_k     & A_k-\lambda I_{m_k} \end{array}\right] \end{equation*}

を得る。ここで,B_1,\cdots,B_kは横長の形状となり,行フルランク(行数に等しい階数)をもつので

(30)   \begin{equation*} \left[\begin{array}{cc} T^TB & T^TAT-\lambda I_n \end{array}\right] = T^T \left[\begin{array}{cc} B & A-\lambda I_n \end{array}\right] \left[\begin{array}{cc} I_n & 0 \\ 0 & T \end{array}\right] \end{equation*}

に,シルベスターの公式を適用して

(31)   \begin{equation*} {\rm rank}\, [\begin{array}{cc} T^TB & T^TAT-\lambda I_n \end{array}] = {\rm rank}\, [\begin{array}{cc} B & A-\lambda I_n \end{array}] =n \end{equation*}

を得る。この階数がnより小さくなるのは,式(29)から明らかなように,\lambdaA_kの固有値に等しい場合である。したがって,\lambdaAの固有値について調べれば十分である。すなわち,条件C4を得た。

証明 <条件C4\Leftrightarrow条件C5>
条件C4は,wn次元ベクトルとして

(32)   \begin{equation*} w^T\, [\begin{array}{cc} B & A-\lambda I_n \end{array}] =0 \Rightarrow w =0 \end{equation*}

すなわち

(33)   \begin{equation*} B^Tw=0,\ A^Tw=\lambda w \Rightarrow w=0 \end{equation*}

と等価である。

証明 <条件C5\Rightarrow条件C2>
条件C5より

(34)   \begin{equation*} w^TA^{i-1}B=0\quad (i=0,1,2,\cdots) \end{equation*}

を得て

(35)   \begin{equation*} w^T\, [\begin{array}{ccccc} B & AB & \cdots & A^{n-1}B \end{array}] =0 \Rightarrow w =0 \end{equation*}

が成り立つ。これは条件C2を意味する。
——————————————————————————————————————–

それでは,可安定性が成り立つための条件について考えよう。可制御であれば可安定であるので,問題は,不可制御で可安定の場合をどのように特徴付けるかである。不可制御であれば条件C0が不成立なので,B_k=0となり,A_kの固有値だけは状態フィードバックによってどうすることもできない。このとき可安定であるためには,A_kが安定行列であることを前提とするしかない。これは条件S0そのものであり,可安定性の条件としての妥当性が得られた。

このことを踏まえて,可安定性の条件はつぎのように表される。

——————————————————————————————————————–
定理 可安定性が成り立つための等価な条件

定義 状態フィードバックにより安定化可能

条件S1: {\rm rank}\, \left[\begin{array}{cc} B & A-\lambda I_n \end{array}\right] =n (\lambdaAのすべての不安定固有値)

条件S2: B^Tw=0,\ A^Tw=\lambda w \Rightarrow w=0\quad (〃)

証明 定義\Leftrightarrow条件S0だから,条件S1\Leftrightarrow条件S2との等価性を示せばよい。そのために,B_k=0のとき,式(29)の行列のランクが落ちるのは

(36)   \begin{equation*} {\rm rank}\, \left. [\begin{array}{cc} B_k & A_k-\lambda I_{m_k} \end{array}] \right|_{B_k=0} = {\rm rank}\, [\begin{array}{cc} 0 & A_k-\lambda I_{m_k} \end{array}] \le m_k \end{equation*}

(37)   \begin{equation*} {\rm rank}\, \left. \left[\begin{array}{c} B_k^T \\ A_k^T-\lambda I_{m_k} \end{array}\right] \right|_{B_k=0} = {\rm rank}\, \left[\begin{array}{c} 0 \\ A_k^T-\lambda I_{m_k} \end{array}\right] \le m_k \end{equation*}

の部分で,\lambdaA_kの固有値に一致した場合である。したがって,条件S0を保証するためには,\lambdaA_kのすべての不安定固有値を入れて,式(37)を調べ,列フルランクとなればよい。すなわち,条件C4\Leftrightarrow条件C5において,\lambdaA_kのすべての不安定固有値という条件をつければよい。
——————————————————————————————————————–

上の可制御性と可安定性のさまざまな条件のうち,理論展開では条件C5とS2が,数値計算では条件C0,C4と,条件S0,S1がよく用いられる(条件C5\Leftrightarrow条件C4と条件S2\Leftrightarrow条件S1による判定法は,PBH(Popov-Blevitch-Hautus)法と呼ばれる。)。

例えば,可制御性は状態フィードバックにより不変であることを示そう。

——————————————————————————————————————–
定理 (A,B)が可制御対ならば(A-BF,B)も可制御対である。
証明 条件C5より,任意の\lambdaに対して

(38)   \begin{equation*} B^Tw=0,\ A^Tw=\lambda w \Rightarrow\ w=0 \end{equation*}

成り立つので,これを仮定して

(39)   \begin{equation*} B^Tw=0,\ (A-BF)^Tw=\lambda w \Rightarrow\ w=0 \end{equation*}

を示せばよい。実際,B^Tw=0より

(40)   \begin{equation*} (A-BF)^Tw=A^Tw-F^TB^Tw=A^Tw=\lambda w \end{equation*}

を得るが,このようなwは,仮定(38)より,w=0でなければならない。
——————————————————————————————————————–