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

私たちが制御したい対象が物である場合,その振る舞いは物理法則により記述される。一般に,動的な振る舞いを表す物理法則は微分方程式で表される。微分方程式というと尻込みする人もいるかもしれないが,最初は,ニュートンの運動第2法則m\ddot{x}=Fに代表されるような,簡単な線形の微分方程式から始めればよい。ただし,これ1本で済むというわけにはいかない。通常は複数本の微分方程式が出てきて,やっかいなことにそれらがお互いに関連し合っている。どんな問題を解く場合も対象の性質をよく分析することが第一歩であるから,連立された微分方程式を扱える統一的表現法があると望ましい。そこで,対象を記述するいくつかの微分方程式を,1階の連立微分方程式としてまとめて,行列表現することがよく行われる。これを制御対象の「状態方程式」と呼ぶ。また,その解変数を「状態変数」という。

さて,制御を行うには対象の動的な振る舞い(状態変数の時間変化)を観測することが必要であり,そのためにセンサが用いられる。ここで,現実には,すべての状態変数にセンサを配置できるわけではなく,通常は状態変数の一部がセンサにより観測される。センサにより観測される変数と状態変数の間の関係式を,制御対象の「出力方程式」と呼ぶ。

本章では,動的な振る舞いが連立線形微分方程式で記述されるいくつかの制御対象に対して,状態方程式と出力方程式をペアにした「状態空間表現」を求める。これを図示した「ブロック線図」は,計算機シミュレーションを行う際にも有用である。

1.1 状態空間表現の導出例
1.1.1 ペースメーカ

高齢化社会の到来に伴い,より優れた福祉・医療機器の開発が工学分野の大きなテーマの一つとなっている。図1.1に示すのは,心臓のペースメーカの簡単な原理図である。これは,まず左側の閉回路でコンデンサへの充電を行い,つぎにスイッチを切り替えてできる右側の閉回路で放電を行うという動作を周期的に繰り返すことにより,心臓のペースメーカの役割を果たそうとするものである。ここでは,状態方程式を導く最初の例として,このようなRC回路における充電と放電について考える。

そのために,キルヒホッフの電圧則より,左側閉回路と右側閉回路の回路方程式を考えると,それぞれ

(1)   \begin{equation*} R_1\underbrace{C\dot{v}(t)}_{i(t)}+v(t)=u(t) \end{equation*}

(2)   \begin{equation*} R_2\underbrace{C\dot{v}(t)}_{i(t)}+v(t)=0 \end{equation*}

となる。ここで,Cはコンデンサの容量,R_1は充電回路の抵抗,R_2は心臓の等価抵抗,v(t)は時刻tにおけるコンデンサ両端の電圧,u(t)は時刻tにおけるバッテリの電圧である。


図1.1 心臓のペースメーカ

式(1)は,すでに,v(t)に関する1階の線形微分方程式であるので,両辺をR_1Cで割って,つぎの状態方程式を得る。この解変数v(t)状態変数と呼ぶ。

(3)   \begin{equation*} \dot{v}(t)=-\frac{1}{R_1C}v(t)+\frac{1}{R_1C}u(t) \end{equation*}

状態方程式(3)を図1.2のように図示し,これを状態方程式に基づくブロック線図と呼ぶ。この描き方のポイントは,式(3)の右辺を表すのに加え合わせ記号○を用いることと,また\dot vを積分してvを得て右辺と左辺を関連付けていることである。なお,加え合わせにおけるプラス符号は省略することが多い。


図1.2 ペースメーカの充電回路のブロック線図

このブロック線図から,外部より与えられる入力変数u(t)が,状態変数v(t)の微分値に影響を与え,v(t)が外部に取り出されることが見てとれる。状態変数は1個であるので,式(3)で表される動的システムを1次システム(first-order system)または1次系^*と呼ぶ。

同様に,式(2)から得られる状態方程式は

(4)   \begin{equation*} \dot{v}(t)=-\frac{1}{R_2C}v(t) \end{equation*}

であり,これによるブロック線図は図1.3のように示される。


図1.3 ペースメーカの放電回路のブロック線図

微分方程式(4)の解が

(5)   \begin{equation*} \displaystyle{v(t)=e^{-\frac{1}{R_2C}t}v(0)} \end{equation*}

