制御系CAD

制御を意識したモノつくりを目指して

出力FB型SM制御(p=m)

投稿日時: 2021年11月1日投稿者: cacsd

出力数＝入力数の場合…Homework

これまでは、SM制御則の線形制御部を状態フィードバックの形で求めていましたが、以下では出力フィードバックに置き換えることを検討していきます。

[1] 制御対象の状態空間表現として次式を考えます。

$\displaystyle{ \begin{array}{ll} (1.1) & \dot{x}(t)=Ax(t)+Bu(t)+f(t,x,u)\\ (1.2) & y(t)=Cx(t)\\ & (x(t)\in{\rm\bf R}^n, u(t)\in{\rm\bf R}^m, y(t)\in{\rm\bf R}^p, 1\le \boxed{m=p} < n) \end{array} }$

ただし、 $f(t,x,u)\in{\rm\bf R}^n$ はモデル誤差、非線形要素、外乱などの影響を表し、次のマッチング条件を満たすとします。

$\displaystyle{ \begin{array}{ll} (2.1) & f(t,x,u)=B\xi(t,x,u)\\ (2.2) & ||\xi(t,x,u)||<k_1||u||+\alpha(t,y)\quad(k_1<1)\\ & (f(t,x,y)\in{\rm\bf R}^n, \xi(t,x,u)\in{\rm\bf R}^m) \end{array} }$

これに対し、次のスイッチング関数を定義します。

$\displaystyle{(3)\quad \boxed{s(t)=Fy(t)}\quad(s(t)\in{\rm\bf R}^m, F\in{\rm\bf R}^{m\times m}) }$

●状態方程式(1.1)は次のSM標準形をとるように座標変換されているとします(Note C22-1参照)。

$\displaystyle{(4.1)\quad \begin{array}{l} \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot x_2(t) \end{array}\right] }_{\dot{x}(t)} = \underbrace{ \left[\begin{array}{cc} A_{11} & A_{12} \\ A_{21} & A_{22} \\ \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} 0\\ B_2 \end{array}\right] }_{B} \underbrace{(u(t)+\xi(t,x,u))}_{u'(t)}\\ (x_1(t)\in{\rm\bf R}^{n-m}, x_2(t)\in{\rm\bf R}^m, B_2\in{\rm\bf R}^{m\times m}) \end{array} }$

これに対応して、観測方程式(1.2)も次のように表しておきます。

$\displaystyle{(4.2)\quad y(t)= \underbrace{ \left[\begin{array}{cc} C_{1} & C_{2} \\ \end{array}\right] }_{C} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x(t)}\quad (y(t)\in{\rm\bf R}^m, C_2\in{\rm\bf R}^{m\times m}) }$

すなわち、以下では行列 $A,B,C$ は(4.1)と(4.2)のような分割をもつとします。

●スイッチング関数(3)を

$\displaystyle{(5)\quad s(t)=Fy(t)=F \underbrace{ \left[\begin{array}{cc} C_{1} & C_{2} \\ \end{array}\right] }_{C} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x(t)} =\underbrace{ \left[\begin{array}{cc} FC_{1} & FC_{2} \\ \end{array}\right] }_{S=FC} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x(t)} }$

と表すと、等価制御が存在する前提条件 ${\rm rank}\,SB=m$ は

$\displaystyle{(6)\quad {\rm rank}\,\underbrace{FCB}_{SB}={\rm rank}\,FC_2B_2=m }$

となります。以下では ${\rm rank}\,CB=m$ を仮定します。このとき、正方行列 $F,C_2,B_2$ はすべて正則となります。

[2] (4.1)に対して、座標変換

$\displaystyle{(7)\quad \begin{array}{l} \underbrace{ \left[\begin{array}{c} x_1(t)\\ y(t) \end{array}\right] }_{\bar{x}(t)} = \underbrace{ \left[\begin{array}{cc} I_{n-m}& 0 \\ C_1 & C_2 \\ \end{array}\right] }_{T_y} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x(t)}\\ \Leftrightarrow \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x(t)} = \underbrace{ \left[\begin{array}{cc} I_{n-m}& 0 \\ -C_2^{-1}C_1 & C_2^{-1} \end{array}\right] }_{T_y^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ y(t) \end{array}\right] }_{\bar{x}(t)} \end{array} }$

を行うために、まず(7)を(4.1)に代入して（ $M=C_2^{-1}C_1$ ）

$\displaystyle{(8)\quad \begin{array}{l} \underbrace{ \left[\begin{array}{cc} I_{n-m} & 0 \\ -M & C_2^{-1} \\ \end{array}\right] }_{T_y^{-1}} \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot y(t) \end{array}\right] }_{\dot{\bar{x}}(t)}\\ = \underbrace{ \left[\begin{array}{cc} A_{11} & A_{12} \\ A_{21} & A_{22} \\ \end{array}\right] }_{A} \underbrace{ \left[\begin{array}{cc} I_{n-m} & 0 \\ -M & C_2^{-1} \\ \end{array}\right] }_{T_y^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ y(t) \end{array}\right] }_{\bar{x}(t)} + \underbrace{ \left[\begin{array}{c} 0\\ B_2 \end{array}\right] }_{B} u'(t) \end{array} }$

を得て、左から $T_y$ をかけて

$\displaystyle{(9)\quad \begin{array}{l} \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot y(t) \end{array}\right] }_{\dot{\bar{x}}(t)} = \underbrace{ \left[\begin{array}{cc} I_{n-m} & 0 \\ C_1 & C_2 \\ \end{array}\right] }_{T_y} \underbrace{ \left[\begin{array}{cc} A_{11}-A_{12}M & A_{12}C_2^{-1} \\ A_{21}-A_{22}M & A_{22}C_2^{-1} \\ \end{array}\right] }_{AT_y^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ y(t) \end{array}\right] }_{\bar{x}(t)}\\ + \underbrace{ \left[\begin{array}{cc} I_{n-m} & 0 \\ C_1 & C_2 \\ \end{array}\right] }_{T_y} \underbrace{ \left[\begin{array}{c} 0\\ B_2 \end{array}\right] }_{B} u'(t) \end{array} }$

すなわち、次式を得ます。

$\displaystyle{(10)\quad \boxed{ %\begin{array}{l} \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot y(t) \end{array}\right] }_{\dot{\bar{x}}(t)} = \underbrace{ \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ \bar{A}_{21} & \bar{A}_{22}\\ \end{array}\right] }_{\bar{A}=T_y A T_y^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ y(t) \end{array}\right] }_{\bar{x}(t)} + \underbrace{ \left[\begin{array}{c} 0\\ \bar{B}_2 \end{array}\right] }_{\bar{B}=T_yB} u'(t) %\end{array} }}$

ただし

$\displaystyle{(10')\quad \left\{\begin{array}{l} \bar{A}_{11}=A_{11}-A_{12}M\quad(M=C_2^{-1}C_1)\\ \bar{A}_{12}=A_{12}C_2^{-1}\\ \bar{A}_{21}=C_2(M\bar{A}_{11}+A_{21}-A_{22}M)\\ \bar{A}_{22}=C_2(M{A}_{12}+A_{22})C_2^{-1}\\ \bar{B}_{2}=C_2B_2 \end{array}\right. }$

また、(4.2)と(3)は、次式のように表せます。

$\displaystyle{(11)\quad \boxed{y(t)= \underbrace{ \left[\begin{array}{cc} 0 & I_{m} \\ \end{array}\right] }_{\bar{C}=CT_y^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ y(t) \end{array}\right] }_{\bar{x}(t)}} }$

$\displaystyle{(12)\quad \boxed{s(t)=Fy(t)= \underbrace{ \left[\begin{array}{cc} 0 & F \\ \end{array}\right] }_{F\bar{C}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ y(t) \end{array}\right] }_{\bar{x}(t)}} }$

ここで、 $\bar{A}_{11}=A_{11}-A_{12}M$ の固有値が $(A,B,C)$ の不変零点に等しく、スライディングモードの振舞いを決めます。しかし、 $M=C_2^{-1}C_1$ ですから、 $M$ には選択の自由度がないことに注意します。

[3] 以下では、(10),(11),(12)を考えますが、煩雑さを避けるために $\bar{A}$ 、 $\bar{B}$ 、 $\bar{C}$ 、 $\bar{x}$ のバーをとって記述します。

●さて、出力FB型のＳＭ制御則として次式を考えます。

$\displaystyle{(13)\quad \boxed{ u(t)=-\gamma Fy(t)\underbrace{-\rho(t,y)\frac{Fy(t)}{||Fy(t)||}}_{\nu_c(t)}\quad(\gamma>\gamma_0) }}$

ここで、 $\gamma$ と $\rho(t,y)$ は設計パラメータ、 $\gamma_0$ は後で示す定数です。

これによる閉ループ系は、(13)を(4.1)に代入して、次式となります。

$\displaystyle{(14)\quad \dot{x}(t)=(A-\gamma BFC)x(t)+B(-\rho(t,y)\frac{Fy(t)}{||Fy(t)||}+\xi(t,x,y)) }$

●このとき、SM制御則(13)を、２次安定性

$\displaystyle{(15)\quad \boxed{ V(x)=x(t)^TPx(t) \Rightarrow \dot{V}(x)\le -x(t)^TL(\gamma)x(t)} }$

が成り立つように決定します（ $P>0$ , $L(\gamma)>0$ ）。ただし、次式を仮定します。

$\displaystyle{(16)\quad B^TP=FC }$
$\displaystyle{(17)\quad \rho(t,y)=\frac{k_1\gamma||Fy(t)||+\alpha(t,y)+\gamma_2}{1-k_1} }$

また、(17)と(13)から次式が得られることに注意します。

$\displaystyle{ \begin{array}{cl} (18.1) &\rho(t,y)=k_1(\gamma||Fy(t)||+\rho(t,y))+\alpha(t,y)+\gamma_2\\ \Downarrow & using\ ||u(t)||\le \gamma||Fy(t)||+\rho(t,y)\\ (18.2) & \ge k_1||u(t)||+\alpha(t,y)+\gamma_2\\ \end{array} }$

●これらの準備の下で、(15)が次のように示されます。

$\displaystyle{ \begin{array}{cl} (19.1) & \dot{V}(x)=2x(t)^TP\dot{x}(t)\\ (19.2) & =2x^T(t)P((A-\gamma BFC)x(t)+B(-\rho(t,y)\frac{Fy(t)}{||Fy(t)||}+B\xi(t,x,u))\\ (19.3) & =x^T(t)\underbrace{(P(A-\gamma BFC)+(A-\gamma BFC)^TP)}_{L(\gamma)}x(t)\\ & -2x^T(t)\underbrace{PB}_{C^TF^T}\rho(t,y)\frac{Fy(t)}{||Fy(t)||}+2x^T(t)\underbrace{PB}_{C^TF^T}\xi(t,x,u))\\ (19.4) & =-x^T(t)L(\gamma)x(t)-2\rho(t,y)||Fy(t)||+2(Fy(t))^T\xi(t,x,u)\\ \Downarrow & by\ (Fy(t))^T\xi(t,x,u)\le||Fy(t)||\cdot||\xi(t,x,u)||\\ (19.5) & \le -x^T(t)L(\gamma)x(t)-2\rho(t,y)||Fy(t)||+2||Fy(t)||\cdot||\xi(t,x,u)||\\ (19.6) & = -x^T(t)L(\gamma)x(t)-2||Fy(t)||(\rho(t,y)-||\xi(t,x,u)||)\\ \Downarrow & by\ (2.2):||\xi(t,x,u)||<k_1||u||+\alpha(t,y)\quad(k_1<1)\\ (19.7) & \le -x^T(t)L(\gamma)x(t)-2||Fy(t)||(\rho(t,y)-k_1||u||-\alpha(t,y))\\ \Downarrow & by\ (18):\rho(t,y)\ge k_1||u(t)||+\alpha(t,y)+\gamma_2\\ (19.8) & \le -x^T(t)L(\gamma)x(t)-2\gamma_2||Fy(t)||\\ (19.9) & \le -x^T(t)L(\gamma)x(t) \end{array} }$

●ここで、次式を示すことができします。

$\displaystyle{(20)\quad \boxed{L(\gamma)=\underbrace{P(A-\gamma BFC)+(A-\gamma BFC)^TP}_{PA+A^TP-2\gamma (FC)^TFC}<0 \Leftrightarrow \gamma>\gamma_0} }$

実際、左辺は

$\displaystyle{ \begin{array}{cl} (21.1) & L(\gamma)=PA+A^TP-2\gamma (FC)^TFC\\ (21.2) & =\underbrace{\left[ \begin{array}{cc} P_1 & 0 \\ 0 & P_2 \\ \end{array}\right]}_{P} \underbrace{\left[\begin{array}{cc} A_{11} & A_{12} \\ A_{21} & A_{22} \\ \end{array}\right]}_{A} + \underbrace{\left[\begin{array}{cc} A_{11} & A_{12} \\ A_{21} & A_{22} \\ \end{array}\right]^T}_{A^T} \underbrace{\left[\begin{array}{cc} P_1 & 0 \\ 0 & P_2 \\ \end{array}\right]}_{P}\\ &-2\gamma \underbrace{\left[\begin{array}{cc} 0 & F \\ \end{array}\right]^T}_{(FC)^T} \underbrace{\left[\begin{array}{cc} 0 & F \\ \end{array}\right]}_{FC}\\ (21.3) & =\left[ \begin{array}{cc} P_1A_{11}+A_{11}^TP_1 & P_1A_{12}+A_{21}^TP_2 \\ (P_1A_{12}+A_{21}^TP_2)^T & P_2A_{22}+A_{22}^TP_2 \\ \end{array}\right] -2\gamma \left[\begin{array}{cc} 0 & 0 \\ 0 & F^TF \\ \end{array}\right]\\ (21.4) & =\left[\begin{array}{cc} P_1A_{11}+A_{11}^TP_1 & P_1A_{12}+A_{21}^TP_2 \\ (P_1A_{12}+A_{21}^TP_2)^T & P_2A_{22}+A_{22}^TP_2-2\gamma F^TF \\ \end{array}\right]<0 \end{array} }$

となり、これは公式

$\displaystyle{(22)\quad \begin{array}{l} \left[\begin{array}{cc} X & Z \\ Z^T & Y \end{array}\right]<0\\ \ \Leftrightarrow\ X-ZY^{-1}Z^T<0,\ Y<0 \\ \ \Leftrightarrow\ X<0,\ Y-Z^TX^{-1}Z<0 \end{array} }$

を用いて、以下と等価になります。

