状態方程式から離散化した逐次式へ

状態方程式を離散化した逐次式に変換

状態方程式を離散化した逐次式に変換する手法を説明しますが、途中に出てくる微分方程式の一般解やPade近似などは参考知識なので、具体的な手法は最後の PとQの計算式 になります。

状態方程式の解を求める

微分式のままではシステムの振る舞いを具体的に計算することが難しいので、状態方程式の解からシステムの時間応答を計算できるようにします。

状態方程式は連続時間の微分方程式
$$\frac{d x(t)}{d t}=A x(t)+B u(t) \\$$ 微分方程式の一般解は
$$\boldsymbol{x}(t)=e^{A t} x_{0}+\int_{0}^{t} e^{A(t-\tau)} B u(\tau) d \tau \\$$

離散式の導出

この一般解を使って離散式を導き出します。
時刻tにおける値を𝒙(𝑘) としたとき、Δt 後の t+Δt の状態 𝒙(𝑘+1) は

$$\boldsymbol{x}(k+1)=e^{A \Delta t} \boldsymbol{x}(k)+A^{-1}\left(e^{A \Delta t}-I\right) B \boldsymbol{u}(k)$$ ここで　$P=e^{A \Delta t} \\$　$Q=A^{-1}\left(e^{A \Delta t}-I\right) B$　とすると
$$x(t+\Delta t)=P x(t)+Q u(t)$$
になります。
線形システムで、状態方程式のA,Bで表すシステムのパラメータが一定なら、離散時間 Δt を定めれば、PとQを計算することで離散式に変換できます。

$e^{A \Delta t}$の計算方法

一般にAは行列なので$e^{A \Delta t}$の直接計算は難しいため、このブログでは Pade近似 を使います。
Pade近似式でPとQを計算すれば状態方程式は逐次式に変換することができます。

Pade近似 $$e^{x}=\frac{1+x / 2}{1-x / 2}$$ を適用すると、P Q は以下の式になります。
$$P=e^{A \Delta t}=\left(I+\frac{\Delta t}{2} A\right) \quad\left(I-\frac{\Delta t}{2} A\right)^{-1}=I+A \Delta t\left(I-\frac{\Delta t}{2} A\right)^{-1}$$ $$Q=A^{-1}\left(e^{A \Delta t}-I\right) B=\Delta t\left(I-\frac{\Delta t}{2} A\right)^{-1}B$$

Excelには行列演算の関数 MINVERSE や MMULT があるので、PとQはExcelで容易に計算できます。計算方法は次のCR回路モデルの時間応答の確認で具体的なモデルを例に説明していきます。

最後までお読みいただき、ありがとうございました。

状態方程式を離散化する方法には、いくつかのアプローチがありますが、その精度の比較については、こちら→ リンク が参考になります。