と与えられることはよいであろう(式(4)に代入して確かめよ)。状態方程式(4)は入力変数をもたないが,状態変数の初期値によって,状態変数の時間的振る舞いが現れる。この意味で,1次系(4)は自励系(autonomous system) 自由系(unforced system)と呼ばれる。つぎのシミュレーション例^{**}をみてみよう。


シミュレーション1.1 式(5)で表されるコンデンサ電圧v〔{\rm V}〕の時間的振る舞いを,R_2C=1v(0)=1\,{\rm V}の場合について図1.4に示す。


図1.4 コンデンサ放電時の電圧変化

問1.1 図1.4において,時刻t=1,3,5におけるvの値を

(6)   \begin{equation*} e^{x}=1+x+\frac{1}{2}x^2+\cdots+\frac{1}{k!}x^k+\cdots \end{equation*}

によって近似計算しなさい。


*系はsystemの訳語。ここでは「××システム」を簡潔に「××系」と書く。
**本書では,時間応答のコンピュータによるシミュレーション(simulation)の欄を設けた。最終的には時間応答の数学的理解が大切であるが,まずは,なぜそのような時間的振る舞いが現れるのかを物理的イメージをもって考えながら,典型的な時間応答に親しみをもってほしい。なお,本書の数値計算については演習問題の【4】を参照のこと。

1.1.2 教室のドア

教室で物の動きを実感できるものに,図1.5に示すようなばねとダンパ^*からなる緩衝装置を付けたドアがある。これは,開いたドアをできるだけ速やかに静かに閉めるためのものである。


図1.5 緩衝装置をつけたドア

このドアの運動は回転運動であるが,話しをわかりやすくするため,図1.6に示すような等価な直線運動として調べてみよう。その出発点は,ニュートンの運動第2法則

(7)   \begin{equation*} M\ddot{r}(t)=F(t) \end{equation*}

である。ここで,Mはドアの質量,r(t)は時刻tにおけるドアの変位,F(t)は時刻tにおいてドアに働く力であり

(8)   \begin{equation*} F(t)=-D\dot{r}(t)-Kr(t)+u(t) \end{equation*}

のように表すことができる。ここで,ダンパが第1項の力を,ばねが第2項の力を与える。uは人がドアに与える力である。式(7)と式(8)より

(9)   \begin{equation*} M\ddot{r}(t)=-D\dot{r}(t)-Kr(t)+u(t) \end{equation*}

を得る。


図1.6 ドアの簡単なモデル

これは2階の線形微分方程式であるが,v(t)=\dot{r}(t)を定義すると

(10)   \begin{equation*} \dot{r}(t)=v(t) \end{equation*}

(11)   \begin{equation*} \dot{v}(t)=-\frac{D}{M}v(t)-\frac{K}{M}r(t)+\frac{1}{M}u(t) \end{equation*}

のような1階の連立線形微分方程式で表される。これらを行列表示すると

(12)   \begin{equation*} \underbrace{ \left[\begin{array}{c} \dot{r}(t) \\ \dot{v}(t) \end{array}\right] }_{\dot{x}(t)}= \underbrace{ \left[\begin{array}{cc} 0 & 1 \\ -\frac{K}{M} & -\frac{D}{M} \end{array}\right] }_A \underbrace{ \left[\begin{array}{c} r(t) \\ v(t) \end{array}\right] }_{x(t)}+ \underbrace{ \left[\begin{array}{c} 0 \\ \frac{1}{M} \end{array}\right] }_B u(t) \end{equation*}

のような状態方程式を得る^{**}。ここで,状態変数はr(t)v(t),入力変数はu(t)である。また,図1.7のようなブロック線図が得られる。


図1.7 ドアのブロック線図

さて,2個の状態変数のうち,ドアの変位rc_r倍の電圧y_r,すなわち

(13)   \begin{equation*} y_r(t)=c_r r(t) \end{equation*}

を得るセンサはあるが,ドアの速度を計測するセンサはないものとする。このとき,y_r出力変数と呼ぶ。これは,つぎの出力方程式により表される。

(14)   \begin{equation*} \underbrace{y_r(t)}_{y(t)}= \underbrace{ \left[\begin{array}{cc} c_r & 0 \end{array}\right] }_C \underbrace{ \left[\begin{array}{c} r(t) \\ v(t) \end{array}\right] }_{x(t)} \end{equation*}

以上から,ドアに対して,状態方程式(12)と出力方程式(14)からなる2次系(second-order system)としての状態空間表現を得た。