$\displaystyle{ \begin{array}{cl} (23.1) & \left\{\begin{array}{l} P_1A_{11}+A_{11}^TP_1<0\\ P_2A_{22}+A_{22}^TP_2 -2\gamma F^TF-(P_1A_{12}+A_{21}^TP_2)^T\\ \times(P_1A_{11}+A_{11}^TP_1)^{-1}(P_1A_{12}+A_{21}^TP_2)<0 \end{array}\right.\\ (23.2) & \Leftrightarrow\ \left\{\begin{array}{l} P_1A_{11}+A_{11}^TP_1<0\\ P_2A_{22}+A_{22}^TP_2 -(P_1A_{12}+A_{21}^TP_2)^T\\ \times(P_1A_{11}+A_{11}^TP_1)^{-1}(P_1A_{12}+A_{21}^TP_2)<2\gamma F^TF \end{array}\right.\\ (23.3) & \Leftrightarrow\ \left\{\begin{array}{l} P_1A_{11}+A_{11}^TP_1<0\\ F^{-T}(P_2A_{22}+A_{22}^TP_2 -(P_1A_{12}+A_{21}^TP_2)^T\\ \times(P_1A_{11}+A_{11}^TP_1)^{-1}(P_1A_{12}+A_{21}^TP_2))F^{-1}<2\gamma I \end{array}\right.\\ (23.4) & \Leftrightarrow\ \left\{\begin{array}{l} P_1A_{11}+A_{11}^TP_1<0\quad(*)\\ \underbrace{\begin{array}{l} \frac{1}{2}\bar{\sigma}(F^{-T}(P_2A_{22}+A_{22}^TP_2 -(P_1A_{12}+A_{21}^TP_2)^T\\ \times(P_1A_{11}+A_{11}^TP_1)^{-1}(P_1A_{12}+A_{21}^TP_2))F^{-1}) \end{array} }_{\gamma_0}<\gamma \end{array}\right. \end{array} }$

ここで、(*)は仮定より満足されるので、 $\gamma>\gamma_0$ が $L(\gamma)>0$ となるための条件となることが分かります。

SMCまとめ

投稿日時: 2021年10月29日投稿者: cacsd

２次安定性

[1]　非線形系

$\displaystyle{(1)\quad \dot{x}(t)=f(x(t))}$

に対してリャプノフ関数

$\displaystyle{(5)\quad V(x)=x^T(t)Px(t)\quad(P>0)}$

を考えます。これを微分して(1)を用いると（(1)の解に沿って微分すると言います）

$\displaystyle{(6)\quad %\begin{array}{l} \frac{d}{dt}V(x)=\dot{x}^T(t)Px(t)+x^T(t)P\dot{x}(t)=2x^T(t)P\dot{x}(t)=2x^T(t)Pf(x) %\end{array} }$

となります。このとき２次安定とは

$\displaystyle{(7)\quad { \boxed{\frac{d}{dt}V(x)=2x^T(t)Pf(x)\le -x^T(t)Qx(t)\quad(Q>0)} }}$

が成り立つことを言います。これから次式を示すことができ、漸近安定性が成り立ちます。

$\displaystyle{(8)\quad \boxed{||x(t)||\le \beta||x(0)||e^{-\frac{1}{2}\alpha t}} }$

ただし

$\displaystyle{(9)\quad {\alpha=\sigma_n(P^{-1}Q),\ \beta=\sqrt{\frac{\sigma_1(P)}{\sigma_n(P)} }}$

[2]　制御対象が平衡状態にあることは線形状態方程式 $\dot{x}(t)=Ax(t)+Bu(t)$ において、 $x=0, u=0$ を意味します。そこで平衡状態が乱されて $x\ne0$ となる時刻を $t=0$ にとると、線形状態方程式は次式となります。

$\displaystyle{(1)\quad \dot{x}(t)=Ax(t)\qquad(x(0)\ne0) }$

これに出力方程式

$\displaystyle{(2)\quad y(t)=Cx(t) }$

を考慮する場合も含めると、漸近安定性の判定法は次のようにまとめられます。

【漸近安定性の定義とその等価な条件】
定義DA:　 $\forall x(0)\ne 0: x(t)=\exp(At)x(0)\rightarrow 0\quad(t\rightarrow\infty)$
条件A0:　 $\exp(At)\rightarrow 0\quad(t\rightarrow\infty)$
条件A1:　 ${\rm Re}(\lambda_i(A))<0\quad(i=1,\cdots,n)$
条件A2:　 $\exists P>0: PA+A^TP+I_n=0$
条件A3:　 $\exists P>0: PA+A^TP+C^TC=0$ 　ただし、 $(A,C)$ は可観測対
条件A4:　 $\exists P>0: PA+A^TP<0$

まず(1)が漸近安定とすると、条件A2のリャプノフ方程式から次式を得ます。

$\displaystyle{(3)\quad x^T(t)(PA+A^TP)x(t)=-x(t)^Tx(t) }$

次に(1)に対してリャプノフ関数

$\displaystyle{(4)\quad V(x)=x^T(t)Px(t) }$

を考えるとき、(3)は２次安定性の条件

$\displaystyle{(5)\quad \begin{array}{l} \frac{d}{dt}V(x)=2x^T(t)P\dot{x}(t)=2x^T(t)PAx(t)=x^T(t)(PA+A^TP)x(t)\\ =-x^T(t)x(t) \end{array} }$

を意味します。これから漸近安定な線形系は２次安定であることが分かります。

●条件A2の両辺から $Q^{1/2}\ (Q>0)$ をかけて

$\displaystyle{(6)\quad \underbrace{Q^{1/2}PQ^{1/2}}_{P'}\underbrace{Q^{-1/2}AQ^{1/2}}_{A'}+\underbrace{Q^{1/2}A^TQ^{-1/2}}_{A'^T}\underbrace{PQ^{1/2}Q^{1/2}}_{P'}+\underbrace{Q^{1/2}Q^{1/2}}_{Q}=0 }$

を得ます。 $(A',Q^{1/2})$ は可観測対なので条件A3を適用して（ $'$ をはずして）

条件A2′:　 $\boxed{\exists P>0: PA+A^TP+Q=0}$ 　ただし、 $Q>0$

も漸近安定性の条件となることに注意します。

●たとえば、安定行列 $\displaystyle{ A= \left[\begin{array}{cc} 0 & 1 \\ -4 & -1 \end{array}\right] }$

が与えられるとき、リャプノフ方程式 $PA+A^TP+I_n=0$ を満足する $P>0$ は

$\displaystyle{ P=\left[\begin{array}{cc} 2.6250 & 0.1250\\ 0.1250 & 0.6250 \end{array}\right] }$

のように求まります。これはMATLABで次のコマンドを用いて計算できます。

MATLAB

P=lyap(A',eye(2))

ここで、P=lyap(A,Q) はリャプノフ方程式 $PA^T+AP+Q=0$ の解を求める関数であることに注意します。

問題設定

[3]　制御対象の状態方程式として次式を考えます。

$\displaystyle{(1)\quad \begin{array}{l} \dot{x}(t)=Ax(t)+Bu(t)+f(t,x,u)\\ (x(t)\in{\rm\bf R}^n, u(t)\in{\rm\bf R}^m, f(t,x,u)\in{\rm\bf R}^n, m< n, {\rm rank}B=m) \end{array} }$

ここで、 $f(t,x,u)$ はモデル誤差、非線形要素、外乱などの影響を表し、有界かつ未知とします。これに対し、次のスイッチング関数を定義します。

$\displaystyle{(2)\quad \boxed{s(t)=Sx(t)}\quad(s(t)\in{\rm\bf R}^m, {\rm rank}S=m) }$

いま(1)の解が、ある時刻 $t_s$ に対して

$\displaystyle{(3)\quad s(t)=0\quad(t\ge t_s) }$

を満足するとき、スライディングモードが発生していると言います。

●以上の準備の下で、スライディングモード制御問題は次のように記述されます。

問題１: 動的システムの振る舞いが閉じ込められる超平面 ${\cal S}=\{x:Sx=0\}$ を定義するスイッチング関数を決定せよ。
問題２: 有限時刻 $t_s$ において状態を超平面内に拘束し、引き続き留まらせる( $\forall t>t_s: s(t)=0$ )ことのできるスライディングモード制御則（SM制御則）を設計せよ。

●いま、SM制御則を、２次安定性の条件
を、２次安定性

$\displaystyle{(4)\quad %\boxed{ \begin{array}{lll} V(\bar{x})= \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]^T }_{\bar{x}^T(t)} \underbrace{ \left[\begin{array}{cc} P_1 & 0\\ 0 & P_2 \end{array}\right] }_{P} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{\bar{x}(t)}\\ \Rightarrow \dot{V}(\bar{x})\le - \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]^T }_{\bar{x}^T(t)} \underbrace{ \left[\begin{array}{cc} Q_1 & 0\\ 0 & I \end{array}\right] }_{Q} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{\bar{x}(t)} \end{array}}} }$

すなわち

$\displaystyle{(4')\quad \begin{array}{lll} V(\bar{x})= \underbrace{x_1^T(t)P_1x_1(t)}_{V(x_1)}+\underbrace{s^T(t)P_2s(t)}_{V(s)}\\ \Rightarrow \dot{V}(\bar{x})\le -x_1^T(t)Q_1x_1(t)-s^T(t)s(t) \end{array}} }$

が成り立つように決定します（ $P_1>0$ , $P_2>0$ , $Q_1>0$ ）。

以下で示すように、このSM制御則は次式で表されます。

$\displaystyle{(5)\quad { \boxed{\begin{array}{l} u(t)=\underbrace{u_\ell(t)}_{linear\ control}+\underbrace{u_n(t)}_{switching\ component}\\ u_\ell(t)=-\underbrace{(SB)^{-1}(SA-\Phi S)}_{L=L_{eq}+L_{phi}}x(t) \quad(\Phi:stable\ matrix)\\ u_n(t)=-\underbrace{(SB)^{-1}\rho(t,x)}_{L_n}\frac{P_2s(t)}{||P_2s(t)||}\quad(P_2\Phi+\Phi^TP_2=-I) \end{array}}} }$

ここで、 $SB$ は正則行列であることを仮定しています。また設計パラメータは、スイッチング関数を決める行列 $S$ 、安定行列 $\Phi$ 、スカラー $\rho>0$ です。

●２重積分器に対して、 $S$ , $\Phi$ , $\rho$ を適当に選んだときのSM制御則の計算は次のように行われます。

MATLAB

%cDINT_smc.m
%-----
clear all, close all
A=[0 1;
   0 0];
B=[0;1];
C=eye(2); 
D=zeros(2,1);
M=1;
S=[M 1];
Phi=-1;
P2=lyap(Phi,1);
P2S=P2*S;
Leq=inv(S*B)*S*A;
LPhi=-inv(S*B)*Phi*S;
L=Leq+LPhi;
rho=0.5;
Ln=inv(S*B)*rho;
x0=B;
sim('DINT_smc')
%-----
%eof

これを実行して、図1aのシミュレーション結果を得ます。

図1a　２重積分器のSM制御

図1b　２重積分器の線形制御（SM制御のスイッチング項を外した場合）

図1bに、SM制御のスイッチング項を外した場合のシミュレーション結果を示します。これらからスイッチング項が速応性の改善に役立っていることが分かります。

●いま(1)において $f(t,x,u)=0$ 、すなわち

$\displaystyle{(6)\quad \dot{x}(t)=Ax(t)+Bu(t) }$

とし、 $SB$ の正則性を仮定します。 $t\ge t_s$ において、スライディングモード時の制御則は等価制御と呼ばれ、次式のように求められます。

$\displaystyle{(7)\quad \begin{array}{l} s(t)=0\\ \Rightarrow\dot{s}(t)=0\\ \Rightarrow S\dot{x}(t)=SAx(t)+SBu(t)=0\\ \Rightarrow \boxed{u_{eq}(t)=-\underbrace{(SB)^{-1}SA}_{L_{eq}}x(t)} \end{array} }$

これによる閉ループ系は次式となります。

$\displaystyle{(8)\quad \boxed{\dot{x}(t)=\underbrace{(I_n-B(SB)^{-1}S)A}_{A_{eq}}x(t)}\quad(t\ge t_s) }$

●いま(1)において $f(t,x,u)\ne 0$ の場合、適合サイズをもつ行列 $R$ を用いて

$\displaystyle{(9)\quad　{ \boxed{f(t,x,u)=BR\xi(t,x)}} }$

のように表されると仮定します。このとき、等価制御は

$\displaystyle{(10)\quad \begin{array}{l} s(t)=0\\ \Rightarrow\dot{s}(t)=0\\ \Rightarrow S\dot{x}(t)=SAx(t)+SBu(t)+SBR\xi(t,x)=0\\ \Rightarrow u(t)=-(SB)^{-1}(SAx(t)+SBR\xi(t,x)) \end{array} }$

となります。これによる閉ループ系は次式となります。

$\displaystyle{(11)\quad \begin{array}{l} \dot{x}(t)=Ax(t)-B(SB)^{-1}(SAx(t)+SBR\xi(t,x))+BR\xi(t,x)\\ =\underbrace{(I_n-B(SB)^{-1}S)A}_{A_{eq}}x(t)+ \underbrace{(I_n-B(SB)^{-1}S)B}_0R\xi(t,x)\\ =\underbrace{(I_n-B(SB)^{-1}S)A}_{A_{eq}}x(t) \end{array} }$

すなわちマッチング条件とよばれる(9)が成り立つとき、スライディングモード時は $f(t,x,u)$ の影響を受けないことが分かります。

●剛体振子の運動方程式は次式で与えられました。

$\displaystyle{(12)\quad \ddot{\theta}(t)=-\underbrace{\frac{3g}{4\ell}}_{a}\sin\theta(t)+\underbrace{\frac{1}{J}\tau(t)}_{u(t)} }$

$x_1=\theta$ , $x_2=\dot{\theta}$ とおいて、次の非線形状態方程式を得ます。

$\displaystyle{(13)\quad \left[\begin{array}{c} \dot{x}_1(t)\\ \dot{x}_2(t) \end{array}\right]= \left[\begin{array}{cc} x_2(t)\\ -a\sin x_1(t)+u(t) \end{array}\right] }$

これは、次式のように、２重積分器に非線形項が加えられたものとみなすことができます。

$\displaystyle{(14)\quad \left[\begin{array}{c} \dot{x}_1(t)\\ \dot{x}_2(t) \end{array}\right]= \left[\begin{array}{cc} 0 & 1\\ 0 & 0 \end{array}\right] \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right]+ \left[\begin{array}{c} 0\\ 1 \end{array}\right] u(t)+ \left[\begin{array}{c} 0\\ -a\sin x_1(t) \end{array}\right] }$

そこで、図1のSM制御則をそのまま用いたときの様子を調べてみます。図2は $a=15$ とした場合（ $\ell=0.5$ ）のシミュレーション結果です。

図2a　非線形項を考慮した場合のSM制御

図2b　非線形項を考慮した場合の線形制御（SM制御のスイッチング項を外した場合）

図1のSM制御則は非線形項を考慮していないのですが、マッチング条件が満たされているので、非線形項の影響を抑制していることが分かります。これは図2bのスイッチング項がない場合と比べると明らかです。

SM制御則

[4]　制御対象の状態方程式として次式を考えます。

$\displaystyle{(1)\quad \begin{array}{l} \dot{x}(t)=Ax(t)+Bu(t)+f(t,x,u)\\ (x(t)\in{\rm\bf R}^n, u(t)\in{\rm\bf R}^m, f(t,x,u)\in{\rm\bf R}^n) \end{array} }$

以下では、状態方程式とスイッチング関数は次のSMC標準形(regular form)をとるように座標変換されていると仮定します（ ${\rm rank}B={\rm rank}S=m$ ）。

$\displaystyle{(2)\quad \boxed{\underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot x_2(t) \end{array}\right] }_{\dot{x}(t)} = \underbrace{ \left[\begin{array}{cc} A_{11} & A_{12} \\ A_{21} & A_{22} \\ \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} 0\\ B_2 \end{array}\right] }_{B} u(t) + \underbrace{ \left[\begin{array}{l} f_u(t,x)\\ f_m(t,x,u) \end{array}\right] }_{f(t,x,u)}} }$

$\displaystyle{(3)\quad \boxed{s(t)= \underbrace{ \left[\begin{array}{cc} S_1 & S_2 \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x(t)} = \underbrace{S_2 \left[\begin{array}{cc} M & I \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x(t)} \ (M=S_2^{-1}S_1)} }$

これに対して、座標変換

$\displaystyle{(4)\quad \begin{array}{l} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{\bar{x}(t)} = \underbrace{ \left[\begin{array}{cc} I & 0 \\ S_1 & S_2 \\ \end{array}\right] }_{T_s} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x(t)}\\ \Leftrightarrow \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x(t)} = \underbrace{ \left[\begin{array}{cc} I & 0 \\ -S_2^{-1}S_1 & S_2^{-1} \\ \end{array}\right] }_{T_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{\bar{x}(t)} \end{array} }$

を行って、次式を得ます。

$\displaystyle{(5)\quad \boxed{ \begin{array}{l} \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot s(t) \end{array}\right] }_{\dot{\bar{x}}(t)} = \underbrace{ \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ \bar{A}_{21} & \bar{A}_{22}\\ \end{array}\right] }_{\bar{A}=T_s A T_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{\bar{x}(t)} + \underbrace{ \left[\begin{array}{c} 0\\ \bar{B}_2 \end{array}\right] }_{\bar{B}=T_sB} u(t)\\ + \underbrace{ \left[\begin{array}{c} f_u(t,x)\\ S_1f_u(t,x)+S_2f_m(t,x,u) \end{array}\right] }_{T_sf(t,x,u)} \end{array}} }$

ただし

$\displaystyle{(6)\quad \left\{\begin{array}{l} \bar{A}_{11}=A_{11}-A_{12}M\\ \bar{A}_{12}=A_{12}S_2^{-1}\\ \bar{A}_{21}=S_2(M\bar{A}_{11}+A_{21}-A_{22}M)\\ \bar{A}_{22}=S_2(M{A}_{12}+A_{22})S_2^{-1}\\ \bar{B}_{2}=S_2B_2 \end{array}\right. }$

以下では、 $\bar{A}_{11}$ が安定行列となるようにスイッチング関数が選ばれていると仮定します。

●SMC標準形に対する等価制御は、 $f_u=0, f_m=0$ として

$\displaystyle{ \begin{array}{cl} & s(t)=0\ \Rightarrow\ \dot{s}(t)=0\\ \Downarrow & by\ (5)\\ &\bar{A}_{21}x_1(t)+\bar{A}_{22}s(t)+\bar{B}_{2}u(t)=0\\ \Downarrow & \\ (7) & u_{eq}(t)=-\underbrace{\bar{B}_{2}^{-1}}_{([0\ I_m]T_sB)^{-1}} \underbrace{\left[\begin{array}{cc} \bar{A}_{21} & \bar{A}_{22} \\ \end{array}\right]}_{[0\ I_m]T_sAT_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{\bar{x}(t)}\\ \Downarrow & by\ S=[0\ I_m]T_s, x=T_s^{-1}\bar{x}\\ (8) & \boxed{u_{eq}(t)=-\underbrace{(SB)^{-1}SA}_{L_\ell}x(t)}} \end{array} }$

のように得られます。線形制御則 $u_\ell$ は、この等価制御をベースして

$\displaystyle{ \begin{array}{cl} (9) & u_{eq}(t)=-\underbrace{\bar{B}_{2}^{-1}}_{([0\ I_m]T_sB)^{-1}} (\underbrace{\left[\begin{array}{cc} \bar{A}_{21} & \bar{A}_{22} \\ \end{array}\right]}_{[0\ I_m]T_sAT_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{\bar{x}(t)}-\Phi s(t))\\ \Downarrow & by\ S=[0\ I_m]T_s, x=T_s^{-1}\bar{x}\\ (10) & \boxed{u_\ell(t)=-\underbrace{(SB)^{-1}(SA-\Phi S)}_{L=L_{eq}+L_{phi}}x(t)}} \end{array} }$

のように構成します。ここで、 $\Phi$ は適当に設定された安定行列です。

●このとき閉ループ系は次式で与えられます。

$\displaystyle{(11)\quad \begin{array}{l} \left[\begin{array}{c} \dot{x}_1(t)\\ \dot{s}(t) \end{array}\right] = \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ 0 & \Phi \\ \end{array}\right] \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]\\ + \left[\begin{array}{c} f_u(t,x)\\ \bar{B}_{2}u_n(t)+S_2(Mf_u(t,x)+f_m(t,x,u)) \end{array}\right] \end{array} }$

これからスライディングモード $s=0$ を達成するために $\Phi$ の影響を与えることが分かります。簡単な例で確かめると $\Phi$ はダンピングを与え、速くスライディングモードを達成するためのパラメータであることが分かります。

さて、 $\Phi$ は安定行列なので

$\displaystyle{(12)\quad \begin{array}{l} \boxed{P_2\Phi+\Phi^TP_2=-I} \end{array} }$

を満たす $P_2>0$ を選ぶことができます。これを用いて、スイッチング制御則 $u_n$ を

$\displaystyle{(13)\quad u_n(t)=-\underbrace{\bar{B}_{2}^{-1}}_{(SB)^{-1}}\rho(t,x)\frac{P_2s(t)}{||P_2s(t)||} }$

すなわち

$\displaystyle{(14)\quad \boxed{u_n(t)=-\underbrace{(SB)^{-1}\rho(t,x)}_{L_n}\frac{P_2s(t)}{||P_2s(t)||}} }$

と定めます。線形制御則(10)とスイッチング制御則(14)を用いて、２次安定性を示すことができます。

スイッチング関数（１）

[5]　スライディングモードの場合、 $s(t)=0$ なので

$\displaystyle{(1)\quad s(t)=0 \Rightarrow Sx(t)=S_1x_1(t)+S_2x_2(t)=0\Rightarrow x_2(t)=-\underbrace{S_2^{-1}S_1}_Mx_1(t)}}$

となって、 $x_2$ の振る舞いは $x_1$ によって決まります。一方、 $x_1$ の振る舞いは

$\displaystyle{(2)\quad \dot{x}_1(t)=\underbrace{(A_{11}-A_{12}M)}_{\bar{A}_{11}}x_1(t)}$

によって決まります。これは

$\displaystyle{(3)\quad \dot{x}_1(t)=A_{11}x_1(t)+A_{12}x_2(t)}$

に対する安定化状態フィードバック

$\displaystyle{(4)\quad x_2(t)=-Mx_1(t)}$

による閉ループ系とみなすことができます。

したがって、 $\bar{A}_{11}$ が安定行列となるようにスイッチング関数が選ぶためには、(3)に対する状態フィードバックによる安定化問題を解けばよいことが分かります。

●そこで、次の２次形式評価関数を最小にするように $M$ を決定する方法が提案されています。

$\displaystyle{(5)\quad \boxed{J=\int_0^\infty \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right]^T }_{x^T(t)} \underbrace{ \left[\begin{array}{cc} Q_{11} & Q_{12} \\ Q_{21} & Q_{22} \end{array}\right] }_{Q} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x(t)} dt} }$

これは

$\displaystyle{(6)\quad J=\int_0^\infty(x_1^T(t)\underbrace{(Q_{11}-Q_{12}Q_{22}^{-1}Q_{21})}_{\hat{Q}}x_1(t)}$
$\displaystyle{+\underbrace{(Q_{12}Q_{22}^{-1}x_1(t)+x_2(t))^T}_{\hat{u}^T}\underbrace{Q_{22}}_{\hat{R}}\underbrace{(Q_{12}Q_{22}^{-1}x_1(t)+x_2(t))}_{\hat{u}})dt}$

と書けます。ここで、(3)と(4)をそれぞれ

$\displaystyle{(7)\quad \dot{x}_1(t)=\underbrace{(A_{11}-A_{12}Q_{22}^{-1}Q_{21})}_{\hat{A}}x_1(t)+\underbrace{A_{12}}_{\hat{B}}\underbrace{(Q_{22}^{-1}Q_{21}x_1(t) +x_2(t))}_{\hat{u}(t)}}$

$\displaystyle{(8)\quad \underbrace{Q_{22}^{-1}Q_{21}x_1(t) +x_2(t)}_{\hat{u}(t)}=-\underbrace{(M-Q_{22}^{-1}Q_{21})}_{\hat{F}}x_1(t)}$

と書き直すと、次のような最適化制御問題の定式化ができていることがわかります。すなわち

$\displaystyle{(9)\quad \boxed{\begin{array}{l} \dot{x}_1(t)=\hat{A}x_1(t)+\hat{B}\hat{u}(t)\\ \hat{A}=A_{11}-A_{12}Q_{22}^{-1}Q_{21}\\ \hat{B}=A_{12} \end{array}} }$

を２次系式評価関数

$\displaystyle{(10)\quad \boxed{\begin{array}{l} J=\int_0^\infty(x_1^T(t)\hat{Q}x_1(t)+\hat{u}^T(t)\hat{R}\hat{u}(t))dt\\ \hat{Q}=Q_{11}-Q_{12}Q_{22}^{-1}Q_{21}\\ \hat{R}=Q_{22} \end{array}} }$

を最小にするように

$\displaystyle{(11)\quad \boxed{\begin{array}{l} \hat{u}(t)=-\hat{F}x_1(t)\\ \hat{F}=M-Q_{22}^{-1}Q_{21} \end{array}} }$

を求める問題です。この解は、リッカチ方程式

$\displaystyle{(12)\quad \boxed{\Pi\hat{A}+\hat{A}^T\Pi-\Pi\hat{B}\hat{R}^{-1}\hat{B}^T\Pi+\hat{Q}=0}}$

を解いて

$\displaystyle{(13)\quad \boxed{\hat{F}=\hat{R}^{-1}\hat{B}^T\Pi}}$

を求め、次式から $M$ を決めればよいことがわかります。

$\displaystyle{(14)\quad \boxed{M=\hat{F}+Q_{22}^{-1}Q_{21}}}$

●２重積分器に対して、関数swflqrによるスイッチング関数Sの決定は次のように行われます。

MATLAB

%cDINT_smc1.m
%-----
clear all, close all
A=[0 1;
   0 0];
B=[0;1];
C=eye(2); 
D=zeros(2,1);
Mx=0.5; Tx=0.5;
Q=diag([1/Mx,Tx/Mx].^2);
S=swflqr(A,B,Q)
Phi=-1;
P2=lyap(Phi,1);
P2S=P2*S;
Leq=inv(S*B)*S*A;
LPhi=-inv(S*B)*Phi*S;
L=Leq+LPhi;
rho=0.5;
Ln=inv(S*B)*rho;
x0=B;
sim('DINT_smc')
%-----
%eof

これを実行して図3のシミュレーション結果を得ます。

図3　関数swflqrによるスイッチング関数Sの決定

ここで、関数swflqrは、Edwards and Spurgeon によるものです。

スイッチング関数（２）

[6]　等価制御

$\displaystyle{(1)\quad u_{eq}(t)=-(SB)^{-1}SAx(t) }$

による閉ループ系は次式で表されました。

$\displaystyle{(2)\quad \dot{x}(t)=\underbrace{(I_n-B(SB)^{-1}S)A}_{A_{eq}}x(t)\quad(t\ge t_s) }$

ここで

$\displaystyle{(3)\quad SA_{eq}=S(I_n-B(SB)^{-1}S)A=0 }$

また、 $A$ 、 $B$ 、 $S$ は、SMC標準形で与えられるとすると

$\displaystyle{(4)\quad \begin{array}{l} A_{eq}= \underbrace{ \left[\begin{array}{cc} I & 0 \\ -M & 0 \\ \end{array}\right] }_{I_n-B(SB)^{-1}S} \left[\begin{array}{cc} A_{11} & A_{12} \\ A_{21} & A_{22} \\ \end{array}\right]= \left[\begin{array}{cc} A_{11} & A_{12} \\ -MA_{11} & -MA_{12} \\ \end{array}\right]\\ = \left[\begin{array}{cc} I & 0 \\ -M & I \\ \end{array}\right] \left[\begin{array}{cc} \bar{A}_{11} & A_{12} \\ 0 & 0 \\ \end{array}\right] \left[\begin{array}{cc} I & 0 \\ -M & I \\ \end{array}\right]^{-1} \end{array} }$

