線形制御系設計

    線形システムの表現
    線形システムの特性解析
    離散時間システム
    制御器設計
    追従制御 (FF)
    外乱抑圧制御 (FB)
    観測器設計

    モーションコントロール

    加速度制御
    剛体の姿勢制御

    モータドライブ

    駆動制御法
    磁極位置推定法

    システム同定

    ホワイトボックスモデル
    パラメトリックモデル

安定逆系の設計

不安定零点を持つシステムの逆系は不安定となるため実装ができないが,未来の情報(非因果的情報)を利用することで安定的に入出力の逆計算が可能となる。 これは追従制御問題のような予め参照軌道が与えられる場合に,不安定零点の影響を考慮した入力決定を可能とする。

伝達関数表現を用いた安定逆系

安定逆系を実現する方法について,初期の段階では状態空間表現を使用して議論が行われていた[1-3]。十河氏は周波数領域で設計されたシステムに対して安定逆系の設計を行う方法の整備が不十分であることを指摘し,伝達関数表現に対して安定逆系を設計する方法を与えた[4]。同氏は文献[5]にて両側ラプラス変換[6]を非因果的な畳み込み積分を表現するための数学ツールとして利用して,文献[4]では両側ラプラス変換を用いた手法が従来の状態空間表現を利用した手法と等価であることを示し,直感的かつ直接的な安定逆系の設計を可能とした。

両側ラプラス変換