シミュレーション 式(12)において, M=1\,{\rm kg}, K=1\,{\rm N/m}r(0)=0.5\,{\rm m}v(0)=0\,{\rm m/s}u(t)=0\,{\rm N}のとき,D=5,2,0.1\,{\rm N\cdot s/m}の三つの場合について,ドア開度r〔{\rm m}〕の時間的振る舞いを図1.8に示す。


図1.8 ドア開度の時間的振る舞い

問1.2 図1.8の三つの時間応答に対応して,ドアはそれぞれどのように閉まるか説明しなさい。

*ばねとダンパの特性値を調整するためのねじを回すことにより行われる。
**本書では,\underbrace{\Delta}_{\circ}のように書いて,△を○で定義・表記する(△は○に等しいとする)。

1.1.3 直流モータ


代表的なアクチュエータとしてモータがある。例えば図1.9に示すのは,ロボットアームを駆動する直流モータである。


図1.9 直流モータ

このモデルは図1.10のように表される。


図1.10 直流モータのモデル

このとき,つぎが成り立つ。

(15)   \begin{equation*} J\dot{\omega}(t)=K_Ti(t)-\tau_L(t) \end{equation*}

(16)   \begin{equation*} L\dot{i}(t)=-Ri(t)-K_E\omega(t)+v(t) \end{equation*}

ここで,式(15)は機械系としての運動方程式であるが,電流による発生トルクの項K_T\,iを含む。K_Tはトルク定数と呼ばれる。また,式(16)は電気系としての回路方程式であるが,角速度\omega=\dot{\theta}による逆起電力の項K_E\,\omegaを含む。K_Eは逆起電力定数と呼ばれる。このように,モータは機械系と電気系の混合系という特徴をもつ。式(15)と式(16)に

(17)   \begin{equation*} \dot{\theta}(t)=\omega(t) \end{equation*}

を加えたものを行列表示すると

(18)   \begin{eqnarray*} %\begin{equation} %\renewcommand{\arraystretch}{1.3} %\begin{array}{rl} \underbrace{ \left[\begin{array}{ccc} 1 & 0 & 0 \\ 0 & J & 0 \\ 0 & 0 & L \end{array}\right] }_{E} \left[\begin{array}{c} \dot{\theta}(t) \\ \dot{\omega}(t) \\ \dot{i}(t) \end{array}\right] & = \underbrace{ \left[\begin{array}{ccc} 0 & 1 & 0 \\ 0 & 0 & K_T \\ 0 & -K_E & -R \end{array}\right] }_{A} \left[\begin{array}{c} \theta(t) \\ \omega(t) \\ i(t) \end{array}\right] \nonumber\\ & + \underbrace{ \left[\begin{array}{cc} 0 & 0 \\ 0 & -1 \\ 1 & 0 \end{array}\right] }_{B} \left[\begin{array}{cc} v(t) \\ \tau_L(t) \end{array}\right] %\end{array} %\end{equation} \end{eqnarray*}

となる^*。この左から,E^{-1}をかけて

(19)   \begin{equation*} \underbrace{ \left[\begin{array}{c} \dot{\theta}(t) \\ \dot{\omega}(t) \\ \dot{i}(t) \end{array}\right] }_{\dot{x}(t)} = \underbrace{ \left[\begin{array}{ccc} 0 & 1 & 0 \\ 0 & 0 & \frac{K_T}{J} \\ 0 & -\frac{K_E}{L} & -\frac{R}{L} \end{array}\right] }_{E^{-1}A} \underbrace{ \left[\begin{array}{c} \theta(t) \\ \omega(t) \\ i(t) \end{array}\right] }_{x(t)} + \underbrace{ \left[\begin{array}{cc} 0 & 0 \\ 0 & -\frac{1}{J} \\ \frac{1}{L} & 0 \end{array}\right] }_{E^{-1}B} \underbrace{ \left[\begin{array}{cc} v(t) \\ \tau_L(t) \end{array}\right] }_{u(t)} \end{equation*}

のような状態方程式を得る。状態方程式(19)は二つの入力変数v,\,\tau_Lをもち,vは操作できるが,\tau_Lは操作できない外乱であることに注意してほしい。

問1.3 式(19)を用いて,直流モータのブロック線図を描きなさい。

さて,この直流モータに対しては,角度\thetac_\theta倍の電圧y_\thetaと,角加速度\omegac_\omega倍の電圧y_\omegaが測れるものとすると,出力方程式は