となって、 $\bar{A}_{11}$ が安定行列なので、 $A_{eq}$ は $m$ 個の零固有値を持ちます。一方、 $A_{eq}$ の非零固有値 $\lambda_i\ (i=1,\cdots,n-m)$ に対応する固有ベクトルを $v_i$ とすると、

$\displaystyle{(5)\quad SA_{eq}=0\Rightarrow SA_{eq}v_i=0\Rightarrow \lambda_iSv_i=0 \Rightarrow Sv_i=0 }$

すなわち

$\displaystyle{(6)\quad { \boxed{S[v_1\cdots v_{n-m}]=0 %　\Rightarrow [v_1\cdots v_{n-m}]^TS^T=0}} }$

が成り立ちます。したがって、まず等価制御による閉ループ系の $A_{eq}$ が望ましい固有値・固有ベクトルをもつように設計して、(5)から $S$ を決定することが考えられます。

●２重積分器に対して、関数swfvplによるスイッチング関数Sの決定は次のように行われます。

MATLAB

%cDINT_smc2.m
%-----
clear all, close all
A=[0 1;
   0 0];
B=[0;1];
C=eye(2); 
D=zeros(2,1);
lambda=-1;
nocomp=0;
specpos=rand(2,1)
specent=rand(2,1)
S=swfvpl(A,B,lambda,nocomp,specpos,specent)
Phi=-1;
P2=lyap(Phi,1);
P2S=P2*S;
Leq=inv(S*B)*S*A;
LPhi=-inv(S*B)*Phi*S;
L=Leq+LPhi;
rho=0.5;
Ln=inv(S*B)*rho;
x0=B;
sim('DINT_smc')
%-----
%eof

これを実行して図4のシミュレーション結果を得ます。

図4　関数swfvplによるスイッチング関数Sの決定

ここで、関数swfvplは、Edwards and Spurgeon によるものです。

積分動作の導入

[6] 定値外乱を受ける制御対象

$\displaystyle{(1)\quad \begin{array}{l} \dot{x}(t)=Ax(t)+Bu(t)+w\\ y(t)=Cx(t)\\ (x(t)\in{\rm\bf R}^n, u(t)\in{\rm\bf R}^m, y(t)\in{\rm\bf R}^p, w\in{\rm\bf R}^n, m=p) \end{array} }$

の出力を、次のコマンド（定値目標）

$\displaystyle{(2)\quad r\in{\rm\bf R}^m %\begin{array}{l} %\dot{r}(t)=\Gamma(r(t)-r_c)\\ %(r(t),r_c\in{\rm\bf R}^m) %\end{array} %\Gammaは安定行列 }$

に追従させることを考えます。そのために、積分動作

$\displaystyle{(3)\quad \begin{array}{l} \dot{x}_I(t)=y(t)-r\\ %\dot{x}_I(t)=r(t)-y(t)\\ (x_I(t)\in{\rm\bf R}^m) \end{array} }$

を考え、次の拡大系を構成します。

偏差系E3：
$\displaystyle{(4)\quad \underbrace{ \frac{d}{dt} \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] }_{\dot{x}_{E3}(t)} = \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] }_{x_{E3}(t)} + \underbrace{ \left[\begin{array}{c} 0 \\ I_m \end{array}\right] }_{B_{E3}} {\dot u}(t) }$

を得ます。この偏差系E3が標準形となっていることに注意して、SMCを設計します。

(4)を、改めて次のように書きます。

$\displaystyle{(5)\quad \begin{array}{l} \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot x_2(t) \end{array}\right] }_{\dot{x}_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} A & B \\ 0 & 0 \end{array}\right] }_{A_{E3}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} + \underbrace{ \left[\begin{array}{c} 0\\ I_m \end{array}\right] }_{B_{E3}} {\dot u}(t)\\ (x_1(t)=x(t)-x_\infty, x_2(t)=u(t)-u_\infty) \end{array} }$

スイッチング関数として、次式を考えます。

$\displaystyle{(6)\quad s(t)= \underbrace{ \left[\begin{array}{cc} S_1 & S_2 \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} = \underbrace{S_2 \left[\begin{array}{cc} M & I \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} \ (M=S_2^{-1}S_1) }$

(15)に対して、座標変換

$\displaystyle{(7)\quad \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} I & 0 \\ S_1 & S_2 \\ \end{array}\right] }_{T_s} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)}\\ \Leftrightarrow \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} I & 0 \\ -S_2^{-1}S_1 & S_2^{-1} \\ \end{array}\right] }_{T_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} }$

を行って、次式を得ます。

$\displaystyle{(8)\quad \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot s(t) \end{array}\right] }_{\dot{x}'_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ \bar{A}_{21} & \bar{A}_{22} \\ \end{array}\right] }_{T_sA_{E3}T_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} + \underbrace{ \left[\begin{array}{cc} 0\\ S_2 \end{array}\right] }_{T_sB_{E3}} {\dot u}(t) }$

$\displaystyle{(9)\quad \left\{\begin{array}{l} \bar{A}_{11}=A-BM\\ \bar{A}_{12}=BS_2^{-1}\\ \bar{A}_{21}=S_1(A-BM)\\ \bar{A}_{22}=S_1BS_2^{-1} \end{array}\right. }$

以下では、 $\bar{A}_{11}$ が安定行列となるようにスイッチング関数が選ばれていると仮定します。

このとき、SM制御則を、２次安定性

$\displaystyle{(10)\quad \begin{array}{lll} V(s)=s(t)^TP_2s(t)& \Rightarrow \dot{V}(s)\le -s^T(t)s(t)&(t\le t_s)\\ V(x_1)=x_1^T(t)P_1x_1(t)& \Rightarrow \dot{V}(x_1)\le -x_1^T(t)Q_1x_1(t)&(t> t_s) \end{array} }$

が成り立つように決定します（ $P_1>0$ , $P_2>0$ , $Q_1>0$ ）。

偏差系E3に対するSMCは次式で与えられます。

$\displaystyle{(11)\quad {\dot u}(t)={\dot u}_\ell(t)+{\dot u}_n(t) }$

$\displaystyle{(12)\quad {\dot u}_\ell(t)=-(SB_{E3})^{-1}(SA_{E3}-\Phi S)x_{E3}(t) }$

$\displaystyle{(13)\quad {\dot u}_n(t)=-(SB_{E3})^{-1}\rho\frac{P_2s(t)}{||P_2s(t)||}}=-(SB_{E3})^{-1}\rho\,{\rm sgn}(P_2Sx_{E3}(t)) }$

これらを積分して、制御対象(1)に対する積分動作をもつSMCを導出します。

$\displaystyle{(14)\quad {u}(t)={u}_\ell(t)+{u}_n(t) }$

まず(12)は次式のように書けます。

$\displaystyle{(15)\quad {\dot u}_\ell(t)=- \boxed{\underbrace{ (SB_{E3})^{-1}(SA_{E3}-\Phi S)S_E^{-1} }_{\left[\begin{array}{cc} F & F_I \end{array}\right]}} \underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right] }_{x_{E2}(t)} }$

これを積分して

$\displaystyle{(16)\quad \boxed{u_\ell(t)=-Fx(t)+F_I\int_0^t(r-y(\tau))d\tau} }$

次に(13)は次式のように書けます。

$\displaystyle{(17)\quad {\dot u}_n(t) =-(SB_{E3})^{-1}\rho\, {\rm sgn}( \boxed{\underbrace{ P_2SS_E^{-1} }_{\left[\begin{array}{cc} G & G_I \end{array}\right]}} \underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right] }_{x_{E2}(t)}) }$

これを積分すれば

$\displaystyle{(18)\quad \boxed{u_n(t)=(SB_{E3})^{-1}\rho\,\int_0^t\underbrace{{\rm sgn}(-G{\dot x}(\tau)+G_I(r-y(\tau)))}_{0,\pm 1}d\tau} }$

●２重積分器に対して、次のように行われます。

MATLAB

%cDINT_smci.m
%-----
clear all, close all
A=[0 1;
   0 0];
B=[0;1];
C=[1 0]; 
D=0;
n=size(A,1); m=size(B,2); p=size(C,1); 
AE=[A B;zeros(m,n+m)];
BE=[zeros(n,m);eye(m)];
SE=[A B;C D];
QE=eye(n+m);
S=-swflqr(AE,BE,QE)
Phi=-0.5;
Leq=inv(S*BE)*S*AE;
LPhi=-inv(S*BE)*Phi*S;
L=(Leq+LPhi)/SE;
F=L(1:2)
FI=L(3)
P2=lyap(Phi,1)
P2S=P2*S/SE;
G=P2S(1:2)
GI=P2S(3) 
rho=1;
Ln=inv(S*BE)*rho
x0=zeros(2,1);
sim('DINT_smci')
%-----
%eof

これを実行して図5のシミュレーション結果を得ます。

図5　積分動作を導入したSM制御

補遺：2次安定系のロバスト性…Homework

いま漸近安定な線形系

$\displaystyle{(1)\quad \dot{x}(t)=Ax(t) }$

において次のようなモデルの不確かさがあるとします。

$\displaystyle{(2)\quad \dot{x}(t)=(A+\Delta A)x(t)=Ax(t)+\underbrace{\Delta Ax(t)}_{\xi(t,x)}\qquad(x(0)\ne0) }$

条件A2’より、適当な $P>0$ と $Q>0$ に対して

$\displaystyle{(3)\quad PA+A^TP+Q=0 }$

が成り立ちます。このとき

$\displaystyle{(4)\quad \boxed{||\xi(t,x)||\le\frac{1}{2}\frac{\sigma_n(Q)}{\sigma_1(P)}||x(t)||} }$

ならば、リャプノフ関数

$\displaystyle{(5)\quad V(x)=x^T(t)Px(t) }$

に対して

$\displaystyle{(6)\quad \frac{d}{dt}V(x)<0 }$

となって、リャプノフ安定となります。

●実際、

$\displaystyle{(7)\quad \begin{array}{l} \frac{d}{dt}V(x)=x^T(t)P\dot{x}(t)+\dot{x}^T(t)Px(t)\\ =x^T(t)P(Ax(t)+\xi(t,x))+(x^T(t)A^T+\xi^T(t,x))Px(t)\\ =x^T(t)PAx(t)+x^T(t)A^TPx(t)+x^T(t)P\xi(t,x)+\xi^T(t,x)Px(t)\\ =x^T(t)(PA+A^TP)x(t)+2\xi^T(t,x)Px(t)\\ =-x^T(t)Qx(t)+2\xi^T(t,x)Px(t)\\ \le -x^T(t)Qx(t)+ 2||\xi(t,x)||\underbrace{||Px(t)||}_{\sqrt{x^T(t)P^2x(t)}} \end{array} }$

ここで、

$\displaystyle{(8)\quad \begin{array}{l} \sigma_n(Q)x^T(t)x(t)\le x^T(t)Qx(t) \Rightarrow \sigma_n(Q)||x(t)||^2\le x^T(t)Qx(t)\\ x^T(t)P^2x(t)\le\sigma_1(P^2)x(t)^Tx(t) \Rightarrow \sqrt{x^T(t)P^2x(t)}\le\sigma_1(P)||x(t)|| \end{array} }$

に注意して、(4)の下では

$\displaystyle{(9)\quad \begin{array}{l} \frac{d}{dt}V(x) \le -\sigma_n(Q)||x(t)||^2+ 2||\xi(t,x)||\sigma_1(P)x(t)^Tx(t) \\ =(||\xi(t,x)||-\frac{1}{2}\frac{\sigma_n(Q)}{\sigma_1(P)}||x(t)||)2\sigma_1(P)||x(t)||<0 \end{array} }$

となります。

補遺：偏差系の安定化SM制御

投稿日時: 2021年10月29日投稿者: cacsd

補遺：偏差系の安定化SM制御…Homework

レギュレータ問題は平衡状態の安定化、追値問題は他の平衡状態への安定化を意味します。LQI制御のところでみたように、偏差系を導入すれば、追値問題は偏差系の安定化問題となります。この観点から積分動作導入による追従SM制御を見直してみます。

[1] 定値外乱を受ける制御対象

$\displaystyle{(1)\quad \begin{array}{l} \dot{x}(t)=Ax(t)+Bu(t)+w\\ y(t)=Cx(t)\\ (x(t)\in{\rm\bf R}^n, u(t)\in{\rm\bf R}^m, y(t)\in{\rm\bf R}^p, w\in{\rm\bf R}^n, m=p) \end{array} }$

の出力を、次のコマンド（定値目標）

$\displaystyle{(2)\quad r\in{\rm\bf R}^m %\begin{array}{l} %\dot{r}(t)=\Gamma(r(t)-r_c)\\ %(r(t),r_c\in{\rm\bf R}^m) %\end{array} %\Gammaは安定行列 }$

に追従させることを考えます。そのために、積分動作

$\displaystyle{(3)\quad \begin{array}{l} \dot{x}_I(t)=y(t)-r\\ %\dot{x}_I(t)=r(t)-y(t)\\ (x_I(t)\in{\rm\bf R}^m) \end{array} }$

を考え、次の拡大系を構成します。

$\displaystyle{(4)\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} }$

定常状態では

$\displaystyle{(5)\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_{\infty}$ , $x_{I\infty}$ , $u_{\infty}$ は定数ベクトル)。まず、(4)から(5)を引いて、つぎの偏差系を得ます。

偏差系E1：
$\displaystyle{(6)\quad \underbrace{ \frac{d}{dt} \left[\begin{array}{c} x(t)-x_\infty \\ x_I(t)-x_{I\infty} \end{array}\right] }_{\dot{x}_{E1}(t)} = \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_{E1}(t)} + \underbrace{ \left[\begin{array}{c} B \\ 0 \end{array}\right] }_{B_{E1}} (u(t)-u_\infty) }$

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

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

を得ます。さらに、(1)の状態方程式と観測方程式をまとめた

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

から、(5)すなわち

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

を引いて、つぎの関係式が成り立ちます。

$\displaystyle{(10)\quad \underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right] }_{x_{E2}(t)} = \underbrace{ \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right] }_{S_E} \underbrace{ \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] }_{x_{E3}(t)} }$

これを用いて、偏差系E2に座標変換を行えば

偏差系E3：
$\displaystyle{(11)\quad \underbrace{ \frac{d}{dt} \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] }_{\dot{x}_{E3}(t)} = \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] }_{x_{E3}(t)} + \underbrace{ \left[\begin{array}{c} 0 \\ I_m \end{array}\right] }_{B_{E3}} {\dot u}(t) }$

を得ます。ここで、つぎの関係式を用いました。

$\displaystyle{(12)\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_E} = \underbrace{ \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right] }_{S_E} \underbrace{ \left[\begin{array}{cc} A & B \\ 0 & 0 \end{array}\right] }_{A_{E3}} }$

$\displaystyle{(13)\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_E} \underbrace{ \left[\begin{array}{c} 0 \\ I_m \end{array}\right] }_{B_{E3}} }$

$\displaystyle{(14)\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_E} = \underbrace{ \left[\begin{array}{cc} C & 0 \end{array}\right] }_{C_{E3}} }$

以下では、この偏差系E3が標準形となっていることに注意して、SMCを設計します。

(11)を、改めて次のように書きます。

$\displaystyle{(15)\quad \begin{array}{l} \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot x_2(t) \end{array}\right] }_{\dot{x}_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} A & B \\ 0 & 0 \end{array}\right] }_{A_{E3}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} + \underbrace{ \left[\begin{array}{c} 0\\ I_m \end{array}\right] }_{B_{E3}} {\dot u}(t)\\ (x_1(t)=x(t)-x_\infty, x_2(t)=u(t)-u_\infty) \end{array} }$

スイッチング関数として、次式を考えます。

$\displaystyle{(16)\quad s(t)= \underbrace{ \left[\begin{array}{cc} S_1 & S_2 \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} = \underbrace{S_2 \left[\begin{array}{cc} M & I \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} \ (M=S_2^{-1}S_1) }$

(15)に対して、座標変換

$\displaystyle{(17)\quad \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} I & 0 \\ S_1 & S_2 \\ \end{array}\right] }_{T_s} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)}\\ \Leftrightarrow \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} I & 0 \\ -S_2^{-1}S_1 & S_2^{-1} \\ \end{array}\right] }_{T_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} }$

を行って、次式を得ます。

$\displaystyle{(18)\quad \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot s(t) \end{array}\right] }_{\dot{x}'_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ \bar{A}_{21} & \bar{A}_{22} \\ \end{array}\right] }_{T_sA_{E3}T_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} + \underbrace{ \left[\begin{array}{cc} 0\\ S_2 \end{array}\right] }_{T_sB_{E3}} {\dot u}(t) }$

$\displaystyle{(19)\quad \begin{array}{l} \bar{A}_{11}=A-BM\\ \bar{A}_{12}=BS_2^{-1}\\ \bar{A}_{21}=S_1(A-BM)\\ \bar{A}_{22}=S_1BS_2^{-1} \end{array} }$

以下では、 $\bar{A}_{11}$ が安定行列となるようにスイッチング関数が選ばれていると仮定します。

このとき、スライディングモード制御則（SM制御則、SMC則）

$\displaystyle{(20)\quad {\dot u}(t)=\underbrace{{\dot u}_\ell(t)}_{linear\ control}+\underbrace{{\dot u}_n(t)}_{switching\ component}} }$

を、２次安定性

$\displaystyle{(21)\quad \begin{array}{lll} V(\bar{x})= \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]^T }_{\bar{x}^T(t)} \underbrace{ \left[\begin{array}{cc} P_1 & 0\\ 0 & P_2 \end{array}\right] }_{P} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{\bar{x}(t)}\\ \Rightarrow \dot{V}(\bar{x})\le - \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]^T }_{\bar{x}^T(t)} \underbrace{ \left[\begin{array}{cc} Q_1 & 0\\ 0 & I \end{array}\right] }_{Q} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{\bar{x}(t)} \end{array}} }$

すなわち

