離散時間系の差分方程式の導出


状態変数モデル

システムを数式化する際に内部状態を使用した数学モデルを状態変数モデルといい、入力、内部状態および出力を表す変数を、それぞれ入力変数、状態変数、出力変数と呼ぶ。

出力変数は状態変数から導出されるため、以下では入力変数、状態変数で表される差分方程式を導出する。


1次系の連続時間系状態方程式、$t\geq0$

\[ \dot{x}(t) = ax(t) + bu(t) \ \] において、初期値 $ x(0) = x_0 $として両辺をラプラス変換すると \[ \ sX(s)-x_0 = aX(s)+bU(s) \ \] となり、これから $ X(s) $ を求めると \[ \ X(s) = (s-a)^{-1}x_0+(s-a)^{-1}bU(s) \ \] が得られる。次に、ラプラス逆変換して $ x(t) $ を求めると \[ \ x(t) = e^{at}x_0+ \displaystyle \int_0^t e^{a(t-\tau)}bu(\tau)d\tau \ \] a、b を行列 $ A $、$ B $ で表された状態方程式 \[ \dot{x}(t) = Ax(t) + Bu(t) \ \] における $ x(t) $ はこれらの過程から類推して、同じ形式で書ける。 \[ \ x(t) = e^{At}x_0+ \displaystyle \int_0^t e^{A(t-\tau)}Bu(\tau)d\tau \ \] ここで、行列 $ A $ に対する指数関数表現は \[ \ e^{At} = I+At+\frac{A^2}{2!}t^2+・・・ \ \] で定義される。$ e^{At} $ は $ x_0 $ から $ t $ 秒後に $ x(t) $ の状態に遷移する道順を決める為、状態遷移行列と呼ばれ $ \Phi(t) $ と記述される。


次に、$ \Phi(t) $ をラプラス逆変換から求めてみる。 状態方程式において、$ u(t) = 0 $ とおいた \[ \dot{x}(t) = Ax(t) \ \] の両辺をラプラス変換すると \[ \ sX(s)-x_0 = AX(s) \ \] \[ \ (sI-A)X(s) = x_0 \ \] となり、両辺に $ (sI-A)^{-1} $ を左からかけると \[ \ X(s) = (sI-A)^{-1}x_0 \ \] が得られる。これを上記の $ x(t) $ の式と比較すると \[ \ \Phi(t) = e^{At} = \mathcal{L}^{-1}[(sI-A)^{-1}] \ \] であることがわかる。


状態方程式の離散化

入力信号 $ u(t) $、状態変数 $ x(t) $ をサンプリング周期 $ T $ 毎に観測するものとすれば、$ k $ 番目のサンプル値を $ u(kT) $、$ x(kT) $ と書くと、離散時間系の状態方程式は \[ \ x((k+1)T) = Px(kT)+Q(kT) \ \] の差分方程式で表せる。 ここで、$ P $、$ Q $ は次のように求めることができる。 $ t = kT $ の状態 $ x(kT) $ を初期値としてサンプリング周期 $ T $ の間の状態の変化を求めると \[ \ x((k+1)T) = \Phi(t)x(kT)+ \displaystyle \int_{kT}^{(k+1)T} \Phi((k+1)T-\tau)Bu(k)d\tau \ \ = \Phi(t)x(kT)+ \displaystyle \left( \int_0^T e^{A\eta}d\eta \right) Bu(kT) \ \ = \Phi(t)x(kT)+A^{-1} \displaystyle \left[ \Phi(T)-I \right] Bu(kT) \ \]


上記より離散時間系の状態方程式, $k=0,1,\cdot\cdot\cdot$

\[ x((k+1)T) = \Phi(T)x(kT) + \Psi(T)u(kT) \ \]

但し、$ \Phi(T) = \mathcal{L}^{-1}[(sI-A)^{-1}] $、$ \Psi(T) = A^{-1}( \Phi(T)-I )B $、$ I = \left(\begin{array}{ccc} 1 & 0 \\ 0 & 1 \end{array}\right) $

上記の行列係数を計算することにより、サンプリング時間 $T$ の離散時間系の差分方程式が得られる。


補足

行列 $ A $ に対する指数関数表現 \[ \ e^{At} = I+At+\frac{A^2}{2!}t^2+・・・ \ \] 両辺を微分すると \[ \ \frac{d}{dt}e^{At} = A+A^2t+\frac{A^3}{2!}t^2+・・・ = A\displaystyle \left( I+At+\frac{A^2}{2!}t^2+・・・ \right) = Ae^{At} \ \] 上式を変形して \[ \ e^{At} = A^{-1}\frac{d}{dt}e^{At} \ \] 両辺を区間 0 から T で積分すると \[ \ \displaystyle \int_0^T e^{At}dt = A^{-1} \left( e^{AT} - I \right) = A^{-1} \left[ \Phi(T)-I \right] \ \]


状態変数モデルの参考書


ホーム

2022-09-12