(20)   \begin{equation*} \underbrace{ \left[\begin{array}{c} y_\theta(t) \\ y_\omega(t) \end{array}\right] }_{y(t)}= \underbrace{ \left[\begin{array}{ccc} c_\theta & 0 & 0 \\ 0 & c_\omega & 0 \end{array}\right] }_C \underbrace{ \left[\begin{array}{c} \theta(t) \\ \omega(t) \\ i(t) \end{array}\right] }_{x(t)} \end{equation*}

となる。このように,直流モータに対して,状態方程式(19)と出力方程式(20)からなる3次系(third-order system)としての状態空間表現を得た。つぎに一つのシミュレーション例を示す。

シミュレーション1.3 式(19)において,J=1.5\times10^{-5}\,{\rm kg\cdot m^2}K_T=0.03\,{\rm N\cdot m/A}L=0.001\,{\rm H}R=1\,{\rm \Omega}K_E=0.05\,{\rm V\cdot s/rad}のとき,\tau_L(t)=0\,{\rm N}v(t)=1\,{\rm V}x(0)=0の場合の\omega〔{\rm rad/s}〕i〔{\rm A}〕の時間的振る舞いを図1.11に示す


図1.11 直流モータの時間応答

ところで,私たちは物理的な感覚として,機械的な動きと電気的な動きでは速さが格段に違うことを知っている。直流モータは機械系と電気系の混合系であることを述べたが,制御目的は位置制御や速度制御のように機械系に関わるのが普通であるので,状態変数としては\theta\omegaだけでよさそうである。式(16)をみると,直流モータの電気的時定数(iの時定数^{**})は

(21)   \begin{equation*} T_e=\frac{L}{R} \end{equation*}

で与えられ,上の例ではT_e=0.001\,{\rm s}である。ところが,図1.11からわかるように,\omegaの時定数は約0.01\,{\rm s}である。したがって,電流は角速度に比べて10倍速く落ち着くので,式(16)の左辺を零とおいてみよう。すなわち

(22)   \begin{equation*} 0=-Ri(t)-K_E\omega(t)+v(t) \end{equation*}

これからiを求めて,式(15)に代入してみると

(23)   \begin{equation*} \dot{\omega}(t)=-\frac{K_TK_E}{JR}\omega(t) +\frac{K_T}{JR}v(t) -\frac{1}{J}\tau_L(t) \end{equation*}

を得る。ここで,\omegaの時定数

(24)   \begin{equation*} T_m=\frac{RJ}{K_TK_E} \end{equation*}

は直流モータの機械的時定数と呼ばれている。上の例で計算してみるとT_m=0.01\,{\rm s}である。したがって,もし,直流モータの電気的時定数が機械的時定数に比べて十分小さい場合(経験則はT_e<0.1\,T_m)は,式(17)と式(23)を合わせて,つぎの状態方程式をもつ2次系としてよい。

(25)   \begin{equation*} \underbrace{ \left[\begin{array}{c} \dot{\theta}(t) \\ \dot{\omega}(t) \end{array}\right] }_{\dot{x}(t)}= \underbrace{ \left[\begin{array}{cc} 0 & 1 \\ 0 & -\frac{1}{T_m} \end{array}\right] }_A \underbrace{ \left[\begin{array}{c} \theta(t) \\ \omega(t) \end{array}\right] }_{x(t)}+ \underbrace{ \left[\begin{array}{cc} 0 & 0 \\ \frac{1}{T_mK_E} & -\frac{1}{J} \end{array}\right] }_B \underbrace{ \left[\begin{array}{c} v(t)\\ \tau_L(t) \end{array}\right] }_{u(t)} \end{equation*}

式(19)と比較すると,状態空間表現の次数を1だけ減らしたことになる。
これは,モデルの低次元化の一例である。

低次元化の過程を図1.12~図1.14に示しておく。


図1.12 式(1.19)に基づく低次元化前のブロック線図


図1.13 式(1.22)を用いた低次元化中のブロック線図


図1.14 式(1.22)を用いた低次元化中のブロック線図

*式(18)は,式(19)のように物理パラメータどうしの演算を含まず,それらの変動の影響を考察するのに便利な形式であり,ディスクリプタ形式の状態方程式と呼ばれる。
**ここでは,2.1.3項で学ぶ時定数の知識を前提にしている。

1.2 状態空間表現へのモデリング