$\displaystyle{(21')\quad \begin{array}{lll} V(\bar{x})= \underbrace{x_1^T(t)P_1x_1(t)}_{V(x_1)}+\underbrace{s^T(t)P_2s(t)}_{V(s)}\\ \Rightarrow \dot{V}(\bar{x})\le -x_1^T(t)Q_1x_1(t)-s^T(t)s(t) \end{array}} }$

が成り立つように決定します（ $P_1>0$ , $P_2>0$ , $Q_1>0$ ）。

[2] 可到達性の検討

等価制御は

$\displaystyle{(22)\quad \begin{array}{l} s(t)=0\Rightarrow\dot{s}(t)=0 \Rightarrow 0=\bar{A}_{21}x_1(t)+\bar{A}_{22}s(t)+S_2{\dot u}(t)\\ \Rightarrow {\dot u}_{eq}(t)=-\underbrace{S_2^{-1}}_{(SB_{E3})^{-1}} \underbrace{\left[\begin{array}{cc} \bar{A}_{21} & \bar{A}_{22} \\ \end{array}\right]x'_{E3}(t)}_{SA_{E3}x_{E3}(t)}} \end{array} }$

のように得られます。(20)の第１項 ${\dot u}_\ell$ は、この等価制御をベースして

$\displaystyle{(23)\quad \begin{array}{l} {\dot u}_\ell(t)=-\underbrace{S_2^{-1}}_{(SB_{E3})^{-1}} \underbrace{(\left[\begin{array}{cc} \bar{A}_{21} & \bar{A}_{22} \\ \end{array}\right]x'_{E3}(t)-\Phi s(t))}_{(SA_{E3}-\Phi S)x_{E3}(t)}} \end{array}} }$

のように構成します（ $\Phi$ は安定行列）。このとき閉ループ系は次式で与えられます。

$\displaystyle{(24)\quad \begin{array}{l} \left[\begin{array}{c} \dot{x}_1(t)\\ \dot{s}(t) \end{array}\right] = \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ 0 & \Phi \\ \end{array}\right] \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] + \left[\begin{array}{c} 0\\ S_2{\dot u}_n(t) \end{array}\right] \end{array }$

すなわち

$\displaystyle{(25)\quad \begin{array}{l} \dot{x}_1(t)=\bar{A}_{11}x_1(t)+\bar{A}_{12}s(t)\\ \dot{s}(t)=\Phi s(t)+S_2{\dot u}_n(t) \end{array} }$

ここで、 $\Phi$ は安定行列なので

$\displaystyle{(26)\quad \begin{array}{l} P_2\Phi+\Phi^TP_2=-I \end{array} }$

を満たす $P_2>0$ を選ぶことができます。これを用いて

$\displaystyle{(27)\quad {\dot u}_n(t)=-\underbrace{S_2^{-1}}_{(SB_{E3})^{-1}}\rho\frac{P_2s(t)}{||P_2s(t)||}} }$

と選びます（ $\rho>0$ は定数）。このとき次式が成り立ちます。

$\displaystyle{(28)\quad \begin{array}{l} \dot{V}(s)=2s^T(t)P_2\dot{s}(t)\\ =2s^T(t)P_2(\Phi s(t)-\rho\frac{P_2s(t)}{||P_2s(t)||})\\ =s^T(t)(P_2\Phi+\Phi^TP_2)s(t)+2s^T(t)P_2(-\rho\frac{P_2s(t)}{||P_2s(t)||})\\ = -||s(t)||^2-2\rho||P_2s(t)||\\ \le -||s(t)||^2 \end{array} }$

[3] スライディングモードの検討

$\bar{A}_{11}$ は安定行列なので

$\displaystyle{(29)\quad \begin{array}{l} P_1\bar{A}_{11}+\bar{A}_{11}^TP_1=-Q_1<0 \end{array} }$

を満たす $P_1>0$ を選ぶことができます。

$\displaystyle{(30)\quad \begin{array}{l} \dot{V}(x_1)=2x_1^T(t)P_1\dot{x}_1(t)\\ =2x_1^T(t)P_1( \bar{A}_{11}x_1(t)+\bar{A}_{12}\underbrace{s(t)}_{0})\\ =x_1^T(t)(P_1\bar{A}_{11}+\bar{A}_{11}^TP_1)x_1(t)\\ =-x_1^T(t)Q_1x_1(t) \end{array} }$

[4] 積分動作をもつSMC

上で求めた偏差系E3に対するSMCは次式で与えられました。

$\displaystyle{(31)\quad {\dot u}(t)={\dot u}_\ell(t)+{\dot u}_n(t) }$

$\displaystyle{(32)\quad {\dot u}_\ell(t)=-(SB_{E3})^{-1}(SA_{E3}-\Phi S)x_{E3}(t) }$

$\displaystyle{(33)\quad {\dot u}_n(t)=-(SB_{E3})^{-1}\rho\frac{P_2s(t)}{||P_2s(t)||}}=-(SB_{E3})^{-1}\rho\,{\rm sgn}(P_2Sx_{E3}(t)) }$

これらを積分して、制御対象(1)に対する積分動作をもつSMCを導出します。

$\displaystyle{(34)\quad {u}(t)={u}_\ell(t)+{u}_n(t) }$

まず(32)は(10)を用いて次式のように書けます。

$\displaystyle{(35)\quad {\dot u}_\ell(t)=- \underbrace{ (SB_{E3})^{-1}(SA_{E3}-\Phi S)S_E^{-1} }_{\left[\begin{array}{cc} F & F_I \end{array}\right]} \underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right] }_{x_{E2}(t)} }$

これを積分して

$\displaystyle{(36)\quad u_\ell(t)=-Fx(t)+F_I\int_0^t(r-y(\tau))d\tau }$

次に(33)は(10)を用いて次式のように書けます。

$\displaystyle{(37)\quad {\dot u}_n(t) =-S_2^{-1}\rho\, {\rm sgn}( \underbrace{ P_2SS_E^{-1} }_{\left[\begin{array}{cc} G & G_I \end{array}\right]} \underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right] }_{x_{E2}(t)}) }$

これを積分すれば

$\displaystyle{(38)\quad u_n(t)=S_2^{-1}\rho\,\int_0^t\underbrace{{\rm sgn}(-G{\dot x}(\tau)+G_I(r-y(\tau)))}_{0,\pm 1}d\tau }$

ここで、 $\dot{x}_1(t)=\frac{d}{dt}(x(t)-x_\infty)=\dot{x}(t)$ に注意し、(25)の第１式を用いて $\dot{x}(t)$ は

$\displaystyle{(39)\quad \begin{array}{l} \underbrace{\dot{x}_1(t)}_{\dot{x}(t)}=\underbrace{(A-BM)}_{\bar{A}_{11}}x_1(t)+\underbrace{BS_2^{-1}}_{\bar{A}_{12}}\underbrace{(S_1x_1(t)+S_2x_2(t))}_{s(t)}\\ =(A-BM)x_1(t)+BS_2^{-1}S_2(Mx_1(t)+x_2(t))\\ =Ax_1(t)+Bx_2(t)\\ =Ax(t)+Bu(t)-\left[\begin{array}{cc} A & B \end{array}\right] \left[\begin{array}{c} x_\infty \\ u_\infty \end{array}\right]\\ =Ax(t)+Bu(t)-\left[\begin{array}{cc} A & B \end{array}\right]S_E^{-1} \left[\begin{array}{c} -w \\ r \end{array}\right]\\ =Ax(t)+Bu(t)-\left[\begin{array}{cc} I & 0 \end{array}\right] \left[\begin{array}{c} -w \\ r \end{array}\right]\\ =Ax(t)+Bu(t)+w \end{array} }$

となって元の状態方程式となりますが、 $w$ を無視し、 $u$ の近似値を使うことも一手段かもしれません。

[5] 数値例(1)

$\displaystyle{(101)\quad \begin{array}{l} \dot{x}(t)=\underbrace{0}_{a}x(t)+\underbrace{1}_{b}u(t)+w\\ y(t)=\underbrace{1}_{c}x(t)\\ \end{array} }$

$\displaystyle{(109)\quad \left[\begin{array}{c} -w \\ r \end{array}\right] = \underbrace{ \left[\begin{array}{cc} 0 & 1 \\ 1 & 0 \end{array}\right] }_{S_E} \left[\begin{array}{c} x_\infty \\ u_\infty \end{array}\right] \Rightarrow \left[\begin{array}{c} x_\infty \\ u_\infty \end{array}\right]= \left[\begin{array}{c} r \\ -w \end{array}\right] }$

$\displaystyle{(111)\quad \underbrace{ \frac{d}{dt} \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] }_{\dot{x}_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} 0 & 1 \\ 0 & 0 \end{array}\right] }_{A_{E3}} \underbrace{ \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] }_{x_{E3}(t)} + \underbrace{ \left[\begin{array}{c} 0 \\ 1 \end{array}\right] }_{B_{E3}} {\dot u}(t) }$

$\displaystyle{(116)\quad s(t)= \underbrace{ \left[\begin{array}{cc} s_1 & s_2 \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} = \underbrace{s_2 \left[\begin{array}{cc} m & 1 \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} \ (m=s_2^{-1}s_1) }$

$\displaystyle{(117)\quad \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} 1 & 0 \\ s_1 & s_2 \\ \end{array}\right] }_{T_s} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)}\\ \Leftrightarrow \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} 1 & 0 \\ -s_2^{-1}s_1 & s_2^{-1} \\ \end{array}\right] }_{T_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} }$

$\displaystyle{(118)\quad \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot s(t) \end{array}\right] }_{\dot{x}'_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} -m & s_2^{-1} \\ -s_1m & m \\ \end{array}\right] }_{T_sA_{E3}T_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} + \underbrace{ \left[\begin{array}{cc} 0\\ s_2 \end{array}\right] }_{T_sB_{E3}} {\dot u}(t) }$

$\displaystyle{(123)\quad \begin{array}{l} {\dot u}_\ell(t)= -\underbrace{s_2^{-1}}_{(SB_{E3})^{-1}} \underbrace{(\left[\begin{array}{cc} -s_1m & m \\ \end{array}\right]x'_{E3}(t)-\Phi s(t))}_{(SA_{E3}-\Phi S)x_{E3}(t)}}\\ =-s_2^{-1}(\left[\begin{array}{cc} 0 & ms_2 \\ \end{array}\right]- \left[\begin{array}{cc} \Phi s_1 & \Phi s_2 \\ \end{array}\right])x_{E3}(t)\\ =-s_2^{-1}\left[\begin{array}{cc} -\Phi s_1 & ms_2-\Phi s_2 \\ \end{array}\right]x_{E3}(t)\\ =-\left[\begin{array}{cc} -\Phi m & m-\Phi \\ \end{array}\right] \underbrace{ \left[\begin{array}{cc} 0 & 1 \\ 1 & 0 \end{array}\right] }_{S_E^{-1}} \underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right] }_{x_{E2}(t)}\\ =-\underbrace{(m-\Phi)}_{f}{\dot x}(t)+\underbrace{m(-\Phi)}_{f_I}(r-y(t))\\ \Rightarrow u_\ell(t)=-fx(t)+f_I\int_0^t (r-y(\tau))d\tau \end{array} }$

$\displaystyle{(126)\quad \begin{array}{l} p_2\Phi+\Phi^Tp_2=-1\Rightarrow p_2=-\frac{1}{2}\Phi^{-1} \end{array} }$

$\displaystyle{(127)\quad \begin{array}{l} {\dot u}_n(t) =-s_2^{-1}\rho\, {\rm sgn}(\underbrace{-\frac{1}{2}\Phi^{-1}}_{p_2} \underbrace{ \left[\begin{array}{cc} s_1 & s_2 \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{cc} 0 & 1 \\ 1 & 0 \end{array}\right] }_{S_E^{-1}} \underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right] }_{x_{E2}(t)})\\ =\rho\, {\rm sgn}(\frac{1}{2}\Phi^{-1} \left[\begin{array}{cc} 1 & m \end{array}\right] \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right])\\ =\rho\, {\rm sgn}(-\underbrace{\frac{1}{2}(-\Phi^{-1})}_{g}{\dot x}(t)+\underbrace{\frac{1}{2}m(-\Phi^{-1})}_{g_I} (r-y(t)))\\ \Rightarrow u_n(t)=\rho\,\int_0^t\underbrace{{\rm sgn}(-g{\dot x}(\tau)+g_I(r-y(\tau)))}_{0,\pm 1}d\tau \end{array} }$

$m=1$ , $\Phi=-0.5$ の場合のシミュレーション結果を次に示します。

[5] 数値例(2)

$\displaystyle{(201)\quad \begin{array}{l} \dot{x}(t)= \underbrace{\left[\begin{array}{cc} 0 & 1 \\ 0 & 0 \end{array}\right]}_{A}x(t) +\underbrace{\left[\begin{array}{cc} 0 \\ 1 \end{array}\right]}_{B}u(t)+ \left[\begin{array}{cc} 0 \\ w \end{array}\right]\\ y(t)=\underbrace{\left[\begin{array}{cc} 1 & 0 \end{array}\right]}_{C}x(t)\\ \end{array} }$

$\displaystyle{(209)\quad \left[\begin{array}{c} 0 \\ -w \\ r \end{array}\right] = \underbrace{ \left[\begin{array}{ccc} 0 & 1 & 0 \\ 0 & 0 & 1 \\ 1 & 0 & 0 \end{array}\right] }_{S_E} \left[\begin{array}{c} x_{1\infty} \\ x_{2\infty} \\ u_\infty \end{array}\right] \Rightarrow \left[\begin{array}{c} x_{1\infty} \\ x_{2\infty} \\ u_\infty \end{array}\right]= \left[\begin{array}{c} r \\ 0 \\ -w \end{array}\right] }$

$\displaystyle{(211)\quad \underbrace{ \frac{d}{dt} \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] }_{\dot{x}_{E3}(t)} = \underbrace{ \left[\begin{array}{ccc} 0 & 1 & 0 \\ 0 & 0 & 1 \\ 0 & 0 & 0 \end{array}\right] }_{A_{E3}} \underbrace{ \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] }_{x_{E3}(t)} + \underbrace{ \left[\begin{array}{c} 0 \\ 0 \\ 1 \end{array}\right] }_{B_{E3}} {\dot u}(t) }$

$\displaystyle{(216)\quad \begin{array}{l} s(t)= \underbrace{ \left[\begin{array}{ccc} s_{11} & s_{12} & s_2 \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} = \underbrace{s_2 \left[\begin{array}{ccc} m_1 & m_2 & 1 \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)}\\ (m_1=s_2^{-1}s_{11},m_2=s_2^{-1}s_{12}) \end{array} }$

$\displaystyle{(217)\quad \begin{array}{l} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} = \underbrace{ \left[\begin{array}{ccc} 1 & 0 & 0\\ 0 & 1 & 0\\ s_{11} & s_{12} & s_2 \\ \end{array}\right] }_{T_s} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)}\Leftrightarrow\\ \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} = \underbrace{ \left[\begin{array}{ccc} 1 & 0 & 0\\ 0 & 1 & 0\\ -s_2^{-1}s_{11} & -s_2^{-1}s_{12} & s_2^{-1} \\ \end{array}\right] }_{T_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} \end{array} }$

$\displaystyle{(218)\quad \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot s(t) \end{array}\right] }_{\dot{x}'_{E3}(t)} = \underbrace{ \left[\begin{array}{ccc} 0 & 1 & 0\\ -m_1 & -m_2 & s_2^{-1} \\ -s_{12}m_1 & -s_{12}m_2 & m_2 \end{array}\right] }_{T_sA_{E3}T_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} + \underbrace{ \left[\begin{array}{cc} 0\\ 0\\ s_2 \end{array}\right] }_{T_sB_{E3}} {\dot u}(t) }$

$\displaystyle{(223)\quad \begin{array}{l} {\dot u}_\ell(t)= -\underbrace{s_2^{-1}}_{(SB_{E3})^{-1}} \underbrace{(\left[\begin{array}{ccc} -s_{12}m_1 & -s_{12}m_2 & m_2 \\ \end{array}\right]x'_{E3}(t)-\Phi s(t))}_{(SA_{E3}-\Phi S)x_{E3}(t)}}\\ =-s_2^{-1}(\left[\begin{array}{ccc} 0 & 0 & m_2s_2 \\ \end{array}\right]- \left[\begin{array}{ccc} \Phi s_{11} & \Phi s_{12} & \Phi s_2 \\ \end{array}\right])x_{E3}(t)\\ =-s_2^{-1}\left[\begin{array}{ccc} -\Phi s_{11} & -\Phi s_{12} & m_2s_2-\Phi s_2 \\ \end{array}\right]x_{E3}(t)\\ =-\left[\begin{array}{ccc} -\Phi m_1 & -\Phi m_2 & m_2-\Phi \\ \end{array}\right] \underbrace{ \left[\begin{array}{ccc} 0 & 0 & 1 \\ 1 & 0 & 0 \\ 0 & 1 & 0 \end{array}\right] }_{S_E^{-1}} \underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right] }_{x_{E2}(t)}\\ =-\underbrace{ \left[\begin{array}{cc} -\Phi m_2 & m_2-\Phi \\ \end{array}\right] }_{F}{\dot x}(t)+\underbrace{m_1(-\Phi)}_{F_I}(r-y(t))\\ \Rightarrow u_\ell(t)=-Fx(t)+F_I\int_0^t (r-y(\tau))d\tau \end{array} }$

$\displaystyle{(226)\quad \begin{array}{l} p_2\Phi+\Phi^Tp_2=-1\Rightarrow p_2=-\frac{1}{2}\Phi^{-1} \end{array} }$

$\displaystyle{(237)\quad \begin{array}{l} {\dot u}_n(t) =-s_2^{-1}\rho\, {\rm sgn}(\underbrace{-\frac{1}{2}\Phi^{-1}}_{p_2} \underbrace{ \left[\begin{array}{ccc} s_{11} & s_{12} & s_2 \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{ccc} 0 & 0 & 1 \\ 1 & 0 & 0 \\ 0 & 1 & 0 \end{array}\right] }_{S_E^{-1}} \underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right] }_{x_{E2}(t)})\\ =s_2^{-1}\rho\, {\rm sgn}(\frac{1}{2}\Phi^{-1} \left[\begin{array}{ccc} s_{12} & s_2 & s_{11} \end{array}\right] \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right])\\ =\rho\, {\rm sgn}(\frac{1}{2}\Phi^{-1} \left[\begin{array}{ccc} m_2 & 1 & m_1 \end{array}\right] \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right])\\ =\rho\, {\rm sgn}(-\underbrace{\frac{1}{2}(-\Phi^{-1}) \left[\begin{array}{cc} m_2 & 1 \end{array}\right] }_{G}{\dot x}(t)+\underbrace{\frac{1}{2}m_1(-\Phi^{-1})}_{G_I} (r-y(t)))\\ \Rightarrow u_n(t)=\rho\,\int_0^t\underbrace{{\rm sgn}(-G{\dot x}(\tau)+G_I(r-y(\tau)))}_{0,\pm 1}d\tau \end{array} }$

$m_1=m_2=1$ , $\Phi=-0.5$ の場合のシミュレーション結果を次に示します。

演習…Flipped Classroom

MATLAB

%cCIP_smci2.m
%-----
 clear all, close all
 global mc m ell g th0
 mc=1; m=0.1; ell=0.2; g=9.8;
 ths=0 %input('ths   = <0,180> ')/180*pi;
 th0=3/180*pi %input('th(0) = <0,180> ')/180*pi;
%-----
 A=[zeros(2,2) eye(2);zeros(2,4)];
 A(3,1)=0; 
 A(3,2)=-3*m*g/(m+4*mc); 
 A(4,1)=0; 
 A(4,2)=3*(m+mc)*g/((m+4*mc)*ell);  
 B=zeros(4,1);
 B(3)=4/(m+4*mc);
 B(4)=-3/((m+4*mc)*ell); 
%-----
 Tcart=1; Mr=0.5;
 Tpend=1/4*2*pi*sqrt(4/3*ell/g); Mth=3/180*pi;
 Tamp=0.01; Mamp=10; 
 CM=eye(2,4); DM=zeros(2,1);
 CS=[1 0]; C=CS*CM; D=CS*DM;
 n=size(A,1); m=size(C,1); r=size(B,2);
%-----
 AE=[A B;zeros(m,n+m)];
 BE=[zeros(n,m);eye(m)];
 SE=[A B;C D];
 CE=eye(n+m);
 QE=diag([1/Mr 1/Mth Tcart/Mr Tpend/Mth 1/Mamp].^2);
 S=swflqr(AE,BE,QE)
%-----
 Phi=-0.5;
 P2=0.5*inv(-Phi);
 Check=P2*Phi+Phi*P2
 P2S=P2*S/SE;
 G=P2S(1:4)
 GI=P2S(5) 
 Leq=inv(S*BE)*S*AE;
 LPhi=-inv(S*BE)*Phi*S;
 L=(Leq+LPhi)/SE;
 F=L(1:4)
 FI=L(5)
 rho=1;
 Ln=inv(S*BE)*rho
%-----
 x0=[0;th0*0;0;0];
 sim('CIP_smci2_2015a.mdl')
%-----
%eof

SCILAB

応用例

投稿日時: 2021年10月29日投稿者: cacsd

応用例…Homework

[1]　ある航空機の線形状態方程式として、次を考えます。

$\displaystyle{(1)\quad \begin{array}{l} \dot{x}(t)=Ax(t)+Bu(t)\\ x(t)= \left[\begin{array}{ll} \theta & pitch\ angle\\ q & pitch\ rate\\ \alpha & angle\ of\ atack\\ \eta & elevator\ deflection\\ \delta & flap\ deflection \end{array}\right],\quad u(t)= \left[\begin{array}{ll} \eta_c & elevator\ command\\ \delta_c& flap\ comannd \end{array}\right]\\ A=\left[\begin{array}{rrrrr} 0 & 1 & 0 & 0 & 0\\ 0 &-1.99 & -13.41 & -18.95 & -3.60\\ 0 & 1 & -1.74 & -0.08 & -0.59\\ 0 & 0 & 0 & -20 & 0\\ 0 & 0 & 0 & 0 & -20 \end{array}\right],\quad B=\left[\begin{array}{rr} 0 & 0 \\ 0 & 0 \\ 0 & 0 \\ 20 & 0\\ 0 & 20 \end{array}\right] \end{array} }$

$A$ 行列の固有値は次のように求められます。