両側ラプラス変換は次のように定義される。
L[g(t)](s)=G(s)g(t)estdt\begin{align} \mathcal{L}[g(t)](s) = G(s) \coloneqq \int^{\infty}_{-\infty} g(t)e^{-st}dt \end{align}
収束域を垂直帯状領域 γ1<Re(s)<γ2\gamma_1 < {\rm Re}(s) < \gamma_2 とする。また,逆変換は次のように記述される。
L1[G(s)](t)={12πjClestG(s)dsif t012πjCrestG(s)dsif t<0={+nRes(estG(s),pn)if t0mRes(estG(s),pm)if t<0\begin{align} \mathcal{L}^{-1}[G(s)](t) &= \left\{ \begin{matrix} \frac{1}{2\pi j}\int_{C_{\rm l}} e^{st} G(s) ds & {\rm if}\ t\geq 0 \\ \frac{1}{2\pi j}\int_{C_{\rm r}} e^{st} G(s) ds & {\rm if}\ t < 0 \end{matrix} \right. \\ &= \left\{ \begin{matrix} +\sum_{\rm n} {\rm Res}(e^{st}G(s), p_{\rm n}) & {\rm if}\ t\geq 0 \\ -\sum_{\rm m} {\rm Res}(e^{st}G(s), p_{\rm m}) & {\rm if}\ t < 0 \end{matrix} \right. \end{align}
ここで,pn,pmp_{\rm n}, p_{\rm m}G(s)G(s) の左半平面の極および右半平面の極を表し,Cl,CrC_{\rm l}, C_{\rm r} はそれぞれ虚軸上の線分 (jR,+jR)(-jR, +jR) に原点を中心とする半径 RR の円の左半分および右半分を付け加えた左回りの閉曲線と右回りの閉曲線とする。RR は十分に大きく,閉曲線が全ての極を含むものとする。ここで,左半平面にのみ極を持つシステムを Gs(s)G_{\rm s}(s),右半平面にのみにのみ極を持つシステムを Gu(s)G_{\rm u}(s) とすれば,以下が成立する[5]
L1[Gs(s)](t)={+nRes(estGs(s),pn)if t00if t<0=Lˉ1[Gs(s)](t)gs(t)L1[Gu(s)](t)={+nRes(estGu(s),pn)if t00if t<0=Lˉ1[Gu(s)](t)gˉu(t)L1[Gu(s)](t)={0if t012πjCrestGu(s)dsif t<0={0if t012πjClestGu(s)dsif t<0=gˉu(t)\begin{align} \mathcal{L}^{-1}[G_{\rm s}(s)](t) &= \left\{ \begin{matrix} +\sum_{\rm n} {\rm Res}(e^{st}G_{\rm s}(s), p_{\rm n}) & {\rm if}\ t\geq 0 \\ 0 & {\rm if}\ t < 0 \end{matrix} \right.\\ &= \bar{\mathcal{L}}^{-1}[G_{\rm s}(s)](t) \equiv g_{\rm s}(t)\\ \mathcal{L}^{-1}[G_{\rm u}(-s)](t) &= \left\{ \begin{matrix} +\sum_{\rm n} {\rm Res}(e^{st}G_{\rm u}(-s), -p_{\rm n}) & {\rm if}\ t\geq 0 \\ 0 & {\rm if}\ t < 0 \end{matrix} \right.\\ &= \bar{\mathcal{L}}^{-1}[G_{\rm u}(-s)](t) \equiv \bar{g}_{\rm u}(t) \\ \mathcal{L}^{-1}[G_{\rm u}(s)](t) &= \left\{ \begin{matrix} 0 & {\rm if}\ t\geq 0 \\ \frac{1}{2\pi j}\int_{C_{\rm r}} e^{st} G_{\rm u}(s) ds & {\rm if}\ t < 0 \end{matrix} \right.\\ &= \left\{ \begin{matrix} 0 & {\rm if}\ t\geq 0 \\ \frac{1}{2\pi j}\int_{-C_{\rm l}} e^{-s't} G_{\rm u}(-s') ds' & {\rm if}\ t < 0 \end{matrix} \right.\\ &= \bar{g}_{\rm u}(-t) \end{align}
ただし,Lˉ\bar{\mathcal{L}} は片側ラプラス変換を表す。また,G(s)=Gs(s)+Gu(s)G(s)=G_{\rm s}(s)+G_{\rm u}(s) のように記述すれば,その逆両側ラプラス変換は次のようになる。
g(t)L1[G(s)](t)={gs(t)if t0gˉu(t)if t<0\begin{align} g(t) &\equiv \mathcal{L}^{-1}[G(s)](t) &= \left\{ \begin{matrix} g_{\rm s}(t) & {\rm if}\ t\geq 0 \\ \bar{g}_{\rm u}(-t) & {\rm if}\ t < 0 \end{matrix} \right.\\ \end{align}
ここで,Yd(s)=G(s)Ud(s)Y_{\rm d}(s)=G(s)U_{\rm d}(s) と記述されるシステムがあれば,後述の実現条件の下で ud(s)u_{\rm d}(s) を求めることができる[5]
ud(t)=L1[G(s)Yd(s)]=g(tτ)yd(τ)dτ=tgs(tτ)yd(τ)dτ+tgˉu(tτ)yd(τ)dτ=tgs(tτ)yd(τ)dτ+tˉgˉu((tˉτˉ))yd(τˉ)dτˉtˉ=t\begin{align} u_{\rm d}(t) &= \mathcal{L}^{-1}\left[ G(s) \cdot Y_{\rm d}(s) \right] \\ &= \int^{\infty}_{-\infty} g(t-\tau)y_{\rm d}(\tau) d\tau \\ &= \int^{t}_{-\infty} g_{\rm s}(t-\tau)y_{\rm d}(\tau) d\tau + \int^{\infty}_{t} \bar{g}_{\rm u}(t-\tau)y_{\rm d}(\tau) d\tau \\ &= \int^{t}_{-\infty} g_{\rm s}(t-\tau)y_{\rm d}(\tau) d\tau + \left. \int^{\bar{t}}_{-\infty} \bar{g}_{\rm u}( - (\bar{t} - \bar{\tau}))y_{\rm d}(-\bar{\tau}) d\bar{\tau} \right|_{\bar{t}=-t} \\ \end{align}
一方で,システムが G(s)=Gs(s)Gu(s)G(s)=G_{\rm s}(s)\cdot G_{\rm u}(s) のように記述される場合には, ud(s)u_{\rm d}(s) は以下のように求めることができる[4]
uu(t)=L1[Gs(s)Gu(s)Yd(s)]=tgs(tτ)v(τ)dτv(t)σgˉu(στ)yd(τ)dτσ=t\begin{align} u_{\rm u}(t) &= \mathcal{L}^{-1}\left[ G_{\rm s}(s) \cdot G_{\rm u}(s) \cdot Y_{\rm d}(s) \right] \\ &= \int^{t}_{-\infty} g_{\rm s}(t-\tau)v(\tau) d\tau \\ v(t) &\equiv \left. \int^{\sigma}_{-\infty}\bar{g}_{\rm u}(\sigma -\tau) y_{\rm d}(-\tau) d\tau \right|_{\sigma=-t} \end{align}

安定逆系の実現条件

不安定零点を持つシステムについて,状態空間表現および同値な伝達関数表現について考える。
状態空間表現
x˙(t)=Ax(t)+Bu(t)y(t)=Cx(t)\begin{align} \dot{\bm{x}}(t) &= \bm{A}\bm{x}(t) + \bm{B}u(t)\\ y(t) &= \bm{C}\bm{x}(t)\\ \end{align}
伝達空間表現
Y(s)=G(s)U(s)\begin{align} Y(s) &= G(s) U(s)\\ \end{align}
ただし,u(t), x(t), y(t)u(t),\ x(t),\ y(t) は入力,内部状態量,および出力を表し,U(s)=L[u(t)](s),Y(s)=L[y(t)](s)U(s)=\mathcal{L}[u(t)](s), Y(s)=\mathcal{L}[y(t)](s) とする。G(s)G(s)虚軸上に零点を持たない (逆系の応答が無限に振動しない) ものとし,相対次数を r\rm r とする。参照軌道 yd(t)y_{\rm d}(t) が与えられた際の内部状態量と入力の参照値をそれぞれ xd(t),u(t)\bm{x}_{\rm d}(t), u(t) と記述すると,これらが満たすべき等式は次のようになる。
x˙d(t)=Axd(t)+Bud(t)yd(t)=Cxd(t)\begin{align} \dot{\bm{x}}_{\rm d}(t) &= \bm{A}\bm{x}_{\rm d}(t) + \bm{B}u_{\rm d}(t)\\ y_{\rm d}(t) &= \bm{C}\bm{x}_{\rm d}(t)\\ \end{align}
また,内部状態および入力は以下の条件を満たすものとする。
ud(t)0,xd(t)0(t±)\begin{align} u_{\rm d}(t)\rightarrow 0, \quad \bm{x}_{\rm d}(t)\rightarrow \bm{0} \quad (t\rightarrow \pm\infty) \end{align}
ここで,i[0,r] yd(i)(t)L1L^{\forall {\rm i} \in [0, {\rm r}]\ } y^{({\rm i})}_{\rm d}(t)\in L_{1} \cap L_{\infty} (時間全体で絶対値積分可能かつ有界)であるならば,ud(t), xd(t)u_{\rm d}(t),\ x_{\rm d}(t) は次のように求めることができる[4]
ud(t)=L1[G1(s)Yd(s)]xd(t)=L1[(sIA)1BUd(s)]\begin{align} u_{\rm d}(t) &= \mathcal{L}^{-1}\left[ G^{-1}(s) \cdot Y_{\rm d}(s) \right] \\ x_{\rm d}(t) &= \mathcal{L}^{-1}\left[ (s\bm{I} - \bm{A})^{-1} \bm{B} U_{\rm d}(s) \right] \end{align}

安定逆系の簡易実装

文献[4]では逆系を安定なシステムと不安定なシステムの積として表現することで安定逆系を与えたが,解が複雑であり数値計算の負荷が大きい。システムが安定零点のみを持つ場合にはその逆系の設計は容易であることを踏まえ,逆系を安定なシステムと不安定なシステムの和として表現する文献[5]の方法が実装の観点からは望ましい。また,逆系の整理では(1) 元のシステムの相対次数と同じ次数を持つ多項式とプロパーな有理関数の和形式への分解,(2) 有理関数の安定システムと不安定システムへの分解を行う必要があるが,前者は参照軌道に対する微分可能性に関する制約を設けて段階的な逆系の設計を許容すれば回避することができる。

有理関数形式で表現される伝達関数 G(s)G(s) を持つ非最小位相システムについて考える。
Y(s)=G(s)U(s)G(s)=A(s)B(s)\begin{align} Y(s) &= G(s) U(s)\\ G(s)&=\frac{A(s)}{B(s)} \end{align}
ただし,U(s),Y(s)U(s), Y(s) は入力と出力を表す。ここに内部状態変数 X(s)X(s) を導入すれば,次のように記述される。
X(s)=A(s)U(s)Y(s)=B1(s)X(s)\begin{align} X(s) &= A(s)U(s)\\ Y(s) &= B^{-1}(s) X(s) \end{align}
追従制御問題を想定して与えられた参照軌道 yd(t)y_{\rm d}(t) から入力 u(t)u(t) を求めることを考えれば,安定逆系を設計する対象は B1(s),A(s)B^{-1}(s), A(s) となる。B1(s)B^{-1}(s) の逆系 B(s)B(s)ss の多項式となり時間領域では微分を含む作用となるが,B(s)B(s) の次数を rr として i[0,r] yd(i)(0)=0, yd(i)(t)L1L^{\forall {\rm i} \in [0, {\rm r}]\ } y^{({\rm i})}_{\rm d}(0)=0, \ y^{({\rm i})}_{\rm d}(t)\in L_{1} \cap L_{\infty} となるように選択すれば実装が可能となる。B(s)B(s) が実係数多項式 k=0rbksk\sum_{{\rm k}=0}^{\rm r} b_{\rm k}s^{\rm k} として与えられる場合,参照軌道 yd(t)y_{\rm d}(t) が与えられた際の内部状態参照値を xd(t)x_{\rm d}(t) とすれば次のように表現される。
xd(t)=L1[k=0rbkskYd(s)]=k=0rbkyd(k)(t)\begin{align} x_{\rm d}(t) = \mathcal{L}^{-1}\left[\sum_{{\rm k}=0}^{\rm r} b_{\rm k}s^{\rm k} Y_{\rm d}(s)\right] = \sum_{{\rm k}=0}^{\rm r} b_{\rm k}y^{({\rm k})}_{\rm d}(t) \end{align}
一方で,A(s)A(s) は不安定零点を含むため逆系 A1(s)A^{-1}(s) をそのまま安定逆系として扱うことはできない。そこで,A1(s)A^{-1}(s) を安定極を持つシステム As(s)A^{\ast}_{\rm s}(s) と不安定極を持つシステム Au(s)A^{\ast}_{\rm u}(s) に分離し,両側ラプラス変換が簡単な形式で表現する。
A1(s)=As(s)+Au(s)\begin{align} A^{-1}(s) &= A^{\ast}_{\rm s}(s) + A^{\ast}_{\rm u}(s) \end{align}
以上より,入力 u(t)u(t) は次のように決定される。
xd(t)=k=0rbkyd(k)(t)u(t)=0tas(tτ)xd(τ)dτ+σaˉu((στ))xd(τ)dτσ=t\begin{align} x_{\rm d}(t) &= \sum_{{\rm k}=0}^{\rm r} b_{\rm k}y^{({\rm k})}_{\rm d}(t)\\ u(t) &= \int^{t}_{0}a_{\rm s}(t -\tau) x_{\rm d}(\tau) d\tau + \left. \int^{\sigma}_{-\infty}\bar{a}_{\rm u}(-(\sigma - \tau)) x_{\rm d}(-\tau) d\tau \right|_{\sigma=-t}\\ \end{align}
ただし,as(t)L1[As(s)](t), aˉu(t)L1[Au(s)](t)a_{\rm s}(t)\equiv \mathcal{L}^{-1}[A^{\ast}_{\rm s}(s)](t),\ \bar{a}_{\rm u}(t)\equiv \mathcal{L}^{-1}[A^{\ast}_{\rm u}(-s)](t) とする。

非因果的入力の簡易導出

解の存在を仮定して議論を進めれば,不安定極のみを持つシステムに対する非因果的入力は簡単に導出することができる。 以下の議論では,非因果的入力を生成するための制約として次の条件を設ける。
limtx(t)=0\begin{align} \lim_{t\rightarrow\infty} \bm{x}(t) = \bm{0} \end{align}
状態空間表現
次の状態空間表現について考える。
x˙(t)=Ax(t)+Bu(t)\begin{align} \dot{\bm{x}}(t) &= \bm{A}\bm{x}(t) + \bm{B}\bm{u}(t) \end{align}
ただし,A\bm{A} は全ての固有値が負のシステム行列,B\bm{B} は入力行列を表す。このシステムの状態 x(t)\bm{x}(t) が有界であると仮定すれば,上記の微分方程式の解は次のように記述される。
x(t)=eA(tt0)x(t0)+t0teA(tτ)Bu(τ)dτ\begin{align} \bm{x}(t) = e^{\bm{A}(t-t_{0})}\bm{x}(t_{0}) + \int^{t}_{t_{0}} e^{\bm{A}(t-\tau)}\bm{B}\bm{u}(\tau) d\tau \end{align}
ただし,t0t_{0} は任意の時刻とする。無限遠の時刻を起点として解を整理すれば,以下の得る。
x(t)=limt0(eA(tt0)x(t0)+t0teA(tτ)Bu(τ)dτ)=teA(tτ)Bu(τ)dτ=teA(tτ)Bu(τ)dτ(limγeAγx(γ)=0)\begin{align} \bm{x}(t) &= \lim_{t_{0}\rightarrow\infty} \left(e^{\bm{A}(t-t_{0})}\bm{x}(t_{0}) + \int^{t}_{t_{0}} e^{\bm{A}(t-\tau)}\bm{B}\bm{u}(\tau) d\tau \right)\\ &= \int^{t}_{\infty} e^{\bm{A}(t-\tau)}\bm{B}\bm{u}(\tau) d\tau = - \int^{\infty}_{t} e^{\bm{A}(t-\tau)}\bm{B}\bm{u}(\tau) d\tau \\ &\left(\because \lim_{\gamma\rightarrow\infty} e^{-\bm{A}\gamma}\bm{x}(\gamma) = 0\right) \end{align}
時間反転した関数 uˉ(t)u(t)\bar{\bm{u}}(t)\equiv \bm{u}(-t) と 変数 στ[,t]\sigma\equiv-\tau \in [-\infty, -t] を導入すれば,以下の解を得る。
x(t)=teA(t+σ)Bu(σ)dσ=teA(tσ)Buˉ(σ)dσ=tˉeA(tˉσ)Buˉ(σ)dσtˉ=t\begin{align} \bm{x}(t) &= \int^{-t}_{-\infty} e^{\bm{A}(t+\sigma)}\bm{B}\bm{u}(-\sigma) d\sigma\\ &= \int^{-t}_{-\infty} e^{-\bm{A}(-t-\sigma)}\bm{B}\bar{\bm{u}}(\sigma) d\sigma\\ &= \left. \int^{\bar{t}}_{-\infty} e^{-\bm{A}(\bar{t}-\sigma)}\bm{B}\bar{\bm{u}}(\sigma) d\sigma \right|_{\bar{t}=-t} \end{align}
伝達関数表現
次の入出力伝達関数を持つシステムについて考える。
X(s)=G(s)U(s)\begin{align} \bm{X}(s) = \bm{G}(s) \bm{U}(s) \end{align}
ただし,G(s)G(s) の極は全て右半平面に存在するものとする。このシステムの応答 x(t)\bm{x}(t) が有界であると仮定すれば,畳み込みを用いて次のように記述される。
x(t)=x(t0)+t0tg(tτ)u(τ)dτg(t)L1[G(s)](t)\begin{align} \bm{x}(t) &= \bm{x}(t_{0}) + \int^{t}_{t_{0}} \bm{\bm{g}}(t-\tau)\bm{u}(\tau) d\tau \\ \bm{\bm{g}}(t) &\equiv \mathcal{L}^{-1} \left[\bm{G}(s)\right](t) \end{align}
ただし,t0t_{0} は任意の時刻とする。無限遠の時刻を起点として解を整理すれば,以下の得る。
x(t)=limt0(x(t0)+tg(tτ)u(τ)dτ)=tg(tτ)u(τ)dτ=tg(tτ)u(τ)dτ\begin{align} \bm{x}(t) &= \lim_{t_{0}\rightarrow\infty} \left(\bm{x}(t_{0}) + \int^{t}_{\infty} \bm{\bm{g}}(t-\tau)\bm{u}(\tau) d\tau \right)\\ &= \int^{t}_{\infty} \bm{\bm{g}}(t-\tau)\bm{u}(\tau) d\tau = - \int^{\infty}_{t} \bm{g}(t-\tau)\bm{u}(\tau) d\tau \end{align}
時間反転した関数 gˉ(t)g(t), uˉ(t)u(t)\bar{\bm{g}}(t)\equiv \bm{\bm{g}}(-t),\ \bar{\bm{u}}(t)\equiv \bm{u}(-t) と 変数 στ[,t]\sigma\equiv-\tau \in [-\infty, -t] を導入すれば,以下の解を得る。
x(t)=tg(t+σ)u(σ)dσ=tgˉ(tσ)uˉ(σ)dσ=tˉgˉ(tˉσ)uˉ(σ)dσtˉ=t\begin{align} \bm{x}(t) &= \int^{-t}_{-\infty} \bm{\bm{g}}(t + \sigma)\bm{u}(-\sigma) d\sigma\\ &= \int^{-t}_{-\infty} \bar{\bm{g}}(-t - \sigma)\bar{\bm{u}}(\sigma) d\sigma\\ &= \left. \int^{\bar{t}}_{-\infty} \bar{\bm{g}}(\bar{t} - \sigma) \bar{\bm{u}}(\sigma) d\sigma \right|_{\bar{t}=-t} \end{align}

参考文献

[1]L.R. Hunt, G. Meyer, and R. Su, "Noncausal inverses for linear systems," in IEEE Transactions on Automatic Control, vol. 41, no. 4, pp. 608-611, Apr. 1996.
[2]S. Devasia, Degang Chen, and B. Paden, "Nonlinear inversion-based output tracking," in IEEE Transactions on Automatic Control, vol. 41, no. 7, pp. 930-942, Jul. 1996.
[3]木下 浩二,十河 拓也,足立 紀彦, "勾配法による反復学習制御とStable Inversionとの関連性について," in 計測自動制御学会論文集, vol. 36, no. 12, pp. 40-46, Dec. 2000.
[4]T. Sogo, "On the equivalence between stable inversion for nonminimum phase systems and reciprocal transfer functions defined by the two-sided Laplace transform," in Automatica, vol. 46, no. 1, pp. 122-126, Jul. 2010.
[5]十河 拓也, "モデルマッチング非因果的解の計算法とその予見フィードフォワード制御への応用," in 計測自動制御学会論文集, vol. 42, no. 1, pp. 40-46, Jan. 2006.
[6]B. van der Pol and H. Bremmer, "Operational Calculus Based on the Two-Sided Laplace Integral (Third edition)," in Cambridge University Press, 1987.