前節で述べた内容は,モデリング(modeling)という枠組でとらえられる。これは,何らかの対象に対して,目的に応じた科学的分析を通してそのモデルを得ようとするものである。私たちが扱う対象は,入力と出力をもつ動的システム^*(dynamical system)であり,そのモデリングを行うのは動的な振る舞いの制御のためであるので,制御対象(controlled object)と呼ばれる。ここでの科学的分析とは,おもに物理法則から線形定係数常微分方程式^{**}で表される運動方程式や回路方程式の導出と,これに含まれるパラメータ値の決定をいう。そして,モデルの一つとして,状態空間表現(state-space description)という数理モデル^{***}(mathematical model)を考え,その構築例を示してきた。

本書では,状態空間表現を次式のように書く。

(26)   \begin{equation*} \dot{x}(t)=Ax(t)+Bu(t) \end{equation*}

(27)   \begin{equation*} y(t)=Cx(t)+Du(t) \end{equation*}

ここで,x(t)n個の状態変数(state variables)をもつベクトル,u(t)m個の入力変数(input variables)をもつベクトル,y(t)p個の出力変数(output variables)をもつベクトルである。また,ABCDは,それぞれサイズn\times nn\times mp\times np\times mの実行列であり,それぞれA行列,\mbox{\boldmathB}行列C行列D行列と呼ばれる。状態方程式(state equation)(26)と出力方程式}(output equation)(27)からなる状態空間表現を,m入力p出力をもつn次系(nth-order system)と参照する。Du(t)の項は直達項(direct part)を表し,存在しないこともある。特に,1次系のときは,\dot{x}(t)=ax(t)+bu(t)などと書く。

*動的システムは,微分方程式・差分方程式のどちらで記述されるかによって連続時間系・離散時間系,重ね合わせの原理が成り立つか否かによって線形系・非線形系,常微分方程式か偏微分方程式かによって集中定数系・分布定数系,係数パラメータの時間依存性によって時変系・時不変系,入出力が確率過程であるか否かによって決定系・確率系などに分類される。
**非線形系の場合の取り扱いは7章で述べる。1~6章までは線形時不変系のみを扱う。
***他の数理モデルとして伝達関数表現がある。状態空間表現と伝達関数表現の間の相互関係については8章で述べる。
****他のアプローチとして,入力と出力の時系列データからモデリングを行うシステム同定がある。

1.3 状態空間表現の座標変換

状態空間表現を見やすくする一つの手段として,座標変換(coordinate transformation)があるので,これについて説明しよう。
いま,n次系

(28)   \begin{equation*} \dot{x}(t)=Ax(t)+Bu(t) \end{equation*}

(29)   \begin{equation*} y(t)=Cx(t) \end{equation*}

に対して,つぎの座標変換を行いたい。

(30)   \begin{equation*} x'(t)=Tx(t) \end{equation*}

ただし,Tは正則とする。式(30)を式(28)に代入すると