$\displaystyle{(2)\quad \lambda(A):\left\{\begin{array}{cl} -1.864\pm j 3.660 & short\ period\ mode\\ 0 & pitch\ attitude\ mode\\ -20 & elevator\ actuator\ mode\\ -20 & flap\ actuator\ mode \end{array}\right. }$

状態変数として、pitch angle $\theta$ の代わりに、fligt path angle $\gamma=\theta-\alpha$ を用いると、次式となります。

$\displaystyle{(3)\quad \begin{array}{l} \dot{x}(t)=Ax(t)+Bu(t)\\ x(t)= \left[\begin{array}{ll} \gamma=\theta-\alpha & fligt\ path\ angle\\ q & pitch\ rate\\ \alpha & angle\ of\ atack\\ \eta & elevator\ deflection\\ \delta & flap\ deflection \end{array}\right],\quad u(t)= \left[\begin{array}{ll} \eta_c & elevator\ command\\ \delta_c& flap\ comannd \end{array}\right]\\ A=\left[\begin{array}{rrrrr} 0 & 0 & 1.74 & 0.08 & 0.59\\ 0 &-1.99 & -13.41 & -18.95 & -3.60\\ 0 & 1 & -1.74 & -0.08 & -0.59\\ 0 & 0 & 0 & -20 & 0\\ 0 & 0 & 0 & 0 & -20 \end{array}\right],\quad B=\left[\begin{array}{rr} 0 & 0 \\ 0 & 0 \\ 0 & 0 \\ 20 & 0\\ 0 & 20 \end{array}\right] \end{array} }$

また、評価変数は、pitch angle $\theta$ とfligt path angle $\gamma$ とすると、その出力方程式は次式のように与えられます。

$\displaystyle{(4)\quad \begin{array}{l} y(t)=Hx(t)\\ H= \left[\begin{array}{ccccc} 1 & 0 & 1 & 0 & 0\\ 1 & 0 & 0 & 0 & 0 \end{array}\right] \end{array} }$

制御目的は、次図のように、pitch angle $\theta$ とfligt path angle $\gamma$ の間の非干渉化を達成し、独立して設定値に合せることです。

図1 航空機の操縦法

[3] 積分動作導入による追従制御

この場合の制御則は次のように構成されます。

$\displaystyle{(20)\quad { \boxed{\begin{array}{l} u(t)=u_\ell(t)+u_n(t)\\ u_\ell(t)=-\underbrace{(SB_E)^{-1}(SA_E-\Phi S)}_{L=L_{eq}+L_\Phi}\left[\begin{array}{c} x_r(t)\\ x(t) \end{array}\right]\\ -\underbrace{(SB_E)^{-1}(\Phi S_r+S_1B_r)}_{L_r} r(t) +\underbrace{(SB_E)^{-1}S_r}_{L_{\dot r}} \dot{r}(t)\\ u_n(t)=-\underbrace{(SB_E)^{-1}\rho(t,x)}_{L_n}\frac{P_2(s(t)-S_rr(t))}{||P_2(s(t)-S_rr(t))||} \end{array}}} }$

●ここでは、等価制御による閉ループ系の $A_{eq}$ が望ましい固有値・固有ベクトルをもつように設計して、次式からスイッチング関数を決める行列 $S$ を決定します。

$\displaystyle{(21)\quad { S[v_1\cdots v_{n-m}]=0 %　\Rightarrow [v_1\cdots v_{n-m}]^TS^T=0} }$

$A_{eq}$ の望ましい固有値を

$\displaystyle{(22)\quad \lambda(A_{eq}):\left\{\begin{array}{cl} -5.6\pm j 4.2 & short\ period\ mode\\ -1 & pitch\ attitude\ mode\\ -0.4 & elevator\ deflection\ mode\\ -0.7 & flap\ deflextion\ mode \end{array}\right. }$

とし、また望ましいモード分布を次により指定します。

$\displaystyle{(23)\quad {\tt specpos}= \left[\begin{array}{ccccc} 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0\\\hline 1 & 1 & 1 & 1 & 1\\ 1 & 0 & 1 & 1 & 1\\ 0 & 1 & 1 & 0 & 1\\\hline 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0\\ \end{array}\right],\quad {\tt specent}= \left[\begin{array}{rrrrr} -1 & -1 & -1 & -1 & -1\\ -1 & -1 & -1 & -1 & -1\\\hline 0 & 0 & 1 & 0 & 1\\ 1 & -1 & 0 & 1 & 0\\ -1 & 1 & 0 & -1 & 0\\\hline -1 & -1 & -1 & -1 & -1\\ -1 & -1 & -1 & -1 & -1\\ \end{array}\right] }$

このとき次のような行列 $S$ を得ます。

$\displaystyle{(24)\quad S=\left[\begin{array}{rrrrr|rr} 0.0531 & 0.0137 & -0.1436 & -0.0260 & -0.1373 & 0.0500 & 0\\ -0.0072 & -0.0612 & 0.1635 & 0.0035 & 0.1661 & 0 & 0.0500 \end{array}\right] }$

このときSMC則における行列 $L$ などは次のように計算されます。

$\displaystyle{(25)\quad L=\left[\begin{array}{rrrrrrr} -1.0616 & -0.2743 & 2.9379 & 0.6060 & 2.4605 & -0.4927 & -0.0900\\ 0.1440 & 1.2236 & -3.3390 & -0.2296 & -3.2769 & 0.0671 & 0.0142 \end{array}\right] }$

$\displaystyle{(26)\quad L_{r}=\left[\begin{array}{rr} -2.9499 & 0.0120\\ 0.4000 & 2.9391 \end{array}\right] }$

$\displaystyle{(27)\quad L_{\dot{r}}=\left[\begin{array}{rr} -0.1448 & 0.0013\\ 0.0196 & 0.1439 \end{array}\right] }$

図2 積分動作導入によるSMC制御系応答

MATLAB

%cAITCRAFT_smci.m
%-----
 clear all, close all
 A0= [0 0 1.74 0.08 0.59;
     0 -1.99 -13.41 -18.95 -3.6;
     0 1 -1.74 -0.08 -0.59;
     0 0 0 -20 0;
     0 0 0 0 -20];
 B0=[0 0;0 0;0 0;20 0;0 20]; 
 H=[1 0 1 0 0;1 0 0 0 0];
 A=[zeros(2,2) -H ;
    zeros(5,2) A0];
 B=[zeros(2,2); B0];
 Gamma=diag([-0.9,-0.7]);
%----- 
 lambda=[-5.6 4.2 -1 -0.4 -0.7];
 nocomp=1;
 specpos=[zeros(2,5);[1 1 1 1 1;1 0 1 1 1;0 1 1 0 1];zeros(2,5)]
 specent=[-ones(2,5);[0 0 1 0 1;1 -1 0 1 0;-1 1 0 -1 0];-ones(2,5)]
 S=swfvpl(A,B,lambda,nocomp,specpos,specent);
 SS=S(:,6:7)\S*0.05
 F=(SS*B)\(SS*A)
 eig(A-B*F)
%-----
 [nn,mm]=size(B);
 S1=SS(:,1:nn-mm);
 S2=SS(:,nn-mm+1:nn);
 Br=eye(nn-mm,mm);
 B2=B(nn-mm+1:nn,:); 
 Lambda=S2*B2;
 P2=diag([0.025,0.025]);
 Phi=-20/0.025*P2;
 Check=P2*Phi+Phi*P2
 L=inv(Lambda)*(SS*A-Phi*SS)
%-----
%Sr=[-0.1448 0.0013;0.0196 0.1439]
 Sr=zeros(2)
 Lr=inv(Lambda)*(Phi*Sr+S1*Br)
 Ldr=inv(Lambda)*Sr
 Ln=0.1*inv(Lambda)
%-----
 sim('AIRCRAFT_smci_2015a.mdl')
%-----
%eof

SCILAB

追従SMI制御

投稿日時: 2021年10月29日投稿者: cacsd

積分動作導入による追従SM制御…Homework

[1] 制御対象

$\displaystyle{(1)\quad \begin{array}{l} \dot{x}(t)=Ax(t)+Bu(t)\\ y(t)=Cx(t)\\ (x(t)\in{\rm\bf R}^n, u(t)\in{\rm\bf R}^m, y(t)\in{\rm\bf R}^p, m=p) \end{array} }$

の出力を、コマンド（次式の解）

$\displaystyle{(2)\quad \begin{array}{l} %\dot{r}(t)=\Gamma(r(t)-r_c)\\ \frac{d}{dt}(r(t)-r_c)=\Gamma(r(t)-r_c)\\ (r(t),r_c\in{\rm\bf R}^m) \end{array} }$

に追従させることを考えます（ $\Gamma$ は安定行列）。そのために、積分動作

$\displaystyle{(3)\quad \begin{array}{l} \dot{x}_r(t)=r(t)-y(t)\\ (x_r(t)\in{\rm\bf R}^m) \end{array} }$

を導入し、次の拡大系を構成します。ここで、(1)はすでに標準形であるとしています。

$\displaystyle{(4)\quad \begin{array}{l} \left[\begin{array}{c} \dot x_r(t)\\ \dot x(t) \end{array}\right] = \left[\begin{array}{c|cc} 0 & -C_1 & -C_2\\\hline 0 & A_{11} & A_{12} \\ 0 & A_{21} & A_{22} \end{array}\right] \left[\begin{array}{c} x_r(t)\\ x(t) \end{array}\right] + \left[\begin{array}{c} 0\\\hline 0\\ B_2 \end{array}\right] u(t) + \left[\begin{array}{c} I_m \\\hline 0 \\ 0 \end{array}\right] r(t)\\ (x_r(t)\in{\rm\bf R}^m, x(t)\in{\rm\bf R}^n) \end{array} }$

これを、次のように分割し直しても標準形であることには変わりありません。

$\displaystyle{(5a)\quad \begin{array}{l} \left[\begin{array}{c} \dot{x}_1(t)\\ \dot{x}_2(t) \end{array}\right] = \underbrace{ \left[\begin{array}{cc|c} 0 & -C_1 & -C_2\\ 0 & A_{11} & A_{12} \\\hline 0 & A_{21} & A_{22} \end{array}\right] }_{A_E} \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] + \underbrace{ \left[\begin{array}{c} 0\\ 0\\\hline B_2 \end{array}\right] }_{B_E} u(t) + \left[\begin{array}{c} I_m \\ 0 \\\hline 0 \end{array}\right] r(t)\\ (x_1(t)\in{\rm\bf R}^n, x_2(t)\in{\rm\bf R}^m) \end{array} }$

ただし

$\displaystyle{(5b)\quad %\begin{array}{l} \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] = \left[\begin{array}{c} x_r(t)\\ x(t) \end{array}\right] %\end{array} }$

●この積分器による拡大系を安定化できれば、積分器の値 $x_r(t)$ は定値となり、被積分項 $r(t)-y(t)$ の値は零となり、 $y(t)$ は $r(t)$ へ漸近します。そこで、SM制御によって拡大系を安定化し、追従制御系を構成することを考えます。この制御系は特別な $r(t)=0$ の場合を含みますので、まずスイッチング関数として、次式を考えます。

$\displaystyle{(6)\quad s(t)= \underbrace{ \left[\begin{array}{cc} S_1 & S_2 \\ \end{array}\right] }_{S} %\underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] %}_{x(t)} = \underbrace{S_2 \left[\begin{array}{cc} M & I_m \\ \end{array}\right] }_{S} %\underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] %}_{x(t)} \ (M=S_2^{-1}S_1) }$

(5)に対して、座標変換

$\displaystyle{(7)\quad \begin{array}{l} \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] = \underbrace{ \left[\begin{array}{cc} I_n & 0 \\ S_1 & S_2 \\ \end{array}\right] }_{T_s} \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right]\\ \Leftrightarrow \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] = \underbrace{ \left[\begin{array}{cc} I_n & 0 \\ -S_2^{-1}S_1 & S_2^{-1} \\ \end{array}\right] }_{T_s^{-1}} \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] \end{array} }$

を行って、次式を得ます。

$\displaystyle{(8a)\quad \begin{array}{l} %\underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot s(t) \end{array}\right] %}_{\dot{x}'(t)} = \underbrace{ \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ S_2\bar{A}_{21} & S_2\bar{A}_{22}S_2^{-1} \\ \end{array}\right] }_{T_sA_ET_s^{-1}} %\underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] %}_{x'(t)} + \underbrace{ \left[\begin{array}{cc} 0\\ S_2B_2 \end{array}\right] }_{T_sB_E} u(t)\\ + \left[\begin{array}{cc} B_r \\ S_1B_r \end{array}\right] r(t) \end{array} }$

ただし

$\displaystyle{(8b)\quad \left\{\begin{array}{l} \bar{A}_{11}= \left[\begin{array}{cc} 0 & -C_1 \\ 0 & A_{11} \end{array}\right] -\left[\begin{array}{c} -C_2\\ A_{12} \end{array}\right]M\\ \bar{A}_{12}= \left[\begin{array}{c} -C_2\\ A_{12} \end{array}\right]S_2^{-1}\\ \bar{A}_{21}=S_2(M\bar{A}_{11} + \left[\begin{array}{cc} 0 & A_{21} \end{array}\right] -A_{22}M)\\ \bar{A}_{22}=S_2(M \left[\begin{array}{c} -C_2\\ A_{12} \end{array}\right] +A_{22})S_2^{-1}\\ B_r=\left[\begin{array}{cc} I_m \\ 0 \end{array}\right] \end{array}\right. }$

ここで、 $\bar{A}_{11}$ が安定行列となるように行列 $M$ が選ばれているとします。

●特別な $r(t)=0$ の場合のスライディングモードは $s(t)=0$ で表されますが、一般の $r(t)\ne0$ の場合のスライディングモードは

$\displaystyle{(9)\quad {\boxed{s(t)-S_rr(t)=0} \quad\Rightarrow \dot{s}(t)-S_r\dot{r}(t)=0 } }$

で表されるとします。ここで、 $S_r$ は適当に選択された行列です。

●以上の準備の下で、SM制御則

$\displaystyle{(10)\quad u(t)=\underbrace{u_\ell(t)}_{linear\ control}+\underbrace{u_n(t)}_{switching\ component} }$

を、２次安定性

$\displaystyle{(11)\quad \boxed{ \begin{array}{lll} V(\bar{x})= \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t)-S_rr(t) \end{array}\right]^T }_{\bar{x}^T(t)} \underbrace{ \left[\begin{array}{cc} P_1 & 0\\ 0 & P_2 \end{array}\right] }_{P} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t)-S_rr(t) \end{array}\right] }_{\bar{x}(t)}\\ \Rightarrow \dot{V}(\bar{x})\le - \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t)-S_rr(t) \end{array}\right]^T }_{\bar{x}^T(t)} \underbrace{ \left[\begin{array}{cc} Q_1 & 0\\ 0 & I \end{array}\right] }_{Q} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t)-S_rr(t) \end{array}\right] }_{\bar{x}(t)} \end{array}}} }$

すなわち

$\displaystyle{(11')\quad \begin{array}{lll} V(\bar{x})= \underbrace{x_1^T(t)P_1x_1(t)}_{V(x_1)}+\underbrace{(s(t)-S_rr(t))^TP_2(s(t)-S_rr(t))}_{V(s)}\\ \Rightarrow \dot{V}(\bar{x})\le -x_1^T(t)Q_1x_1(t)-(s(t)-S_rr(t))^T(s(t)-S_rr(t)) \end{array}} }$

が成り立つように決定します（ $P_1>0$ , $P_2>0$ , $Q_1>0$ ）。

[2] 可到達性の検討

ここでは、スライディングモード制御則(10)の具体的な表現を求めます。

●等価制御は、(8)において $r(t)=0$ の場合

$\displaystyle{ \begin{array}{cl} (12.1) & s(t)=0\\ \Downarrow &\\ (12.2) & \dot{s}(t)=0\\ \Downarrow &\\ (12.3) & 0=S_2\bar{A}_{21}x_1(t)+S_2\bar{A}_{22}S_2^{-1}s(t)+S_2B_2u(t)\\ \Downarrow &\\ (12.4) & u_{eq}(t)=-\underbrace{(S_2B_2)^{-1}}_{([0\ I]T_sB_E)^{-1}} \underbrace{\left[\begin{array}{cc} S_2\bar{A}_{21} & S_2\bar{A}_{22}S_2^{-1} \\ \end{array}\right]}_{[0\ I]T_sA_ET_s^{-1}} \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]\\ \Downarrow & by\ S=[0\ I]T_s\\ (12.5) & u_{eq}(t)=-(SB_E)^{-1}SA_ET_s^{-1} \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]\\ \Downarrow & by\ (7),(8b)\\ (12.6) & u_{eq}(t)=-(SB_E)^{-1}SA_E \left[\begin{array}{c} x_r(t)\\ x(t) \end{array}\right]} \end{array} }$

のように得られます。(10)の第１項 $u_\ell$ は、この等価制御をベースして、(9)を考慮して

$\displaystyle{(13)\quad \begin{array}{l} u_\ell(t)=-(S_2B_2)^{-1} (\left[\begin{array}{cc} S_2\bar{A}_{21} & S_2\bar{A}_{22}S_2^{-1} \\ \end{array}\right]\left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]-\Phi s(t))\\ -(S_2B_2)^{-1}(\Phi S_r+S_1B_r)r(t)+(S_2B_2)^{-1}S_r\dot{r}(t)\\ \end{array} }$

すなわち

$\displaystyle{(13')\quad { \boxed{\begin{array}{l} u_\ell(t) =-\underbrace{(SB_E)^{-1}(SA_E-\Phi S)}_{L=L_{eq}+L_\Phi}\left[\begin{array}{c} x_r(t)\\ x(t) \end{array}\right]\\ -\underbrace{(SB_E)^{-1}(\Phi S_r+S_1B_r)}_{L_r} r(t) +\underbrace{(SB_E)^{-1}S_r}_{L_{\dot r}} \dot{r}(t)\\ \end{array}}} }$

のように構成します（ $\Phi$ は安定行列）。このとき閉ループ系は次式で与えられます。

$\displaystyle{(14)\quad \begin{array}{l} \left[\begin{array}{c} \dot x_1(t)\\ \dot s(t) \end{array}\right] = \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ S_2\bar{A}_{21} & S_2\bar{A}_{22}S_2^{-1} \\ \end{array}\right] \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]\\ + \left[\begin{array}{cc} 0\\ S_2B_2 \end{array}\right] (-(S_2B_2)^{-1} (\left[\begin{array}{cc} S_2\bar{A}_{21} & S_2\bar{A}_{22}S_2^{-1} \\ \end{array}\right]\left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]-\Phi s(t))\\ -(S_2B_2)^{-1}(\Phi S_r+S_1B_r)r(t)+(S_2B_2)^{-1}S_r\dot{r}(t)+u_n(t)) + \left[\begin{array}{cc} B_r \\ S_1B_r \end{array}\right] r(t)\\ = \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ S_2\bar{A}_{21} & S_2\bar{A}_{22}S_2^{-1} \\ \end{array}\right] \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]\\ + \left[\begin{array}{cc} 0\\ I_m \end{array}\right] \left[\begin{array}{cc} -S_2\bar{A}_{21} & -S_2\bar{A}_{22}S_2^{-1}+\Phi \\ \end{array}\right] \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]\\ + \left[\begin{array}{cc} 0\\ -(\Phi S_r+S_1B_r)r(t)+S_r\dot{r}(t)+S_2B_2u_n(t) \end{array}\right] + \left[\begin{array}{cc} B_r \\ S_1B_r \end{array}\right] r(t)\\ = \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ 0 & \Phi \\ \end{array}\right] \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]+ \left[\begin{array}{cc} B_rr(t)\\ -\Phi S_rr(t)+S_r\dot{r}(t)+S_2B_2u_n(t) \end{array}\right] \end{array} }$

これを変形して

$\displaystyle{(15)\quad \begin{array}{l} \left[\begin{array}{c} \dot{x}_1(t)\\ \dot{s}(t)-S_r\dot{r}(t) \end{array}\right] = \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ 0 & \Phi \\ \end{array}\right] \left[\begin{array}{c} x_1(t)\\ s(t)-S_rr(t) \end{array}\right]\\ + \left[\begin{array}{cc} B_rr(t)\\ -\Phi S_rr(t)+S_r\dot{r}(t)+S_2B_2u_n(t) \end{array}\right] \end{array} }$

すなわち、次式が成り立ちます。