(31)   \begin{equation*} \dot{x}(t)=T^{-1}\dot{x'}(t) \end{equation*}

に注意して

(32)   \begin{equation*} T^{-1}\dot{x'}(t)=AT^{-1}x'(t)+Bu(t) \end{equation*}

%すなわち

(33)   \begin{equation*} \dot{x'}(t)=\underbrace{TAT^{-1}}_{A'}x'(t)+\underbrace{TB}_{B'}u(t) \end{equation*}

となる。また,式(30)を式(29)に代入すると

(34)   \begin{equation*} y(t)=\underbrace{CT^{-1}}_{C'}x'(t) \end{equation*}

となる。この結果を,参照しやすいようにつぎにまとめておく。

定理1.1 n次系\dot{x}(t)=Ax(t)+Bu(t),\ y(t)=Cx(t)に対して,座標変換x'(t)=Tx(t)を行うと,新しいn次系は次式で表される。

(35)   \begin{equation*} \dot{x'}(t)=A'x'(t)+B'u(t) \end{equation*}

(36)   \begin{equation*} y(t)=C'x'(t) \end{equation*}

ただし

(37)   \begin{equation*} A'=TAT^{-1},\ B'=TB,\ C'=CT^{-1} \end{equation*}

例題1.1 直流モータの状態方程式(25)において,\tau_Lを零とおくと

(38)   \begin{equation*} \underbrace{ \left[\begin{array}{c} \dot{\theta}(t)\\ \dot{\omega}(t) \end{array}\right] }_{\dot{x}(t)}= \underbrace{ \left[\begin{array}{cc} 0 & 1 \\ 0 & -\frac{1}{T_m} \end{array}\right] }_A \underbrace{ \left[\begin{array}{c} \theta(t)\\ \omega(t) \end{array}\right] }_{x(t)} + \underbrace{ \left[\begin{array}{c} 0 \\ \frac{1}{T_mK_E} \end{array}\right] }_B \underbrace{ v(t) }_{u(t)} % \end{equation*}

である。これに対して,座標変換

(39)   \begin{equation*} \underbrace{ \left[\begin{array}{c} x'_1(t) \\ x'_2(t) \end{array}\right] }_{x'(t)} = \underbrace{ \left[\begin{array}{cc} 1 & T_m \\ 0 & -T_m \end{array}\right] }_T \underbrace{ \left[\begin{array}{c} \theta(t) \\ \omega(t) \end{array}\right] }_{x(t)} % \end{equation*}

を行うと,新しい状態方程式は

(40)   \begin{equation*} \underbrace{ \left[\begin{array}{c} \dot{x}'_1(t) \\ \dot{x}'_2(t) \end{array}\right] }_{\dot{x}'(t)}= \underbrace{ \left[\begin{array}{cc} 0 & 0 \\ 0 & -\frac{1}{T_m} \end{array}\right] }_{A'} \underbrace{ \left[\begin{array}{c} x'_1(t) \\ x'_2(t) \end{array}\right] }_{x'(t)} + \underbrace{ \left[\begin{array}{c} \frac{1}{K_E} \\ -\frac{1}{K_E} \end{array}\right] }_{B'} \underbrace{ v(t) }_{u(t)} \end{equation*}

となることを示しなさい。

解答 座標変換後のA行列とB行列は,定理1.1を用いて

(41)   \begin{equation*} \underbrace{ \left[\begin{array}{cc} 1 & T_m \\ 0 & -T_m \end{array}\right] \left[\begin{array}{cc} 0 & 1 \\ 0 & -\frac{1}{T_m} \end{array}\right] \left[\begin{array}{cc} 1 & 1 \\ 0 & -\frac{1}{T_m} \end{array}\right] }_{TAT^{-1}} = \underbrace{ \left[\begin{array}{cc} 0 & 0 \\ 0 & -\frac{1}{T_m} \end{array}\right] }_{A'} \end{equation*}

(42)   \begin{equation*} \underbrace{ \left[\begin{array}{cc} 1 & T_m \\ 0 & -T_m \end{array}\right] \left[\begin{array}{cc} 0 \\ \frac{1}{T_mK_E} \end{array}\right] }_{TB} = \underbrace{ \left[\begin{array}{c} \frac{1}{K_E} \\ -\frac{1}{K_E} \end{array}\right] }_{B'} \end{equation*}

のように得られる。

ここで,2次系の状態方程式が,二つの1次系の状態方程式

(43)   \begin{equation*} \dot{x}'_1(t)=\frac{1}{K_E}v(t),\quad \dot{x}'_2(t)=-\frac{1}{T_m}x'_2(t)-\frac{1}{K_E}v(t) \end{equation*}

に分離されており,入力から状態変数への影響の考察をしやすくなっていることに注意してほしい。

1.4 状態空間表現の直列結合

制御対象の状態空間表現を求める際に,図1.15に示すように,二つの部分システムの状態空間表現を求めておいて,これらを直列結合(serial connection)する場合がある。このときの結合システムの状態空間表現を求めることを考える。


図1.15 直列結合(u_1=y_2)

まず,その結果を定理の形で示そう。

定理1.2 二つの状態空間表現

(44)   \begin{equation*} \dot x_1(t)=A_1x_1(t)+B_1u_1(t) \\ \end{equation*}

(45)   \begin{equation*} y_1(t)=C_1x_1(t)+D_1u_1(t) \end{equation*}

および

(46)   \begin{equation*} \dot x_2(t)=A_2x_2(t)+B_2u_2(t) \end{equation*}

(47)   \begin{equation*} y_2(t)=C_2x_2(t)+D_2u_2(t) % \end{equation*}

に対して,u_1=y_2のように直列結合した場合の状態空間表現は

(48)   \begin{equation*} \left[\begin{array}{c} \dot{x}_1(t) \\ \dot{x}_2(t) \end{array}\right]= \left[\begin{array}{cc} A_1 & B_1C_2 \\ 0 & A_2 \end{array}\right] \left[\begin{array}{c} x_1(t) \\ x_2(t) \end{array}\right] + \left[\begin{array}{cc} B_1D_2 \\ B_2 \end{array}\right] u_2(t) \end{equation*}

(49)   \begin{equation*} y_1(t)= \left[\begin{array}{cc} C_1 & D_1C_2 \end{array}\right] \left[\begin{array}{c} x_1(t) \\ x_2(t) \end{array}\right] +D_1D_2u_2(t) \end{equation*}

のように得られる。

証明 \dot x_1(t)=A_1x_1(t)+B_1u_1(t)y_1(t)=C_1x_1(t)+D_1u_1(t)に,u_1(t)=y_2(t)=C_2x_2(t)+D_2u_2(t)を代入して

(50)   \begin{equation*} \dot x_1(t)=A_1x_1(t)+B_1(C_2x_2(t)+D_2u_2(t)) \end{equation*}

(51)   \begin{equation*} y_1(t)=C_1x_1(t)+D_1(C_2x_2(t)+D_2u_2(t)) \end{equation*}

となる。第1式と\dot x_2(t)=A_2x_2(t)+B_2u_2(t)をまとめたものと,第2式から,定理の結果を得る。

例題1.2 2次系の制御対象

(52)   \begin{equation*} \dot x(t)= \left[\begin{array}{cc} 0 & 1 \\ -\omega_n^2 & -2\zeta\omega_n \end{array}\right] x(t)+ \left[\begin{array}{cc} 0 \\ K\omega_n^2 \end{array}\right] u(t) \end{equation*}

(53)   \begin{equation*} y(t)=\left[\begin{array}{cc} 1 & 0 \end{array}\right]x(t) \end{equation*}

に対して(x(t)は2次元ベクトル),1次系のアクチュエータ

(54)   \begin{equation*} \displaystyle{\dot{x}_a(t) =-\frac{1}{T_a}x_a(t) +\frac{K_a}{T_a}u_a(t)} \end{equation*}

(55)   \begin{equation*} y_a(t)=x_a(t) \end{equation*}

を,u=y_aのように直列結合した場合の状態空間表現を求めなさい。

解答 定理1.2を用いて,直列結合の状態空間表現として

(56)   \begin{equation*} \left[\begin{array}{c} \dot{x}(t) \\ \dot{x}_a(t) \end{array}\right]= \left[\begin{array}{cc|c} 0 & 1 & 0\\ -\omega_n^2 & -2\zeta\omega_n & K\omega_n^2 \\\hline 0 & 0 & -\frac{1}{T_a} \end{array}\right] \left[\begin{array}{c} x(t) \\ x_a(t) \end{array}\right]+ \left[\begin{array}{cc} 0 \\ 0 \\\hline \frac{K_a}{T_a} \end{array}\right] u_a(t) \end{equation*}

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

が得られる^*

問1.4 例題1.2の直列結合の状態空間表現を,状態ベクトルが[x_a(t)\ x(t)]^Tとなるように求めなさい。

*ここで,A行列の縦線と横線,B行列の横線は,状態ベクトルの要素xx_aのサイズに適合するように引かれている。


演習問題


【1】 いろいろな計測装置の基礎となる電気回路の一つにブリッジ回路がある。
例えば,図1.16に示すブリッジ回路^*を考えてみよう。この回路方程式は

(58)   \begin{equation*} L\dot{i}(t)=-\frac{R_1+R_3}{R_1+R_2}R_2i(t)-\frac{R_1+R_3}{R_1+R_2}v(t) +u(t)-R_3C\dot{v}(t) \end{equation*}

(59)   \begin{equation*} C\dot{v}(t)=-\left(\frac{R_2}{R_1+R_2}-\frac{R_4}{R_3+R_4}\right)i(t) -\left(\frac{1}{R_1+R_2}+\frac{1}{R_3+R_4}\right)v(t) \end{equation*}

で与えられる。いま,ブリッジ条件

(60)   \begin{equation*} R_1R_4=R_2R_3 \end{equation*}

が成り立つとして,つぎの状態方程式を導出しなさい。

(61)   \begin{equation*} \underbrace{ \left[\begin{array}{c} \dot{i}(t) \\ \dot{v}(t) \end{array}\right] }_{\dot{x}(t)}= \underbrace{ \left[\begin{array}{cc} -\frac{1}{L}(\frac{R_1R_2}{R_1+R_2}+\frac{R_3R_4}{R_3+R_4})\ \ 0 \\ 0\ \ -\frac{1}{C}(\frac{1}{R_1+R_2}+\frac{1}{R_3+R_4}) \end{array}\right] }_A \underbrace{ \left[\begin{array}{c} i(t) \\ v(t) \end{array}\right] }_{x(t)}+ \underbrace{ \left[\begin{array}{c} \frac{1}{L} \\ 0 \end{array}\right] }_B u(t) \end{equation*}

この状態方程式に基づいて,平衡ブリッジ回路のブロック線図を描きなさい。


図1.16 ブリッジ回路

【2】 さまざまな柔軟構造物の制振問題は,重要な制御のテーマである。
その特徴は,図1.17に示す連結台車^{**}にもみられる。この運動方程式は

(62)   \begin{equation*} m_1\ddot{x}_1(t)=-k(x_1(t)-x_2(t))+f_1(t)+w_1(t) \end{equation*}

(63)   \begin{equation*} m_2\ddot{x}_2(t)=-k(x_2(t)-x_1(t))+w_2(t) \end{equation*}

で与えられる。ここで,m_1m_2はそれぞれ台車1と台車2の質量,kはばね定数である。このとき,つぎの状態方程式を導出しなさい。

(64)   \begin{equation*} \begin{array}{rcl} \underbrace{ \left[\begin{array}{c} \dot{x_1}(t) \\ \dot{x_2}(t) \\ \dot{v_1}(t) \\ \dot{v_2}(t) \end{array}\right] }_{\dot{x}(t)} &=& \underbrace{ \left[\begin{array}{cccc} 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \\ -\frac{k}{m_1} & \frac{k}{m_1} & 0 & 0 \\ \frac{k}{m_2} & -\frac{k}{m_2} & 0 & 0 \end{array}\right] }_A \underbrace{ \left[\begin{array}{c} x_1(t) \\ x_2(t) \\ v_1(t) \\ v_2(t) \end{array}\right] }_{x(t)} \\ & & \hspace*{-1zw} + \underbrace{ \left[\begin{array}{ccc} 0 & 0 & 0 \\ 0 & 0 & 0 \\ \frac{1}{m_1} & 0 & \frac{1}{m_1} \\ 0 & \frac{1}{m_2} & 0 \end{array}\right] }_B \underbrace{ \left[\begin{array}{c} w_1(t)\\ w_2(t)\\ f_1(t) \end{array}\right] }_{u(t)} \end{array} \end{equation*}

この状態方程式に基づいて,連結台車のブロック線図を描きなさい。


図1.17 連結台車

【3】 式23で表される直流モータにおいて,一定入力v^*,一定負荷\tau_L^*のもとで,一定角速度\omega^*の平衡状態が達成されているものとする。この平衡状態を基準とする直流モータの時間的振る舞いを表す状態方程式を示しなさい。

【4】 本書におけるすべての数値計算は,対話型の行列計算環境である学生版MATLAB^{***}を用いて行っている。また,すべての時間応答のグラフは,(非線形)微分方程式による対話型シミュレーション環境である学生版SIMULINK^{****}を用いて得ている。時間応答のシミュレーションのためには,状態方程式のブロック線図を描くことが必要となる。例えば,心臓のペースメーカのブロック線図(図1.3)を得たとすると,SIMULINKでは,これを図1.18のようにほぼそのままの構成で,対話型操作により表現する。ブロックIntegratorの初期値とブロックGainの値を設定し,微分方程式のソルバーの種類,サンプリング周期,シミュレーション時間などを設定すれば,ブロックScopeに図1.1の時間応答を直ちにみることができる。時系列データの処理やグラフ化はMATLABで行える。

MATLABとSIMULINKが手元にあれば,シミュレーション1.3と同一条件下で,直流モータの低次元化後の状態方程式25による角速度の応答を,低次元化前の状態方程式19によるものと比較しなさい。

図1.18 SIMULINKによる微分方程式のブロック表現

*高橋・有本:回路網とシステム理論,コロナ社 (1974)のpp.65\sim66から引用。
**B.Wie, D.2.Bernstein : Benchmark Problems for Robust Control Design, ACC Proc. pp.2047\sim2048 (1992) から引用。
***The Student Edition of MATLAB-Version\,5 User’s Guide, Prentice Hall (1997)
****The Student Edition of SIMULINK-Version\,2 User’s Guide, Prentice Hall (1998)