$\displaystyle{ \begin{array}{cl} (15.1) & \dot{x}_1(t)=\bar{A}_{11}x_1(t)+\bar{A}_{12}(s(t)-S_rr(t))+\underbrace{(B_r+\bar{A}_{12}S_r)}_{B_r'}r(t)\\ (15.2) & \dot{s}(t)-S_r\dot{r}(t)=\Phi(s(t)-S_rr(t))+S_2B_2u_n(t) \end{array} }$

(15.2)に基づいてスライディングモード $s-S_rr=0$ を達成する $u_n$ を明らかにします。

$\Phi$ は安定行列なので

$\displaystyle{(16)\quad \begin{array}{l} P_2\Phi+\Phi^TP_2=-I_m \end{array} }$

を満たす $P_2>0$ を選ぶことができます。これを用いて

$\displaystyle{(17)\quad \begin{array}{l} u_n(t)=-(S_2B_2)^{-1}\rho(u_\ell,y)\frac{P_2(s(t)-S_rr(t))}{||P_2(s(t)-S_rr(t))||} \end{array} }$

すなわち

$\displaystyle{(17')\quad { \boxed{ \begin{array}{l} u_n(t) =-\underbrace{(SB_E)^{-1}\rho(u_\ell,y)}_{L_n}\frac{P_2(s(t)-S_rr(t))}{||P_2(s(t)-S_rr(t))||} \end{array} }} }$

と選びます（ $\rho>0$ ）。このとき次式が成り立ち、(11)が示されます。

$\displaystyle{(18)\quad \begin{array}{l} \dot{V}(s)=2(s(t)-S_rr(t))^TP_2(\dot{s}(t)-S_r\dot{r}(t))\\ =2(s(t)-S_rr(t))^TP_2 (\Phi (s(t)-S_rr(t))-\rho(u_\ell,y)\frac{P_2(s(t)-S_rr(t))}{||P_2(s(t)-S_rr(t))||})\\ =(s(t)-S_rr(t))^T\underbrace{(P_2\Phi+\Phi^TP_2)}_{-I_m}(s(t)-S_rr(t))\\ +2(s(t)-S_rr(t))^TP_2(-\rho(u_\ell,y)\frac{P_2(s(t)-S_rr(t))}{||P_2(s(t)-S_rr(t))||})\\ \le -||s(t)-S_rr(t)||^2-2\rho(u_\ell,y)||P_2(s(t)-S_rr(t))||\\ \le -||s(t)-S_rr(t)||^2 \end{array} }$

これより(15.2)における２次安定性が成り立つことが分かります。２次安定であれば $s-S_rr$ は0に漸近します。

[3] スライディングモードの検討

次に $s-S_rr=0$ のもとで、(15.1)から $x_1$ の振舞いについて調べます。

$\bar{A}_{11}$ は安定行列なので

$\displaystyle{(19)\quad \begin{array}{l} P_1\bar{A}_{11}+\bar{A}_{11}^TP_1=-Q_1<0 \end{array} }$

を満たす $P_1>0$ を選ぶことができます。そこで、(15.1)に対するリャプノフ関数 $V(x_1)=x_1^T(t)P_1x_1(t)$ の微分について調べます。

$\displaystyle{(20)\quad \begin{array}{l} \dot{V}(x_1)=2x_1^T(t)P_1\dot{x}_1(t)\\ =2x_1^T(t)P_1( \bar{A}_{11}x_1(t)+\bar{A}_{12}\underbrace{(s(t)-S_rr(t))}_{0}+B_r'r(t))\\ =x_1^T(t)\underbrace{(P_1\bar{A}_{11}+\bar{A}_{11}^TP_1)}_{-Q_1}x_1(t)+2x_1^T(t)P_1B_r'r(t)\\ =-x_1^T(t)Q_1x_1(t)+2x_1^T(t)P_1B_r'r(t)\\ \le -\sigma_n(Q_1)x_1^T(t)x_1(t)+2||x_1^T(t)P_1||\,||B_r'r(t)||\\ \le -\sigma_n(Q_1)||x_1(t)||^2+2\sqrt{x_1^T(t)P_1^2x_1(t)}\,||B_r'r(t)||\\ \le -\sigma_n(Q_1)||x_1(t)||^2+2\sqrt{\sigma_1(P_1^2)x_1^T(t)x_1(t)}\,||B_r'r(t)||\\ = -\sigma_n(Q_1)||x_1(t)||^2+2\sigma_1(P_1)||x_1(t)||\,||B_r'r(t)||\\ = -\sigma_1(P_1)||x_1(t)||(\frac{\sigma_n(Q_1)}{\sigma_1(P_1)}||x_1(t)||-2||B_r'r(t)||) \end{array} }$

したがって、次式が成り立ち、リャプノフ安定性が示されます。

$\displaystyle{(21)\quad||B_r'r(t)||<\frac{1}{2}\frac{\sigma_n(Q_1)}{\sigma_1(P_1)}||x_1(t)|| \Rightarrow \dot{V}(x_1)<0 }$

これから、ロバスト性での議論と同様に、 $x_1$ が発散することはないことが分かります。

演習C61…Flipped Classroom

MATLAB

%cCIP_smci.m
%-----
 clear all, close all
 global mc m ell g th0
 mc=1; m=0.1; ell=0.2; g=9.8;
 ths=0; %input('ths   = <0,180> ')/180*pi;
 th0=3/180*pi; %input('th(0) = <0,180> ')/180*pi;
%-----
 A=[zeros(2,2) eye(2);zeros(2,4)];
 A(3,1)=0; 
 A(3,2)=-3*m*g/(m+4*mc); 
 A(4,1)=0; 
 A(4,2)=3*(m+mc)*g/((m+4*mc)*ell);  
 B=zeros(4,1);
 B(3)=4/(m+4*mc);
 B(4)=-3/((m+4*mc)*ell); 
 A0=A; B0=B; 
 H=[1 0 0 0];
 A=[0 -H ;
    zeros(4,1) A0];
 B=[0; B0];
 CM=[1 0 0 0;
     0 1 0 0];
 Gamma=-5; 
%-----
 lambda=[-0.8 0.5 -1.5 1];
 nocomp=2;
 specpos=rand(4,4);
 specent=rand(4,4);
 S=swfvpl(A,B,lambda,nocomp,specpos,specent)      
 SS=S;
 F=(SS*B)\(SS*A)
 eig(A-B*F)
%-----
[nn,mm]=size(B);
 S1=SS(:,1:nn-mm);
 S2=SS(:,nn-mm+1:nn);
 Br=eye(nn-mm,mm);
 B2=B(nn-mm+1:nn,:); 
 Lambda=S2*B2;
 Phi=-0.1;
 P2=0.5*inv(-Phi);
 Check=P2*Phi+Phi*P2
 P2S=P2*S;
 Leq=inv(S*B)*S*A;
 LPhi=-inv(S*B)*Phi*S;
 L=Leq+LPhi;
 rho=2;
 Ln=inv(S*B)*rho;
%-----
 Sr=0;
 Lr=inv(Lambda)*(Phi*Sr+S1*Br)
 Ldr=inv(Lambda)*Sr
%Ln=0.1*inv(Lambda)
%-----
 x0=[0;th0*0;0;0];
 sim('CIP_smci_2015a.mdl')
%-----
%eof

SCILAB

応用例

投稿日時: 2021年10月29日投稿者: cacsd

応用例…Homework

[1]　ある航空機の線形状態方程式として、次を考えます。

$A$ 行列の固有値は次のように求められます。

状態変数として、pitch angle $\theta$ の代わりに、fligt path angle $\gamma=\theta-\alpha$ を用いると、次式となります。

また、評価変数として、pitch angle $\theta$ とfligt path angle $\gamma$ をとると、その出力方程式は次式のように与えられます。

$\displaystyle{(4)\quad \begin{array}{l} y(t)=Hx(t)\\ H= \left[\begin{array}{ccccc} 1 & 0 & 1 & 0 & 0\\ 1 & 0 & 0 & 0 & 0 \end{array}\right] \end{array} }$

制御目的は、次図のように、pitch angle $\theta$ とfligt path angle $\gamma$ の間の非干渉化を達成し、独立して設定値に合せることです。

図1 航空機の操縦法

[2] モデル規範による追従制御

この場合の制御則は次のように構成されます。

$\displaystyle{(5)\quad { \boxed{\begin{array}{l} u(t)=-Fx(t)+Gr(t)+\underbrace{u_\ell(t)}_{linear\ control}+\underbrace{u_n(t)}_{switching\ component}\\ u_\ell(t)=-\underbrace{(SB)^{-1}(SA_m-\Phi S)}_{L=L_{eq}+L_{phi}}x(t) \quad(\Phi:stable\ matrix)\\ u_n(t)=-\underbrace{(SB)^{-1}\rho(t,x)}_{L_n}\frac{P_2s(t)}{||P_2s(t)||}\quad(P_2\Phi+\Phi^TP_2=-I) \end{array}}} }$

●追従すべき規範モデルは

$\displaystyle{(6)\quad \begin{array}{l} \dot{x}(t)=A_mx(t)+B_mu(t)\\ y(t)=Hx(t) \end{array}} }$

ただし

$\displaystyle{(7)\quad \begin{array}{l} A_m=A-BF\\ B_m=BG\\ G=-(H(A-BF)^{-1}B)^{-1} \end{array} }$

いまこの規範モデルにもたせるべき望ましい固有値を

$\displaystyle{(8)\quad \lambda(A_m):\left\{\begin{array}{cl} -5.6\pm j 4.2 & short\ period\ mode\\ -1 & pitch\ attitude\ mode\\ -20 & elevator\ actuator\ mode\\ -20 & flap\ actuator\ mode \end{array}\right. }$

とし、また望ましいモード分布を次により指定します。

$\displaystyle{(9)\quad {\tt specpos}= \left[\begin{array}{ccc|cc} 1 & 1 & 1 & 0 & 0\\ 1 & 0 & 1 & 0 & 0\\ 0 & 1 & 1 & 0 & 0\\\hline 0 & 0 & 0 & 1 & 0\\ 0 & 0 & 0 & 0 & 1\\ \end{array}\right],\quad {\tt specent}= \left[\begin{array}{rrr|rr} 0 & 0 & 1 & -1 & -1\\ 1 & -1 & 0 & -1 & -1\\ -1 & 1 & 0 & -1 & -1\\\hline -1 & -1 & -1 & 1 & -1\\ -1 & -1 & -1 & -1 & 1\\ \end{array}\right] }$

ここで、 ${\tt specpos}$ は望ましい固有ベクトルの要素を指定する場合１、指定しない場合０で表した行列、 ${\tt specエent}$ は望ましい固有ベクトルの要素が非零指定の場合１、零の場合０、非零任意の場合－１で表した行列です。

このとき以下に示すプログラムを使って、次が得られます。

$\displaystyle{(10)\quad \begin{array}{l} F=WV^{-1} =\left[\begin{array}{rrrrr} -1.9338 & -0.5744 & -2.1188 & 0.4750 & 0.1066\\ 1.9571 & 0.2253 & 3.1273 & -0.0694 & -0.0515 \end{array}\right]\\ W=\left[\begin{array}{rrrrr} 0.8445 & 1.9935 & -0.2742 & 0 & 0\\ -2.7055 & 1.2780 & 0.1767 & 0 & 0 \end{array}\right]\\ V=\left[\begin{array}{rrrrr} 0 & -0.0000 & 0.6667 & 0.0030 & -0.0313\\ 1.0000 & -9.5000 & -0.3333 & 0.9964 & 0.1493\\ -0.9286 & 1.0000 & -0.3333 & -0.0528 & 0.0238\\ -1.8252 & -2.2364 & 0.2886 & 1.0000 & -0.0649\\ 2.9860 & -2.6459 & -0.1860 & -0.0825 & 1.0000 \end{array}\right] \end{array} }$

●ここでは、スイッチング関数

$\displaystyle{(11)\quad s(t)= \underbrace{ \left[\begin{array}{cc} M & I \end{array}\right] }_{S} x(t) }$

を構成するために

$\displaystyle{(12)\quad \begin{array}{l} \dot{x}_1(t)=A_{11}x_1(t)+A_{12}\hat{u}(t) \end{array} }$

を２次形式評価関数

$\displaystyle{(13)\quad J=\int_0^\infty(x_1^T(t)Q_{11}x_1(t)+\hat{u}^T(t)Q_{22}\hat{u}(t))dt }$

を最小にする

$\displaystyle{(14)\quad \hat{u}(t)=-Mx_1(t) }$

を求めます。

$\displaystyle{(15)\quad Q_{11}={\rm diag}\{10,5,5\},\ Q_{22}={\rm diag}\{20,20\} }$

と選ぶとき、スイッチング関数を決める行列 $S$ が、次のように得られます。

$\displaystyle{(16)\quad S=\left[\begin{array}{rrr|rr} 0.7025 & 0.4209 & 0.1894 & -1 & 0\\ -0.0810 & 0.0635 & 0.0267 & 0 & -1\\ \end{array}\right] }$

● $P_2\Phi+\Phi^TP_2=-I$ を満たす $\Phi$ と $P_2$ は、たとえば次のように定めます。

$\displaystyle{(17)\quad \underbrace{ \left[\begin{array}{rr} 0.025 & 0\\ 0 & 0.025 \end{array}\right] }_{P_2} \underbrace{ \left[\begin{array}{rr} 20 & 0\\ 0 & 20 \end{array}\right] }_{\Phi} +(*)^T=-I }$

このときSMC則における行列 $L$ は次のように計算されます。

$\displaystyle{(18)\quad L=\left[\begin{array}{rrrrr} -1.2285 & -0.1858 & -2.1624 & 0.0782 & 0.0460\\ 1.8631 & 0.2827 & 3.0804 & -0.1299 & -0.0060 \end{array}\right] }$

図2 モデル規範法によるSMC制御系応答

$\displaystyle{{ \boxed{\begin{array}{l} u(t)=-Fx(t)+Gr(t)+\underbrace{u_\ell(t)}_{linear\ control}+\underbrace{u_n(t)}_{switching\ component}\\ u_\ell(t)=-\underbrace{(SB)^{-1}(SA_m-\Phi S)}_{L=L_{eq}+L_{phi}}x(t) \quad(\Phi:stable\ matrix)\\ u_n(t)=-\underbrace{(SB)^{-1}\rho(t,x)}_{L_n}\frac{P_2s(t)}{||P_2s(t)||}\quad(P_2\Phi+\Phi^TP_2=-I) \end{array}}} }$

MATLAB

%cAIRCRAFT_smcm.m
%-----
 clear all, close all
 A= [0 0 1.74 0.08 0.59;
     0 -1.99 -13.41 -18.95 -3.6;
     0 1 -1.74 -0.08 -0.59;
     0 0 0 -20 0;
     0 0 0 0 -20];
 B=[0 0;0 0;0 0;20 0;0 20];
%----- 
 lambda=[-5.6 4.2 -1 -20 -20];
 nocomp=1;
%固有ベクトルの非零要素の位置 (１：指定、０：指定せず)
 specpos=[ 1  1  1  0  0; 
           1  0  1  0  0; 
           0  1  1  0  0;
           0  0  0  1  0;
           0  0  0  0  1]
%固有ベクトルの要素の値 (０：零、１：非零指定、－１：非零任意)       
 specent=[ 0  0  1 -1 -1;
           1 -1  0 -1 -1;
          -1  1  0 -1 -1;
          -1 -1 -1  1 -1;
          -1 -1 -1 -1  1]
 F=vplace2(A,B,lambda,nocomp,specpos,specent)       
 pl=eig(A-B*F)
%-----
 H=[1 0 1 0 0;1 0 0 0 0];
 Am=A-B*F;
 G=-inv(H*(Am\B));
 Bm=B*G;
%-----
 Q=diag([10,5,5,20,20]);
 S=swflqr(A,B,Q)
%-----
 P2=diag([0.025,0.025]);
 Phi=-20/0.025*P2;
 Check=P2*Phi+Phi*P2
 L=inv(S*B)*(S*Am-Phi*S)
 Ln=inv(S*B)
%-----
 sim('AIRCRAFT_smcm_2015a.mdl')
%-----
%eof

モデル規範による追従SM制御

投稿日時: 2021年10月29日投稿者: cacsd

モデル規範による追従SM制御…Homework

[1]　次の状態方程式で表される制御対象を考えます。

$\displaystyle{(1)\quad \dot{x}(t)=Ax(t)+Bu(t) }$

この状態が理想的なモデル

$\displaystyle{(2)\quad \dot{x}_m(t)=A_mx_m(t)+B_mr(t) }$

の状態を追従するように、すなわち

$\displaystyle{(3)\quad e(t)=x(t)-x_m(t)\rightarrow 0 \quad(t\rightarrow\infty) }$

となるように制御則を決定したいとします。

●追従誤差のダイナミックスは

$\displaystyle{(4)\quad \begin{array}{l} %\bullet \dot{e}(t)=Ax(t)+Bu(t)-(A_mx_m(t)+B_mr(t))\\ %=Ax(t)-Ax_m(t)+Ax_m(t)-A_mx_m(t)+Bu(t)-B_mr(t)\\ %=Ae(t)+(A-A_m)x_m(t)+Bu(t)-B_mr(t)\\ %\bullet \dot{e}(t)=Ax(t)+Bu(t)-(A_mx_m(t)+B_mr(t))\\ =Ax(t)-A_mx(t)+A_mx(t)-A_mx_m(t)+Bu(t)-B_mr(t)\\ =A_me(t)+(A-A_m)x(t)+Bu(t)-B_mr(t) \end{array} }$

となります。いま条件

$\displaystyle{(5)\quad \begin{array}{l} (BB^\dag-I)(A-A_m)=0\Leftrightarrow B\underbrace{B^\dag(A-A_m)}_{F}=A-A_m\\ (BB^\dag-I)B_m=0\Leftrightarrow B\underbrace{B^\dag B_m}_{G}=B_m \end{array} }$

を仮定すると（ $B^\dag$ は $BB^\dag B=B$ を満たす疑似逆行列）、制御則

$\displaystyle{(6)\quad u(t)=-Ke(t)-Fx(t)+Gr(t) }$

の下では、誤差方程式は次式となります。

$\displaystyle{(7)\quad \dot{e}(t)=(A_m-BK)e(t) }$

実際、

$\displaystyle{(8)\quad \begin{array}{l} \dot{e}(t)=A_me(t)+(A-A_m)x(t)+Bu-B_mr(t)\\ =A_me(t)+(A-A_m)x(t)+B(-Ke(t)-Fx(t)+Gr(t))-B_mr(t)\\ =A_me(t)+(A-A_m)x(t)-BKe(t)-BB^\dag(A-A_m)x(t)\\ +BB^\dag B_mr(t)-B_mr(t)\\ =(A_m-BK)e(t)-\underbrace{(BB^\dag-I)(A-A_m)}_{0}x(t)+\underbrace{(BB^\dag -I)B_m}_{0}r(t)\\ \end{array} }$

●以上から、制御則(6)の設計の指針として、望ましいモデル(2)を(1)に対する状態フィードバックと入力変換で得るものとします。すなわち

$\displaystyle{(9)\quad \boxed{\begin{array}{l} A_m=A-BF\\ B_m=BG \end{array}} }$

ここで $G$ は、 $r$ から出力

$\displaystyle{(10)\quad y(t)=Hx(t) }$

までの定常ゲイン

$\displaystyle{(11)\quad %\begin{array}{l} \left\{\begin{array}{l} 0=(A-BF)x(\infty)+BGr(\infty)\\ y(\infty)=Hx(\infty) \end{array}\right. \Rightarrow y(\infty)=\underbrace{-H(A-BF)^{-1}BG}_{I}r(\infty) %\end{array} }$

が単位行列となるように、次式で決めます。

$\displaystyle{(12)\quad \boxed{G=-(H(A-BF)^{-1}B)^{-1}} }$

ちなみに $K$ については、(7)が漸近安定であることが前提となり、望ましいモデル(2)より、追従誤差は速いダイナミックスをもつことが必要であることに注意します。

[2]　誤差方程式(4)すなわち

$\displaystyle{(13)\quad \dot{e}(t)=A_me(t)+(A-A_m)x(t)+Bu(t)-B_mr(t) }$

に対して、スライディングモード制御を適用することを考えます。スイッチング関数を

$\displaystyle{(14)\quad s(t)=Se(t) }$

とするとき、スライディングモード時は

$\displaystyle{(15)\quad Se(t)=0\Rightarrow S\dot{e}(t)=0 }$

すなわち

$\displaystyle{(16)\quad S\dot{e}(t)=S(A_me(t)+(A-A_m)x(t)+Bu(t)-B_mr(t))=0 }$

が成り立ち、等価制御は次式となります。

$\displaystyle{(17)\quad u_{eq}(t)=-(SB)^{-1}S(A_me(t)+(A-A_m)x(t)-B_mr(t)) }$

このとき追従誤差のダイナミックスは次式となります。

$\displaystyle{(18)\quad \begin{array}{l} \dot{e}(t)=A_me(t)+(A-A_m)x(t)+Bu(t)-B_mr(t)\\ =A_me(t)+(A-A_m)x(t)\\ +B(-(SB)^{-1}S(A_me(t)+(A-A_m)x(t)-B_mr(t)))-B_mr(t)\\ =(I-B(SB)^{-1}S)(A_me(t)+(A-A_m)x(t)-B_mr(t))\\ =(I-B(SB)^{-1}S)(A_me(t)+BFx(t)-BGr(t))\\ =(I-B(SB)^{-1}S)A_me(t)+\underbrace{(I-B(SB)^{-1}S)B}_{0}(Fx(t)-Gr(t))\\ =\underbrace{(I-B(SB)^{-1}S)A_m}_{A_{eq}}e(t) \end{array} }$

したがって、スライディングモード制御を適用する場合、全体の制御則は、次式で表されます。

$\displaystyle{(19)\quad { \boxed{\begin{array}{l} u(t)=-Fx(t)+Gr(t)+\underbrace{u_\ell(t)}_{linear\ control}+\underbrace{u_n(t)}_{switching\ component}\\ u_\ell(t)=-\underbrace{(SB)^{-1}(SA_m-\Phi S)}_{L=L_{eq}+L_{phi}}x(t) \quad(\Phi:stable\ matrix)\\ \displaystyle{u_n(t)=-\underbrace{(SB)^{-1}\rho(t,x)}_{L_n}\frac{P_2s(t)}{||P_2s(t)||}\quad(P_2\Phi+\Phi^TP_2=-I)} \end{array}}} }$

演習C51…Flipped Classroom

MATLAB

%cCIP_smcm.m
%-----
 clear all, close all
 global mc m ell g th0
 mc=1; m=0.1; ell=0.2; g=9.8;
 ths=0; %input('ths   = <0,180> ')/180*pi;
 th0=3/180*pi; %input('th(0) = <0,180> ')/180*pi;
%-----
 A=[zeros(2,2) eye(2);zeros(2,4)];
 A(3,1)=0; 
 A(3,2)=-3*m*g/(m+4*mc); 
 A(4,1)=0; 
 A(4,2)=3*(m+mc)*g/((m+4*mc)*ell);  
 B=zeros(4,1);
 B(3)=4/(m+4*mc);
 B(4)=-3/((m+4*mc)*ell); 
 H=[1 0 0 0];
%-----
 lambda=[-0.8*1.5 0.5 -1.5*1.5 1];
 nocomp=2;
 specpos=rand(4,4);
 specent=rand(4,4);
 F=vplace2(A,B,lambda,nocomp,specpos,specent)       
 pl=eig(A-B*F)
%-----
 H=[1 0 0 0];
 Am=A-B*F;
 G=-inv(H*(Am\B));
 Bm=B*G;
%-----
 Mr=0.5; Mth=3/180*pi; 
 Tc=1; Tpen=0.25*2*pi*sqrt(2*ell/g);
 Q=diag([1/Mr^2,1/Mth^2,Tc^2/Mr^2,Tpen^2/Mth^2]);
 S=swflqr(A,B,Q) 
%-----
 Phi=-0.1;
 P2=0.5*inv(-Phi);
 Check=P2*Phi+Phi*P2
 P2S=P2*S;
 Leq=inv(S*B)*S*A;
 LPhi=-inv(S*B)*Phi*S;
 L=Leq+LPhi;
 rho=1;
 Ln=inv(S*B)*rho;
%-----
 x0=[0;th0*0;0;0];
 sim('CIP_smcm_2015a.mdl')
%-----
%eof

固有ベクトル設定問題

投稿日時: 2021年10月29日投稿者: cacsd

[0] 次の２入力２次系

$\displaystyle{(1)\quad \dot{x}= \underbrace{ \left[\begin{array}{cc} 0 & 0\\ 0 & -1 \end{array}\right] }_{A} x+ \underbrace{ \left[\begin{array}{cc} 1 & 1\\ 1 & -1 \end{array}\right] }_{B} u }$

に対して２つの異なる状態フィードバック

$\displaystyle{(2)\quad %\begin{array}{l} u=- \underbrace{ \left[\begin{array}{cc} 1 & 1\\ 1 & -1 \end{array}\right] }_{F_1} x,\quad u=- \underbrace{ \left[\begin{array}{cc} 1 & 0\\ 1 & 2 \end{array}\right] }_{F_2} x %\end{array} }$

による閉ループ系は、それぞれ

$\displaystyle{(3)\quad \dot{x}= \underbrace{ \left[\begin{array}{cc} -2 & 0\\ 0 & -3 \end{array}\right] }_{A-BF_1} x,\quad \dot{x}= \underbrace{ \left[\begin{array}{cc} -2 & 2\\ 0 & -3 \end{array}\right] }_{A-BF_2} x }$

となって、同じ固有値-2,-3を持ちます。このように多入力系では閉ループ系の固有値を指定しても状態フィードバックは一意に定まりません。これは固有ベクトルの指定に任意性があるからです（ちなみに１入力系では一意に定まります）。そこで、上の２つの閉ループ系の固有ベクトルを求めてみると

$\displaystyle{(4)\quad V_1=\left[\begin{array}{cc} 1 & 0\\ 0 & 1 \end{array}\right],\quad V_2=\left[\begin{array}{cc} 1 & 1\\ 0 & 1 \end{array}\right] }$

また応答は次式となります。

$\displaystyle{(5)\quad \dot{x}= \underbrace{ \left[\begin{array}{cc} e^{-2t} & 0\\ 0 & e^{-3t} \end{array}\right] }_{\exp((A-BF_1)t)=V_1\exp(\Lambda t)V_1^{-1}} x(0) }$
$\displaystyle{(6)\quad \dot{x}= \underbrace{ \left[\begin{array}{cc} 1 & 1\\ 0 & 1 \end{array}\right] \left[\begin{array}{cc} e^{-2t} & 0\\ 0 & e^{-3t} \end{array}\right] \left[\begin{array}{cc} 1 & -1\\ 0 & 1 \end{array}\right] }_{\exp((A-BF_2)t)=V_2\exp(\Lambda t)V_2^{-1}} x(0) }$

明らかに(5)ではモード $e^{-2t}$ と $e^{-3t}$ が分離していますが、(6)ではそれらがミックスされています。このように、多入力系の場合の状態フィードバックの設計は、モードをどのように分布させるかがポイントとなります。

[1] 実固有値の場合

$\displaystyle{(7)\quad (A-BF)v=\lambda v \Rightarrow \underbrace{ \left[\begin{array}{cc} \lambda I-A & B \end{array}\right] }_{G(\lambda)} \left[\begin{array}{c} v\\ \xi \end{array}\right] =0 \quad(\xi=Fv) }$

より、 $v$ と $\xi$ は、サイズ $n\times(n+m)$ の行列 $G(\lambda)$ の零化空間の基底行列に対して、適当な $m$ 次元ベクトル $\delta$ を用いて

$\displaystyle{(8)\quad G(\lambda) \left[\begin{array}{c} N(\lambda)\\ M(\lambda)\\ \end{array}\right] =0 \Rightarrow \left[\begin{array}{c} v\\ \xi \end{array}\right]= \left[\begin{array}{c} N(\lambda)\\ M(\lambda)\\ \end{array}\right]\delta }$

のように表されます。そこで望ましい固有ベクトル $v$ を与えて

$\displaystyle{(9)\quad N(\lambda)\delta=v }$

を解いて $\delta$ を求めます。ただし、望ましい固有ベクトル $v$ のすべての成分を指定できるわけではないので

$\displaystyle{(10)\quad \tilde{R}N(\lambda)\delta=\tilde{R}v =\left[\begin{array}{c} \tilde{v}^s\\ \tilde{v}^u\\ \end{array}\right] \Rightarrow \left[\begin{array}{cc} I_q & 0 \end{array}\right] \tilde{R}N(\lambda)\times\delta=\tilde{v}^s }$

の最小２乗解として $\delta=\delta^*$ を定めます。ここで、セレクタ行列 $\tilde{R}$ により、望ましい固有ベクトル $v$ の指定された $q$ 個の成分だけを集めた $\tilde{v}^s$ が得られるとしています。

このとき、割当可能な $v^*$ と $\xi^*$ は次式で与えられます。

$\displaystyle{(11)\quad \left[\begin{array}{c} v^*\\ \xi^* \end{array}\right]= \left[\begin{array}{c} N(\lambda)\\ M(\lambda)\\ \end{array}\right]\delta^* }$

●もし $A-\lambda I_n$ が正則ならば

$\displaystyle{(12)\quad (\lambda I_n-A)\underbrace{(\lambda I_n-A)^{-1}B}_{N(\lambda)}+B\underbrace{(-I_m)}_{M(\lambda)}=0 }$

より、 $N(\lambda)$ と $M(\lambda)$ が特定できます。(8)より

$\displaystyle{(13)\quad \begin{array}{l} v=(\lambda I_n-A)^{-1}B\delta\\ \xi=-\delta \end{array} \quad\Rightarrow\quad \boxed{v=(A-\lambda I_n)^{-1}B\xi} }$

となって、指定した固有ベクトル $v$ を近似する $\xi$ を直接求めることができます。

●具体的に、(1)に対して、固有値-2,-3を持たせる状態フィードバックを、対応する固有ベクトルが(4)の $V_1$ となるように求めてみます。

$\displaystyle{(14)\quad \begin{array}{l} \underbrace{ \left[\begin{array}{cc} 0-(-2) & 0\\ 0 & -1-(-2) \end{array}\right]^{-1} \left[\begin{array}{cc} 1 & 1\\ 1 & -1 \end{array}\right] }_{(A-\lambda I_n)^{-1}B} \xi_1 =\underbrace{ \left[\begin{array}{c} 1\\ 0 \end{array}\right]}_{V_1(:,1)}\\ \Rightarrow \xi_1= \left[\begin{array}{cc} 1/2 & 1/2\\ 1 & -1 \end{array}\right]^{-1} \left[\begin{array}{c} 1\\ 0 \end{array}\right] = \left[\begin{array}{cc} 1 & 1/2\\ 1 & -1/2 \end{array}\right] \left[\begin{array}{c} 1\\ 0 \end{array}\right] = \left[\begin{array}{cc} 1\\ 1 \end{array}\right] \end{array} }$

$\displaystyle{(15)\quad \begin{array}{l} \underbrace{ \left[\begin{array}{cc} 0-(-3) & 0\\ 0 & -1-(-3) \end{array}\right]^{-1} \left[\begin{array}{cc} 1 & 1\\ 1 & -1 \end{array}\right] }_{(A-\lambda I_n)^{-1}B} \xi_2 =\underbrace{ \left[\begin{array}{c} 0\\ 1 \end{array}\right]}_{V_1(:,2)}\\ \Rightarrow \xi_2= \left[\begin{array}{cc} 1/3 & 1/3\\ 1/2 & -1/2 \end{array}\right]^{-1} \left[\begin{array}{c} 0\\ 1 \end{array}\right] = \left[\begin{array}{cc} 3/2 & 1\\ 3/2 & -1 \end{array}\right] \left[\begin{array}{c} 0\\ 1 \end{array}\right] = \left[\begin{array}{cc} 1\\ -1 \end{array}\right] \end{array} }$

状態フィードバックゲイン $F$ は

$\displaystyle{(16)\quad \left[\begin{array}{cc} \xi_1 & \xi_2 \end{array}\right] =FV_1 \quad\Rightarrow\quad F= \left[\begin{array}{cc} \xi_1 & \xi_2 \end{array}\right]V_1^{-1} = \left[\begin{array}{cc} 1 & 1\\ 1 & -1 \end{array}\right] }$

[2] 複素共役固有値の場合

$\displaystyle{(17)\quad \begin{array}{l} (A-BF)(v_R\pm jv_I)=(\lambda_R\pm j\lambda_I)(v_R\pm jv_I)\\ \Rightarrow (A-BF)v_R=\lambda_Rv_R-\lambda_Iv_I,\ (A-BF)v_I =\lambda_Rv_I+\lambda_Iv_R \\ \Rightarrow \underbrace{ \left[\begin{array}{ccc} \lambda_R I-A & \lambda_I I & B \end{array}\right] }_{G(\lambda)} \left[\begin{array}{cc} v_R & v_I\\ -v_I & v_R\\ \xi_R & \xi_I \end{array}\right] =0 \quad(\xi_R=Fv_R,\xi_I=Fv_I) \end{array} }$

より、 $v_R,v_I$ と $\xi_R,\xi_I$ は、サイズ $n\times(2n+m)$ の行列 $G(\lambda)$ の零化空間の基底行列に対して、適当な $m$ 次元ベクトル $\delta_R,\delta_I$ を用いて

$\displaystyle{(18)\quad G(\lambda) \left[\begin{array}{c} N(\lambda)\\ P(\lambda)\\ M(\lambda) \end{array}\right] =0 \Rightarrow \left[\begin{array}{cc} v_R & v_I\\ -v_I & v_R\\ \xi_R & \xi_I \end{array}\right]= \left[\begin{array}{c} N(\lambda)\\ P(\lambda)\\ M(\lambda) \end{array}\right] \left[\begin{array}{cc} \delta_R & \delta_I \end{array}\right] }$

のように表されます。そこで望ましい $v_R,v_I$ を与えて

$\displaystyle{(19)\quad \left[\begin{array}{c} v_R \\ v_I \end{array}\right]= \underbrace{ \left[\begin{array}{cc} N(\lambda) \\ -P(\lambda) \end{array}\right] }_{\alpha} \delta_R= \underbrace{ \left[\begin{array}{cc} P(\lambda)\\ N(\lambda) \end{array}\right] }_{\beta} \delta_I }$

の解として $\delta_R,\delta_I$ を求めますが、 $\alpha$ と $\beta$ が独立ではないので、その求解には工夫が必要です。ここでは、Edwards and Spurgeon に従って、その手順を述べます。

まず、行列 $\alpha$ に対して、 $\alpha^T$ の零化空間の基底行列を $K_\alpha$ で表すと、 $\alpha^TK_\alpha=0$ より $K_\alpha^T\alpha=0$ を得て、 $K_\alpha^T$ の零化空間の基底行列は $\alpha$ 、すなわち ${\cal N}(K_\alpha^T)={\cal R}(\alpha)$ となることに注意します。同様に行列 $\beta$ に対して、 $\beta^T$ の零化空間の基底行列を $K_\beta$ で表すと、 $K_\beta^T$ の零化空間の基底行列は $\beta$ 、すなわち ${\cal N}(K_\beta^T)={\cal R}(\beta)$ となります。そこで、

$\displaystyle{(20)\quad \underbrace{ \left[\begin{array}{c} K_\alpha^T\\ K_\beta^T \end{array}\right] }_{\gamma^T} K_\gamma=0 \Rightarrow {\cal R}(K_\gamma)\subset {\cal N}(K_\alpha^T)\cap {\cal N}(K_\beta^T)={\cal R}(\alpha)\cap {\cal R}(\beta) }$

に注意して、 $\gamma=[K_\alpha,K_\beta]$ の零化空間の基底行列 $K_\gamma$ を用いて、望ましい $v_R,v_I$ を与えて

$\displaystyle{(21)\quad \left[\begin{array}{c} v_R \\ v_I \end{array}\right]= K_\gamma \delta }$

を解いて $\delta$ を求めます。ただし、望ましい $v_R,v_I$ のすべての成分を指定できるわけではないので

$\displaystyle{(22)\quad \tilde{R}K_\gamma\delta=\tilde{R}v =\left[\begin{array}{c} \tilde{v}_R^s\\ \tilde{v}_I^s\\ \tilde{v}_R^u\\ \tilde{v}_I^u\\ \end{array}\right] \Rightarrow \left[\begin{array}{cc} I_{2q} & 0 \end{array}\right] \tilde{R}K_\gamma\times\delta= \left[\begin{array}{c} \tilde{v}_R^s\\ \tilde{v}_I^s\\ \end{array}\right] }$

の最小２乗解として $\delta=\delta^*$ を定めます。ここで、セレクター行列 $\tilde{R}$ により、望ましい $v_R,v_I$ の指定された $2q$ 個の成分だけを集めた $\tilde{v}_R^s,\tilde{v}_I^s$ が得られるとしています。

このとき、割当可能な $v_R^*,v_I^*$ は次式で与えられます。

$\displaystyle{(23)\quad \left[\begin{array}{c} v_R^*\\ v_I^* \end{array}\right]= K_\gamma\delta^* }$

これは(14)の $\delta_R^*,\delta_I^*$ を、それぞれ $\alpha^{\dag}K_\gamma\delta^*$ と $\beta^{\dag}K_\gamma\delta^*$ で定めれば（ $\dag$ は疑似逆行列を表す）、対応する $\xi_R^*,\xi_I^*$ は次式で与えられます。

$\displaystyle{(24)\quad \left[\begin{array}{c} \xi_R^*\\ \xi_I^* \end{array}\right]= \left[\begin{array}{c} M(\lambda)\alpha^{\dag}K_\gamma\delta^*\\ M(\lambda)\beta^{\dag}K_\gamma\delta^*\\ \end{array}\right] }$

[3] 複素共役固有値の場合（別のアプローチ）

(17)は次式のようにも表すことができます。

$\displaystyle{(25)\quad \begin{array}{l} (A-BF)(v_R\pm jv_I)=(\lambda_R\pm j\lambda_I)(v_R\pm jv_I)\\ \Rightarrow (A-BF)v_R=\lambda_Rv_R-\lambda_Iv_I,\ (A-BF)v_I =\lambda_Rv_I+\lambda_Iv_R \\ \Rightarrow %\underbrace{ \left[\begin{array}{cccc} \lambda_R I-A & -\lambda_I I & B & 0\\ \lambda_I I & \lambda_R I-A & 0 & B\\ \end{array}\right] %}_{G(\lambda)} \left[\begin{array}{cc} v_R \\ v_I \\ \xi_R \\ \xi_I \end{array}\right] =0 \quad(\xi_R=Fv_R,\xi_I=Fv_I) \end{array} }$

したがって、[2]において

$\displaystyle{(26)\quad \begin{array}{l} G(\lambda)\leftarrow \left[\begin{array}{cccc} \lambda_R I-A & -\lambda_I I & B & 0\\ \lambda_I I & \lambda_R I-A & 0 & B \end{array}\right]\\ v\leftarrow\left[\begin{array}{c} v_R\\ v_I \end{array}\right],\ \xi\leftarrow\left[\begin{array}{c} \xi_R\\ \xi_I \end{array}\right],\ \delta\leftarrow\left[\begin{array}{c} \delta_R\\ \delta_I \end{array}\right] \end{array} }$

と置き換えて同様の議論ができそうです（[2]の議論と等価であるかは検討中です）。

●もし $A-\lambda_R I_n$ が正則ならば、公式

$\displaystyle{(27)\quad \begin{array}{l} \left[\begin{array}{cc} M_{11} & M_{12} \\ M_{21} & M_{22} \end{array}\right]^{-1} = \left[\begin{array}{cc} \Phi^{-1} & -\Phi^{-1}M_{12}M_{22}^{-1} \\ -M_{22}^{-1}M_{21}\Phi^{-1} & M_{22}^{-1}+M_{22}^{-1}M_{21}\Phi^{-1}M_{12}M_{22}^{-1} \end{array}\right]\\ \Phi\stackrel{\rm def}{=}M_{11}-M_{12}M_{22}^{-1}M_{21}\ (M_{22}\ {\rm nonsingular}) \end{array} }$

を参照すると

$\displaystyle{(28)\quad \left[\begin{array}{cc} \lambda_R I-A & -\lambda_I I\\ \lambda_I I & \lambda_R I-A \end{array}\right]^{-1} }$

を計算できます。そこで

$\displaystyle{(29)\quad \begin{array}{l} \left[\begin{array}{cc} \lambda_R I-A & -\lambda_I I\\ \lambda_I I & \lambda_R I-A \end{array}\right] \underbrace{ \left[\begin{array}{cc} \lambda_R I-A & -\lambda_I I\\ \lambda_I I & \lambda_R I-A \end{array}\right]^{-1} \left[\begin{array}{cc} B & 0\\ 0 & B \end{array}\right] }_{N(\lambda)}\\ + \left[\begin{array}{cc} B & 0\\ 0 & B \end{array}\right] \underbrace{(-I_{2m})}_{M(\lambda)}=0 \end{array} }$

より、 $N(\lambda)$ と $M(\lambda)$ が特定できます。(8)に相当する式より

$\displaystyle{(30)\quad \boxed{ \left[\begin{array}{c} v_R\\ v_I \end{array}\right] = \left[\begin{array}{cc} A-\lambda_R I & \lambda_I I\\ -\lambda_I I & A-\lambda_R I \end{array}\right]^{-1} \left[\begin{array}{cc} B & 0\\ 0 & B \end{array}\right] \left[\begin{array}{c} \xi_R\\ \xi_I \end{array}\right]} }$

となって、指定した固有ベクトルの実部 $v_R$ と虚部 $v_I$ を近似する $\xi_R$ と $\xi_I$ を直接求めることができます。

●具体的に、(1)に対して、固有値 $-1\pm j$ を持たせる状態フィードバックを、対応する固有ベクトルが

$\displaystyle{(31)\quad \underbrace{ \left[\begin{array}{c} 1\\ 0 \end{array}\right]}_{v_R} \pm j \underbrace{ \left[\begin{array}{c} 0\\ 1 \end{array}\right]}_{v_I} }$

となるように求めてみます。

$\displaystyle{(32)\quad \begin{array}{l} \underbrace{ \left[\begin{array}{cccc} 0-(-1) & 0 & 1 & 0\\ 0 & -1-(-1) & 0 & 1\\ -1 & 0 & 0-(-1) & 0\\ 0 & -1 & 0 & -1-(-1) \end{array}\right]^{-1} }_{\left[\begin{array}{cc} A-\lambda_R I & \lambda_I I\\ -\lambda_I I & A-\lambda_R I \end{array}\right]^{-1}}\\ \times\underbrace{ \left[\begin{array}{cccc} 1 & 1 & 0 & 0\\ 1 & -1& 0 & 0\\ 0 & 0 & 1 & 1\\ 0 & 0 & 1 & -1 \end{array}\right] }_{\left[\begin{array}{cc} B & 0\\ 0 & B \end{array}\right]} \left[\begin{array}{c} \xi_R\\ \xi_I \end{array}\right] =\underbrace{ \left[\begin{array}{c} 1\\ 0\\ 0\\ 1 \end{array}\right]}_{\left[\begin{array}{c} v_R\\ v_I \end{array}\right]}\\ \Rightarrow \left[\begin{array}{c} \xi_R\\ \xi_I \end{array}\right]= \left[\begin{array}{cccc} 1/2 & 1/2 & -1/2 & -1/2\\ 0 & 0 & -1 & 1\\ 1/2 & 1/2 & 1/2 & 1/2\\ 1 & -1& 0 & 0 \end{array}\right]^{-1} \left[\begin{array}{c} 1\\ 0\\ 0\\ 1 \end{array}\right]= \left[\begin{array}{c} 1\\ 0\\ -1/2\\ -1/2 \end{array}\right] \end{array} }$

状態フィードバックゲイン $F$ は

$\displaystyle{(33)\quad \begin{array}{l} \left[\begin{array}{cc} \xi_R & \xi_I \end{array}\right] =F \left[\begin{array}{cc} v_R & v_I \end{array}\right]\\ \Rightarrow F= \left[\begin{array}{cc} \xi_R & \xi_I \end{array}\right] \left[\begin{array}{cc} v_R & v_I \end{array}\right]^{-1} = \left[\begin{array}{cc} 1 & -1/2\\ 0 & -1/2 \end{array}\right] \end{array} }$

演習C42…Flipped Classroom
$1^\circ$ 　(16)と(33)の計算を行う次のプログラムを実行せよ。

MATLAB

%test_vplace2.m
%-----
 clear all, close all
 A=[0 0;
    0 -1];
 B=[1 1;
    1 -1];
%-----
 lambda=[-2 -3];
 nocomp=0;
%固有ベクトルの要素の指定の有無 (１：指定、０：任意)
 specpos=[ 1  1; 
           1  1];       
%固有ベクトルの要素の値 (０：零指定、１：非零指定、－１：非零任意)       
 specent=[ 1  0;
           0  1];
 F=vplace2(A,B,lambda,nocomp,specpos,specent)
 [V,pl]=eig(A-B*F)
 VV=V*diag([1/V(2,1),1/V(1,2)])
 disp("---")
%-----
 lambda=[-1 1];
 nocomp=1;
%固有ベクトルの要素の指定の有無 (１：指定、０：任意)
 specpos=[ 1  1; 
           1  1];        
%固有ベクトルの要素の値 (０：零指定、１：非零指定、－１：非零任意)       
 specent=[ 1  0;
           0  1];        
 F=vplace2(A,B,lambda,nocomp,specpos,specent) 
 [V,pl]=eig(A-B*F) 
 VV=V*diag([1/V(1,1),1/V(1,2)])

$2^\circ$ 　１入力系の場合は、specposとspecentは乱数で与えてよい、なぜか？

MATLAB

%test_vplace2a.m
%-----
 clear all, close all
 A=[0 1;
    0 0];
 B=[0;
    1];
%-----
%lambda=[-1 -1];
 lambda=[-2 -3];  
 nocomp=0;
%固有ベクトルの要素の指定の有無 (１：指定、０：任意)
 specpos=rand(2,2)
%固有ベクトルの要素の値 (０：零指定、１：非零指定、－１：非零任意)       
 specent=rand(2,2) 
 F=vplace2(A,B,lambda,nocomp,specpos,specent) 
 [V,pl]=eig(A-B*F)
%-----
 lambda=[-1 1];
 nocomp=1;
%固有ベクトルの要素の指定の有無 (１：指定、０：任意)
 specpos=rand(2,2)
%固有ベクトルの要素の値 (０：零指定、１：非零指定、－１：非零任意)       
 specent=rand(2,2) 
 F=vplace2(A,B,lambda,nocomp,specpos,specent) 
 [V,pl]=eig(A-B*F) 
%-----
%eof

Note C42 固有ベクトル設定プログラム
これはEdwards and Spurgeonによるプログラムspvplを基に作成できます。

スイッチング関数(2)

投稿日時: 2021年10月29日投稿者: cacsd

スイッチング関数(2)

[1] 状態方程式とスイッチング関数は次の標準形をとるように座標変換されていると仮定します。

$\displaystyle{(1)\quad \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot x_2(t) \end{array}\right] }_{\dot x(t)} = \underbrace{ \left[\begin{array}{cc} A_{11} & A_{12} \\ A_{21} & A_{22} \\ \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} 0\\ B_2 \end{array}\right] }_{B} u(t) }$

$\displaystyle{(2)\quad s(t)= \underbrace{ \left[\begin{array}{cc} S_1 & S_2 \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x(t)} = \underbrace{S_2 \left[\begin{array}{cc} M & I \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x(t)} \ (M=S_2^{-1}S_1) }$

このとき、等価制御

$\displaystyle{(3)\quad u_{eq}(t)=-(SB)^{-1}SAx(t) }$

による閉ループ系は次式で表されました。

$\displaystyle{(4)\quad \dot{x}(t)=\underbrace{(I_n-B(SB)^{-1}S)A}_{A_{eq}}x(t)\quad(t\ge t_s) }$

ここで $A$ 、 $B$ 、 $S$ は、(1)と(2)で与えられるとすると

$\displaystyle{(5)\quad \begin{array}{l} A_{eq}= \underbrace{ \left[\begin{array}{cc} I & 0 \\ -M & 0 \\ \end{array}\right] }_{I_n-B(SB)^{-1}S} \left[\begin{array}{cc} A_{11} & A_{12} \\ A_{21} & A_{22} \\ \end{array}\right]= \left[\begin{array}{cc} A_{11} & A_{12} \\ -MA_{11} & -MA_{12} \\ \end{array}\right]\\ = \left[\begin{array}{cc} I & 0 \\ -M & I \\ \end{array}\right] \left[\begin{array}{cc} \bar{A}_{11} & A_{12} \\ 0 & 0 \\ \end{array}\right] \left[\begin{array}{cc} I & 0 \\ -M & I \\ \end{array}\right]^{-1} \end{array} }$

$\displaystyle{(6)\quad SA_{eq}=0\Rightarrow SA_{eq}v_i=0\Rightarrow \lambda_iSv_i=0 \Rightarrow Sv_i=0 }$

すなわち

$\displaystyle{(7)\quad { \boxed{S[v_1\cdots v_{n-m}]=0 %　\Rightarrow [v_1\cdots v_{n-m}]^TS^T=0}} }$

が成り立ちます。したがって、まず等価制御による閉ループ系の $A_{eq}$ が望ましい固有値・固有ベクトルをもつように設計して、(7)から $S$ を決定することが考えられます。

● １入力系の場合は固有ベクトル設定の自由度はありませんが、多入力系の場合は固有ベクトルの向きを設定する必要があり、これは元の状態方程式(状態変数の並び)に対して与えます。そこで(7)が元の状態方程式ではどう変わるかを調べておきます。

元の状態方程式

$\displaystyle{(8)\quad \begin{array}{l} \dot{x}(t)=Ax(t)+Bu(t)+f(t,x,u)\\ (x(t)\in{\rm\bf R}^n, u(t)\in{\rm\bf R}^m, f(t,x,u)\in{\rm\bf R}^n) \end{array} }$

に対する標準形(1)は、適当な直交行列 $T_r$ 座標による座標変換

$\displaystyle{(9)\quad x_{r}(t)=T_rx(t) \Leftrightarrow x(t)=T_r^Tx_r(t) }$

を行って（ $T_r^{-1}=T_r^T$ に注意）

$\displaystyle{(10)\quad \begin{array}{l} \underbrace{ \left[\begin{array}{c} \dot{x}_{r1}(t)\\ \dot{x}_{r2}(t) \end{array}\right] }_{\dot{x}_r(t)} = \underbrace{ \left[\begin{array}{cc} A_{r11} & A_{r12} \\ A_{r21} & A_{r22} \\ \end{array}\right] }_{A_r=T_rAT_r^T} \underbrace{ \left[\begin{array}{c} x_{r1}(t)\\ x_{r2}(t) \end{array}\right] }_{x_r(t)} + \underbrace{ \left[\begin{array}{c} 0\\ \Sigma_1V^T \end{array}\right] }_{B_r=T_rB} u(t)\\ + \underbrace{ \left[\begin{array}{l} f_{ru}(t,x)\\ f_{rm}(t,x,u) \end{array}\right] }_{f_r(t,x,u)=T_rf(t,x,u)} \end{array} }$

のように得られたものでした。また(2)は次式に

$\displaystyle{(11)\quad s(t)= \underbrace{ \left[\begin{array}{cc} S_{r1} & S_{r2} \\ \end{array}\right] }_{S_r=ST_r^T} \underbrace{ \left[\begin{array}{c} x_{r1}(t)\\ x_{r2}(t) \end{array}\right] }_{x_r(t)=T_rx(t)} }$

に、(7)は次式に相当します。

$\displaystyle{(12)\quad S_r[v_{r1}\cdots v_{r,n-m}]=0 }$

いま、元の状態空間での固有ベクトル $v_1\cdots v_{n-m}$ は、座標変換後には

$\displaystyle{(13)\quad Av_i=\lambda_iv_i\Rightarrow T_rAT_r^TT_rv_i=\lambda_iT_rv_i\Rightarrow v_{ri}=T_rv_i }$

に変わることに注意します。(12)から

$\displaystyle{(14)\quad S_rT_rT_r^T[v_{r1}\cdots v_{r,n-m}]=S[v_1\cdots v_{n-m}]=0 }$

を得ます。これから元の状態方程式に対するスイッチング関数と固有ベクトルについても、(12)と同様の関係が成り立つと言えます。すなわち固有ベクトル $v_1\cdots v_{n-m}$ は元の状態空間で与えておいて、(14)から $S$ を決定してよいことがわかります。

演習C41…Flipped Classroom

$1^\circ$ スイッチング関数決定プログラムswfvplを用いて、図1の台車駆動型倒立振子を安定化する次のSMC則を求めよ。

$\displaystyle{{ \boxed{\begin{array}{l} u(t)=\underbrace{u_\ell(t)}_{linear\ control}+\underbrace{u_n(t)}_{switching\ component}\\ u_\ell(t)=-\underbrace{(SB)^{-1}(SA-\Phi S)}_{L=L_{eq}+L_{phi}}x(t) \quad(\Phi:stable\ matrix)\\ u_n(t)=-\underbrace{(SB)^{-1}\rho(t,x)}_{L_n}\frac{P_2s(t)}{||P_2s(t)||}\quad(P_2\Phi+\Phi^TP_2=-I) \end{array}}} }$

図1 台車駆動型倒立振子

MATLAB

%cCIP_smc2.m
%-----
 clear all, close all
 global mc m ell g th0
 mc=1; m=0.1; ell=0.2; g=9.8;
 ths=0; %input('ths   = <0,180> ')/180*pi;
 th0=3/180*pi; %input('th(0) = <0,180> ')/180*pi;
%-----
 A=[zeros(2,2) eye(2);zeros(2,4)];
 A(3,1)=0; 
 A(3,2)=-3*m*g/(m+4*mc); 
 A(4,1)=0; 
 A(4,2)=3*(m+mc)*g/((m+4*mc)*ell);  
 B=zeros(4,1);
 B(3)=4/(m+4*mc);
 B(4)=-3/((m+4*mc)*ell); 
 C=diag([100 180/pi])*eye(2,4); 
%-----
 lambda=[-0.8 0.5 -3];
 nocomp=1;
 specpos=rand(4,3)
 specent=rand(4,3)
 S=swfvpl(A,B,lambda,nocomp,specpos,specent)
%-----
 Phi=-0.1;
 P2=0.5*inv(-Phi);
 Check=P2*Phi+Phi*P2
 P2S=P2*S;
 Leq=inv(S*B)*S*A;
 LPhi=-inv(S*B)*Phi*S;
 L=Leq+LPhi;
 rho=0.6;
 Ln=inv(S*B)*rho;
%-----
 x0=[0;th0;0;0];
 sim('CIP_smc_2015a.mdl')
%-----
%eof

Note C41 スイッチング関数決定プログラム２
Edwards and Spurgeon によって、(26)から、スイッチング関数の行列 $S$ を求めるプログラムswfvplが開発されています。

MATLAB

%test_swfvpl.m
%-----
 clear all, close all
 A=[0 1 0;
    0 0 1;
    0 0 0];
 B=[0;
    0;
    1];
%-----
 lambda=[-2 -3];
 nocomp=0;
%固有ベクトルの要素の指定の有無 (１：指定、０：任意)
 specpos=rand(3,3)
%固有ベクトルの要素の値 (０：零指定、１：非零指定、－１：非零任意)       
 specent=rand(3,3) 
 S=swfvpl(A,B,lambda,nocomp,specpos,specent) 
 Leq=inv(S*B)*S*A;
 [V,pl]=eig(A-B*Leq) 
%-----
 lambda=[-1 1];
 nocomp=1;
%固有ベクトルの要素の指定の有無 (１：指定、０：任意)
 specpos=rand(3,3)
%固有ベクトルの要素の値 (０：零指定、１：非零指定、－１：非零任意)       
 specent=rand(3,3) 
 S=swfvpl(A,B,lambda,nocomp,specpos,specent) 
 Leq=inv(S*B)*S*A;
 [V,pl]=eig(A-B*Leq) 
%-----
%eof

補遺：LQ制御

投稿日時: 2021年10月29日投稿者: cacsd

LQ制御…Homework

[1]　可制御な制御対象

$\displaystyle{(1)\quad %\left\{\begin{array}{ll} \dot{x}(t)=Ax(t)+Bu(t)\quad (x(t)\in{\bf R}^n,u(t)\in{\bf R}^m)\\ %z(t)=Cx(t)&(z(t)\in{\bf R}^q) %\end{array}\right. }$

を安定化する状態フィードバック

$\displaystyle{(2)\quad u(t)=-Fx(t) }$

の決定法を考えます。一つの方法は，閉ループ系

$\displaystyle{(3)\quad %\left\{\begin{array}{l} \dot{x}(t)=(A-BF)x(t) \\ %z(t)=Cx(t) %\end{array}\right. }$

の時間応答に関する評価規範として，２次形式評価関数

$(4)\quad \boxed{\begin{array}{l} \displaystyle{J=\int_0^\infty (z^T(t)Qz(t)+u^T(t)Ru(t))\,dt\quad (Q>0,R>0)}\\ \displaystyle{z(t)=Cx(t)\quad (z(t)\in{\bf R}^q)} \end{array}}$

を設定し，これを最小化する問題を解くことです。ただし、 $(A,C)$ は可観測対とします。これによる状態フィードバックのゲイン行列 $F$ は，リッカチ方程式

$\displaystyle{(5)\quad {\boxed{\Pi A+A^T\Pi-\Pi BR^{-1}B^T\Pi+C^TQC=0}} }$

の解 $\Pi>0$ を用いて，次式で与えられます。

$\displaystyle{(6)\quad {\boxed{F=R^{-1}B^T\Pi}} }$

この証明はNoteに示しています。

Note A52-1 行列による微分　

いま、任意の行列 $X$ の $(i,j)$ 要素を $x_{ij}=[X]_{ij}$ で表すとき、スカラ関数 $f$ を行列変数 $X$ の各要素で微分して得られる行列を $\displaystyle{[\frac{\partial f}{\partial X}]_{ij}=\frac{\partial f}{\partial x_{ij}}}$ で定義します。このとき、行列のトレースについて、次が成り立ちます。

$\displaystyle{(1)\ {\rm tr}AB={\rm tr}BA}$
$\displaystyle{(2)\ \frac{\partial}{\partial X}{\rm tr}AXB=A^TB^T}$
$\displaystyle{(3)\ \frac{\partial}{\partial X}{\rm tr}AX^TB=BA}$
$\displaystyle{(4)\ \frac{\partial}{\partial X}{\rm tr}AXBX^T=A^TXB^T+AXB}$

実際、

$\displaystyle{(1)\ \sum_{i}[AB]_{ii}=\sum_{i}\sum_{j}a_{ij}b_{ji} =\sum_{j}\sum_{i}b_{ji}a_{ij}=\sum_{j}[BA]_{jj}}$

$\displaystyle{(2)\ [\frac{\partial}{\partial X}{\rm tr}AXB]_{ij} =\frac{\partial}{\partial x_{ij}}\sum_{k}[AXB]_{kk} =\frac{\partial}{\partial x_{ij}}\sum_{k}\sum_{i,j}a_{ki}x_{ij}b_{jk}}$
$\displaystyle{=\sum_{k}b_{jk}a_{ki}=[BA]_{ji}=[A^TB^T]_{ij}}$

$\displaystyle{(3)\ [\frac{\partial}{\partial X}{\rm tr}AX^TB]_{ij} =\frac{\partial}{\partial x_{ij}}\sum_{k}[AX^TB]_{kk} =\frac{\partial}{\partial x_{ij}}\sum_{k}\sum_{i,j}a_{ki}x_{ji}b_{jk}}$
$\displaystyle{=\sum_{k}b_{ik}a_{kj}=[BA]_{ij}}$

$\displaystyle{(4)\ \frac{\partial}{\partial X}{\rm tr}AXBX^T =\frac{\partial}{\partial X}{\rm tr}(AXB)X^T +\frac{\partial}{\partial X}{\rm tr}AX(BX^T)}$
$\displaystyle{=AXB+A^TXB^T}$

Note A52-2 LQ制御問題の解法　

可制御かつ可観測な $n$ 次系

$\displaystyle{(1)\quad \dot{x}(t)=Ax(t)+Bu(t),\ y(t)=Cx(t) }$

に対する状態フィードバック

$\displaystyle{(2)\quad u(t)=-Fx(t) }$

による閉ループ系

$\displaystyle{(3)\quad \begin{array}{l} \dot{x}(t)=(A-BF)x(t),\ z(t)=Cx(t)\\ A_F=A-BF:\ stable\ matrix \end{array} }$

に対して、評価関数

$\displaystyle{(4)\quad J=\int_0^\infty(z^T(t)Qz(t)+u^T(t)Ru(t))\,dt }$

を最小化するように $F$ を決める問題を考えます。

閉ループ系における状態の振る舞いは次式で与えられます。

$\displaystyle{(5)\quad x(t)=\exp({A_Ft})x(0) }$

ここで、１次系の場合は初期状態は $x(0)\ne0$ であればよかったのですが、一般の場合はインパルス応答となるように $B$ の列ベクトル $B^{(i)}$ を考えます。各インパルス応答

$\displaystyle{(6)\quad x^{(i)}(t)=\exp({A_Ft})B^{(i)} }$

に対する評価関数 $J$ の総和は

$(7)\quad \begin{array}{l} \displaystyle{\sum\int_0^\infty x^{(i)}\,^T(t)(C^TQC+F^TRF)x^{(i)}(t)\,dt}\\ \displaystyle{=\sum B^{(i)}\,^T\underbrace{\int_0^\infty \exp(A_F^Tt)(C^TQC+F^TRF)\exp(A_Ft)\,dt}_{\Pi}B^{(i)}}\\ \displaystyle{=\sum B^{(i)}\,^T\Pi B^{(i)}={\rm tr}\ \Pi BB^T} \end{array} }$

と書けます。いま $\Gamma$ をラグランジュの未定定数として

$\displaystyle{(8)\quad J'={\rm tr}\ \Pi x(0)x^T(0)+{\rm tr}\ ((\Pi A_F+A_F^T\Pi+C^TQC+F^TRF)\Gamma) }$

を最小化する問題を考えます。ここで、制約条件は、リャプノフ方程式と呼ばれる

$\displaystyle{(9)\quad \Pi A_F+A_F^T\Pi+C^TQC+F^TRF=0 }$

ですが、 $\Pi>0$ は $A_F$ が安定行列を意味することに注意します。

そこで、必要条件として次を得ます。

$\displaystyle{(11)\quad \frac{\partial J'}{\partial\Pi}&=&BB^T+\Gamma A_F^T+A_F\Gamma=0\Rightarrow \Gamma>0 }$
$\displaystyle{(12)\quad \frac{\partial J'}{\partial F}&=&2(RF-B^T\Pi)\Gamma=0\Rightarrow F=R^{-1}B^T\Pi }$
$\displaystyle{(13)\quad \frac{\partial J'}{\partial\Gamma}&=&\Pi A_F+A_F^T\Pi+C^TQC+F^TRF=0 }$

ここで、第２式から得られる $F=R^{-1}B^T\Pi$ を第３式に代入して

$\displaystyle{(14)\quad \begin{array}{l} \Pi (A-BR^{-1}B^T\Pi)+(A-BR^{-1}B^T\Pi)^T\Pi+C^TQ\nonumber\\ +(R^{-1}B^T\Pi)^TRR^{-1}B^T\Pi=0 \end{array} }$

すなわち、リッカチ方程式と呼ばれる $\Pi$ の行列方程式

$\displaystyle{(15)\quad \Pi A+A^T\Pi-\Pi BR^{-1}B^T\Pi+C^TQC=0 }$

を得ます。これから $\Pi>0$ を求めて、 $F$ は

$\displaystyle{(16)\quad F=R^{-1}B^T\Pi }$

のように得られます。このような制御方式をLQ制御と呼びます。

一方、十分性の議論は次のように行われます。まず、被積分項は次のように表すことができます。

$\displaystyle{(17)\quad \begin{array}{lll} &&y^TQy+u^TRu=(u+R^{-1}B^T\Pi x)^TR(u+R^{-1}B^T\Pi x)\nonumber\\ &&-\underbrace{\frac{d}{dt}x^T\Pi x}_{2x^T\Pi\dot{x}} \end{array} }$

実際、右辺に $\dot{x}=Ax+Bu$ を代入し、リッカチ方程式を用いると

$\displaystyle{(18)\quad \begin{array}{lll} &&{\rm RHS}=u^TRu+2x^T\Pi Bu+x^T\Pi BR^{-1}B^T\Pi x-2x^T\Pi(Ax+Bu)\nonumber\\ &&=u^TRu+x^T(\Pi A+A^T\Pi+C^TQC)x-2x^T\Pi Ax={\rm LHS} \end{array} }$

したがって、上記の両辺を積分して

$\displaystyle{(19)\quad J=\int_0^\infty (u+R^{-1}B^T\Pi x)^TR(u+R^{-1}B^T\Pi x)\,dt-\left[x^T\Pi x\right]_0^\infty }$

を得ます。ここで、 $x(t)=\exp(A_Ft)x(0)\rightarrow 0\ (t\rightarrow\infty)$ を前提とするので

$\displaystyle{(20)\quad J=\int_0^\infty (u+R^{-1}B^T\Pi x)^TR(u+R^{-1}B^T\Pi x)\,dt+x^T(0)\Pi x(0) }$

を得ます。これから $u=-R^{-1}B^T\Pi x$ が評価関数を最小化することが分かります。