现代控制理论 ——03 状态方程的解

证明

test

证明

eAt=I+At+12!A2t2+e^{\bold{A}t}=\bold{I}+\bold{A}t+\frac{1}{2!}\bold{A}^2t^2+\cdots

=[100010001]+[λ1000λ2000λn]t+12![λ12000λ22000λn2]t2+=\begin{bmatrix} 1&0&\cdots&0\\ 0&1&\cdots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&1\end{bmatrix}+\begin{bmatrix} \lambda_1&0&\cdots&0\\ 0&\lambda_2&\cdots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&\lambda_n\end{bmatrix}t+\frac{1}{2!}\begin{bmatrix} \lambda_1^2&0&\cdots&0\\ 0&\lambda_2^2&\cdots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&\lambda_n^2\end{bmatrix}t^2+\cdots

=[n=0+1n!λ1ntn000n=0+1n!λ2ntn000n=0+1n!λnntn]=[eλ1t000eλ2t000eλnt]=\begin{bmatrix} \sum^{+\infty}_{n=0}\frac{1}{n!}\lambda_1^nt^n&0&\cdots&0\\ 0&\sum^{+\infty}_{n=0}\frac{1}{n!}\lambda_2^nt^n&\cdots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&\sum^{+\infty}_{n=0}\frac{1}{n!}\lambda_n^nt^n\end{bmatrix}=\begin{bmatrix} e^{\lambda_1t}&0&\cdots&0\\ 0&e^{\lambda_2t}&\cdots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&e^{\lambda_nt}\end{bmatrix}

  1. A\bold{A}m×mm\times{m} 的若尔当块,即A=[λ1000λ10λ00λ1000λ]m×m\bold{A}=\begin{bmatrix} \lambda&1&0&\cdots&0\\ 0&\lambda&1&\cdots&0\\ \vdots&\vdots&\lambda&\ddots&\vdots\\0&0&\cdots&\lambda&1\\ 0&0&\cdots&0&\lambda\end{bmatrix}_{m\times{m}},则eAt=eλt[1tt22!tm1(m1)!01ttm2(m2)!001t0001]m×me^{\bold{A}t}=e^{\lambda t}\begin{bmatrix} 1&t&\frac{t^{2}}{2!}&\cdots&\frac{t^{m-1}}{(m-1)!}\\ 0&1&t&\cdots&\frac{t^{m-2}}{(m-2)!}\\ \vdots&\vdots&\ddots&\ddots&\vdots\\0&0&\cdots&1&t\\ 0&0&\cdots&0&1\end{bmatrix}_{m\times{m}}
若尔当块

形如[λ1000λ10λ00λ1000λ]m×m\begin{bmatrix} \lambda&1&0&\cdots&0\\ 0&\lambda&1&\cdots&0\\ \vdots&\vdots&\lambda&\ddots&\vdots\\0&0&\cdots&\lambda&1\\ 0&0&\cdots&0&\lambda\end{bmatrix}_{m\times{m}}mm 阶若尔当矩阵,1 阶若尔当矩阵为λ\lambda

  1. A\bold{A} 为一个有多个若尔当块的若尔当矩阵(即若当标准型),即A=[A1000A2000An]\bold{A}=\begin{bmatrix} \bold{A}_1&0&\cdots&0\\ 0&\bold{A}_2&\cdots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&\bold{A}_n\end{bmatrix},则eAt=[eA1t000eA2t000eAnt]e^{\bold{A}t}=\begin{bmatrix} e^{\bold{A}_1t}&0&\cdots&0\\ 0&e^{\bold{A}_2t}&\cdots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&e^{\bold{A}_nt}\end{bmatrix}

# 矩阵指数的计算

  1. 定义计算:eAt=n=0+1n!An(t)ne^{\bold{A}t}=\sum^{+\infty}_{n=0}\frac{1}{n!}\bold{A}^n(t)^n。该方法适用于计算机运算。
例题

已知A=[0110]\bold{A}=\begin{bmatrix} 0&1\\-1&0\end{bmatrix},求eAte^{\bold{A}t}

由定义,

eAt=I+At+12!+=[1001]+[0tt0]+12![t200t2]+=[1t22!+tt33!+t+t33!1t22!+]=[costsintsintcost]\begin{aligned}e^{\bold{A}t}&=\bold{I}+\bold{A}t+\frac{1}{2!}+\cdots=\begin{bmatrix}1&0\\0&1\end{bmatrix}+\begin{bmatrix}0&t\\-t&0\end{bmatrix}+\frac{1}{2!}\begin{bmatrix}-t^2&0\\0&-t^2\end{bmatrix}+\cdots\\&=\begin{bmatrix}1-\frac{t^2}{2!}+\cdots&t-\frac{t^3}{3!}+\cdots\\-t+\frac{t^3}{3!}-\cdots&1-\frac{t^2}{2!}+\cdots\end{bmatrix}=\begin{bmatrix} \cos{t}&\sin{t}\\-\sin{t}&\cos{t}\end{bmatrix}\end{aligned}.

  1. 拉氏变换法:利用拉氏变换在频域中求解齐次状态方程的解。

设线性时不变齐次状态方程为 x˙=Ax(t)\bold{\dot{x}}=\bold{Ax}(t)x(0)=x0\bold{x}(0)=\bold{x}_0tt0t\geq{t_0}

作拉氏变换有 sX(s)x(0)=AX(s)s\bold{X}(s)-\bold{x}(0)=\bold{AX}(s),即 (sIA)X(s)=x(0)(s\bold{I}-\bold{A})\bold{X}(s)=\bold{x}(0),那么

X(s)=(sIA)1x(0)\bold{X}(s) =(s\bold{I}-\bold{A})^{-1}\bold{x}(0)

取拉氏逆变换有 x(0)=L1[(sIA)1x(0)]=L1[(sIA)1]x(0)\bold{x}(0)=L^{-1}[(s\bold{I}-\bold{A})^{-1}\bold{x}(0)]=L^{-1}[(s\bold{I}-\bold{A})^{-1}]\bold{x}(0),因此

eAt=L1[(sIA)1](2)e^{\bold{A}t}=L^{-1}[(s\bold{I}-\bold{A})^{-1}] \tag{2}

例题

计算矩阵 A=[0123]\bold{A}=\begin{bmatrix}0&1\\-2&-3\end{bmatrix} 的矩阵指数。

由拉氏变换法,(sIA)=[s12s+3](s\bold{I}-\bold{A})=\begin{bmatrix} s&-1\\2&s+3\end{bmatrix},则(sIA)1=[s+3(s+1)(s+2)1(s+1)(s+2)2(s+1)(s+2)s(s+1)(s+2)](s\bold{I}-\bold{A})^{-1}=\begin{bmatrix}\frac{s+3}{(s+1)(s+2)}&\frac{1}{(s+1)(s+2)}\\\frac{-2}{(s+1)(s+2)}&\frac{s}{(s+1)(s+2)}\end{bmatrix}

eAt=L1[s+3(s+1)(s+2)1(s+1)(s+2)2(s+1)(s+2)s(s+1)(s+2)]=[2etete2t2et+2e2tet+2e2t]e^{\bold{A}t}=L^{-1}\begin{bmatrix}\frac{s+3}{(s+1)(s+2)}&\frac{1}{(s+1)(s+2)}\\\frac{-2}{(s+1)(s+2)}&\frac{s}{(s+1)(s+2)}\end{bmatrix}=\begin{bmatrix}2e^{-t}&e^{-t}-e^{-2t}\\-2e^{-t}+2e^{-2t}&-e^{-t}+2e^{-2t}\end{bmatrix}

  1. 将矩阵化为对角标准型或若尔当标准型。

A=[λ1000λ2000λn]\bold{A}=\begin{bmatrix} \lambda_1&0&\cdots&0\\ 0&\lambda_2&\cdots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&\lambda_n\end{bmatrix} 为对角矩阵,则eAte^{\bold{A}t} 也为对角矩阵(性质 8),即eAt=[eλ1t000eλ2t000eλnt]e^{\bold{A}t}=\begin{bmatrix} e^{\lambda_1t}&0&\cdots&0\\ 0&e^{\lambda_2t}&\cdots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&e^{\lambda_nt}\end{bmatrix}

(1)当矩阵A\bold{A}nn 个特征值 λ1,λ2λn\lambda_1,\lambda_2\dots\lambda_n 均两两互异时,则可确定变换阵 P\bold{P} 及其逆矩阵 P1\bold{P}^{-1} ,使得矩阵A\bold{A} 对角化:A=P[λ1000λ2000λn]P1\bold{A} = \bold{P}\begin{bmatrix}\lambda_1&0&\cdots&0\\ 0&\lambda_2&\cdots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&\lambda_n\end{bmatrix}\bold{P}^{-1},则有

eAt=P[eλ1t000eλ2t000eλnt]P1(3)e^{\bold{A}t}=\bold{P}\begin{bmatrix} e^{\lambda_1t}&0&\cdots&0\\ 0&e^{\lambda_2t}&\cdots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&e^{\lambda_nt}\end{bmatrix}\bold{P}^{-1} \tag{3}

解题步骤
  1. 求解系统矩阵A\bold{A} 的特征值 λ1,λ2λn\lambda_1,\lambda_2\dots\lambda_n 。(特征值两两互异)
  2. 求解特征值对应的特征向量p1,p2pnp_1,p_2\dots p_n,构造变换阵 P\bold{P} 并求解其逆矩阵 P1\bold{P}^{-1}
  3. 求解矩阵指数 eAt=P[eλ1t000eλ2t000eλnt]P1e^{\bold{A}t}=\bold{P}\begin{bmatrix} e^{\lambda_1t}&0&\cdots&0\\ 0&e^{\lambda_2t}&\cdots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\cdots&e^{\lambda_nt}\end{bmatrix}\bold{P}^{-1}
例题

试用化为对角标准型法求解矩阵A=[0123]\bold{A}=\begin{bmatrix}0&1\\-2&-3\end{bmatrix} 的矩阵指数 eAte^{\bold{A}t}

求解特征值λIA=λ12λ+3=(λ+1)(λ+2)|\lambda\bold{I}-\bold{A}|=\begin{vmatrix}\lambda&-1\\2&\lambda+3\end{vmatrix}=(\lambda+1)(\lambda+2),得到特征值为λ1=1\lambda_1=-1λ2=2\lambda_2=-2。继而求解特征向量p1=[11]p_1=\begin{bmatrix}1\\-1\end{bmatrix}p2=[12]p_2=\begin{bmatrix}1\\-2\end{bmatrix}

故变换矩阵 P=[1112]\bold{P}=\begin{bmatrix}1&1\\-1&-2\end{bmatrix},求逆有 P1=[2111]\bold{P}^{-1}=\begin{bmatrix}2&1\\-1&-1\end{bmatrix}

则矩阵指数为 eAt=P[et00e2t]P1=[2ete2tete2t2et+2e2tet+2e2t]e^{\bold{A}t}=\bold{P}\begin{bmatrix} e^{-t}&0\\ 0&e^{-2t}\end{bmatrix}\bold{P}^{-1}=\begin{bmatrix} 2e^{-t}-e^{-2t}&e^{-t}-e^{-2t}\\ -2e^{-t}+2e^{-2t}&-e^{-t}+2e^{-2t}\end{bmatrix}

试用化为对角标准型法求解矩阵A=[01161166115]\bold{A}=\begin{bmatrix}0&1&-1\\-6&-11&6\\-6&-11&5\end{bmatrix} 的矩阵指数 eAte^{\bold{A}t}

求解特征值λIA=λ116λ+116611λ5=(λ+1)(λ+2)(λ+3)|\lambda\bold{I}-\bold{A}|=\begin{vmatrix}\lambda&-1&1\\6&\lambda+11&-6\\6&11&\lambda-5\end{vmatrix}=(\lambda+1)(\lambda+2)(\lambda+3),得到特征值为λ1=1\lambda_1=-1λ2=2\lambda_2=-2λ3=3\lambda_3=-3。继而求解特征向量p1=[101]p_1=\begin{bmatrix}1\\0\\1\end{bmatrix}p2=[124]p_2=\begin{bmatrix}1\\2\\4\end{bmatrix}p3=[169]p_3=\begin{bmatrix}1\\6\\9\end{bmatrix}

故变换矩阵 P=[111026149]\bold{P}=\begin{bmatrix}1&1&1\\0&2&6\\1&4&9\end{bmatrix},求逆有 P1=[35223431321]\bold{P}^{-1}=\begin{bmatrix}3&\frac{5}{2}&-2\\-3&-4&3\\1&\frac{3}{2}&-1\end{bmatrix}

则矩阵指数为

eAt=P[et000e2t000e3t]P1=[111026149][et000e2t000e3t][35223431321]=[3et3e2t+e3t52et4e2t+32e3t2et+3e2te3t6et+6e3t8e2t+9e3t6e2t6e3t3et12e2t+9e3t52et16e2t+272e3t2et+12e2t9e3t]\begin{aligned}e^{\bold{A}t}&=\bold{P}\begin{bmatrix} e^{-t}&0&0\\ 0&e^{-2t}&0\\0&0&e^{-3t}\end{bmatrix}\bold{P}^{-1}=\begin{bmatrix}1&1&1\\0&2&6\\1&4&9\end{bmatrix}\begin{bmatrix} e^{-t}&0&0\\ 0&e^{-2t}&0\\0&0&e^{-3t}\end{bmatrix}\begin{bmatrix}3&\frac{5}{2}&-2\\-3&-4&3\\1&\frac{3}{2}&-1\end{bmatrix}\\&=\begin{bmatrix} 3e^{-t}-3e^{-2t}+e^{-3t}&\frac{5}{2}e^{-t}-4e^{-2t}+\frac{3}{2}e^{-3t}&-2e^{-t}+3e^{-2t}-e^{-3t}\\ -6e^{-t}+6e^{-3t}&-8e^{-2t}+9e^{-3t}&6e^{-2t}-6e^{-3t}\\3e^{-t}-12e^{-2t}+9e^{-3t}&\frac{5}{2}e^{-t}-16e^{-2t}+\frac{27}{2}e^{-3t}&-2e^{-t}+12e^{-2t}-9e^{-3t}\end{bmatrix}\end{aligned}

(2)当 n×nn\times{n} 矩阵A\bold{A}nn 重特征根时,存在线性非奇异变换 P\bold{P} 及其逆矩阵 P1\bold{P}^{-1} ,将矩阵 A\bold{A} 转化为若尔当标准型:A=P[λ100λ0100λ]n×nP1\bold{A} = \bold{P}\begin{bmatrix}\lambda&1&\cdots&0\\ 0&\lambda&\cdots&0\\ \vdots&\vdots&\ddots&1\\ 0&0&\cdots&\lambda\end{bmatrix}_{n\times{n}}\bold{P}^{-1},则有

eAt=Peλt[1tt22!tn1(n1)!01ttn2(n2)!000t0001]n×nP1(4)e^{\bold{A}t}=\bold{P}e^{\lambda t}\begin{bmatrix} 1&t&\frac{t^2}{2!}&\cdots&\frac{t^{n-1}}{(n-1)!}\\ 0&1&t&\cdots&\frac{t^{n-2}}{(n-2)!}\\ \vdots&\vdots&\vdots&\ddots&\vdots\\0&0&0&\cdots&t\\ 0&0&0&\cdots&1\end{bmatrix}_{n\times{n}}\bold{P}^{-1} \tag{4}

拓展到一般情况,矩阵A\bold{A} 同时存在重特征根和单特征根时,以有三重根λ1\lambda_1、两重根λ2\lambda_2 和单根λ3\lambda_3 的矩阵A\bold{A} 为例,若存在变换阵 P\bold{P} 及其逆矩阵 P1\bold{P}^{-1} ,将矩阵 A\bold{A} 转化为若尔当标准型:A=P[λ110λ11λ1λ21λ20λ1]P1\bold{A} = \bold{P}\begin{bmatrix}\lambda_1&1&&&&\bold{0}\\ &\lambda_1&1&&&\\ &&\lambda_1&&&\\&&&\lambda_2&1&\\&&&&\lambda_2&\\\bold{0}&&&&&\lambda_1\end{bmatrix}\bold{P}^{-1},则有

eAt=P[eλ1tteλ1t12t2eλ1t0000eλ1tteλ1t00000eλ1t000000eλ2tteλ2t00000eλ2t000000eλ3t]P1(5)e^{\bold{A}t}=\bold{P}\begin{bmatrix}e^{\lambda_1t}&te^{\lambda_1t}&\frac{1}{2}t^2e^{\lambda_1t}&0&0&0\\ 0&e^{\lambda_1t}&te^{\lambda_1t}&0&0&0\\ 0&0&e^{\lambda_1t}&0&0&0\\0&0&0&e^{\lambda_2t}&te^{\lambda_2t}&0\\0&0&0&0&e^{\lambda_2t}&0\\0&0&0&0&0&e^{\lambda_3t}\end{bmatrix}\bold{P}^{-1} \tag{5}

例题

试求矩阵A=[065102324]\bold{A}=\begin{bmatrix}0&6&-5\\1&0&2\\3&2&4\end{bmatrix} 的矩阵指数。

求解特征值λIA=λ651λ232λ4=(λ1)2(λ2)|\lambda\bold{I}-\bold{A}|=\begin{vmatrix}\lambda&-6&5\\-1&\lambda&-2\\-3&-2&\lambda-4\end{vmatrix}=(\lambda-1)^2(\lambda-2),得到特征值为λ1=λ2=1\lambda_1=\lambda_2=1λ3=2\lambda_3=2。继而求解特征向量和广义特征向量p1=[13757]p_1=\begin{bmatrix}1\\-\frac{3}{7}\\-\frac{5}{7}\end{bmatrix}p2=[122494649]p_2=\begin{bmatrix}1\\-\frac{22}{49}\\-\frac{46}{49}\end{bmatrix}p3=[212]p_3=\begin{bmatrix}2\\-1\\-2\end{bmatrix}

故变换矩阵 P=[11237224915746492]\bold{P}=\begin{bmatrix}1&1&2\\-\frac{3}{7}&-\frac{22}{49}&-1\\-\frac{5}{7}&-\frac{46}{49}&-2\end{bmatrix},求逆有 P1=[26572874111]\bold{P}^{-1}=\begin{bmatrix}2&-6&5\\7&28&-7\\-4&-11&1\end{bmatrix}

则矩阵指数为

eAt=P[ettet00et000e2t]P1=[9et+7tet8e2t22et+28tet+22e2t2et7tet+2e2t4et3tet+4e2t10et12tet+11e2tet+3tete2t8et5tet+8e2t22et20tet22e2t3et+5tet2e2t]\begin{aligned}e^{\bold{A}t}&=\bold{P}\begin{bmatrix} e^{-t}&te^{-t}&0\\ 0&e^{t}&0\\0&0&e^{2t}\end{bmatrix}\bold{P}^{-1}\\&=\begin{bmatrix} 9e^{t}+7te^{t}-8e^{2t}&22e^{t}+28te^{t}+-22e^{2t}&-2e^{t}-7te^{t}+2e^{2t}\\ -4e^{t}-3te^{t}+4e^{2t}&-10e^{t}-12te^{t}+11e^{2t}&e^{t}+3te^{t}-e^{2t}\\-8e^{t}-5te^{t}+8e^{2t}&-22e^{t}-20te^{t}-22e^{2t}&3e^{t}+5te^{t}-2e^{2t}\end{bmatrix}\end{aligned}

  1. 化矩阵指数为矩阵A\bold{A} 的有限项。

该方法将矩阵指数表示为eAt=a0(t)I+a1(t)A++an1An1e^{\bold{A}t}=a_0(t)\bold{I}+a_1(t)\bold{A}+\cdots+a_{n-1}\bold{A}^{n-1}

当特征值两两互异时,

[a0(t)a1(t)an1(t)]=[1λ1λ1n11λ2λ2n11λnλnn1]1[eλ1teλ2teλnt](6)\begin{bmatrix}a_0(t)\\a_1(t)\\\vdots\\a_{n-1}(t)\end{bmatrix}=\begin{bmatrix}1&\lambda_1&\cdots&\lambda_1^{n-1}\\1&\lambda_2&\cdots&\lambda_2^{n-1}\\\vdots&\vdots&\ddots&\vdots\\1&\lambda_n&\cdots&\lambda_n^{n-1}\end{bmatrix}^{-1}\begin{bmatrix}e^{\lambda_1t}\\e^{\lambda_2t}\\\vdots\\e^{\lambda_nt}\end{bmatrix} \tag{6}

当存在重特征值时(以三重根λ1\lambda_1 和二重根λ2\lambda_2,其余根为单根为例),

[a0(t)a1(t)a2(t)a3(t)a4(t)a5(t)an1(t)]=[0013λ1(n1)(n2)2!λ1n3012λ13λ12(n1)1!λ1n21λ1λ12λ13λ1n1012λ23λ22(n1)1!λ2n21λ2λ22λ23λ2n11λ3λ32λ33λ3n11λnλn2λn3λnn1]1[12!t2eλ1t11!teλ1teλ1t11!teλ2teλ2teλ3teλn3t](7)\begin{bmatrix}a_0(t)\\a_1(t)\\a_2(t)\\a_3(t)\\a_4(t)\\a_5(t)\\\vdots\\a_{n-1}(t)\end{bmatrix}=\begin{bmatrix}0&0&1&3\lambda_1&\cdots&\frac{(n-1)(n-2)}{2!}\lambda_1^{n-3}\\0&1&2\lambda_1&3\lambda_1^2&\cdots&\frac{(n-1)}{1!}\lambda_1^{n-2}\\1&\lambda_1&\lambda_1^2&\lambda_1^3&\cdots&\lambda_1^{n-1}\\0&1&2\lambda_2&3\lambda_2^2&\cdots&\frac{(n-1)}{1!}\lambda_2^{n-2}\\1&\lambda_2&\lambda_2^2&\lambda_2^3&\cdots&\lambda_2^{n-1}\\1&\lambda_3&\lambda_3^2&\lambda_3^3&\cdots&\lambda_3^{n-1}\\\vdots&\vdots&\vdots&\vdots&\ddots&\vdots\\1&\lambda_n&\lambda_n^2&\lambda_n^3&\cdots&\lambda_n^{n-1}\end{bmatrix}^{-1}\begin{bmatrix}\frac{1}{2!}t^2e^{\lambda_1t}\\\frac{1}{1!}te^{\lambda_1t}\\e^{\lambda_1t}\\\frac{1}{1!}te^{\lambda_2t}\\e^{\lambda_2t}\\e^{\lambda_3t}\\\vdots\\e^{\lambda_{n-3}t}\end{bmatrix} \tag{7}

证明:Cayley-Hamilton定理
解题步骤
  1. 求解系统矩阵A\bold{A} 的特征值 λ1,λ2λn\lambda_1,\lambda_2\dots\lambda_n

  2. 求解有限项,根据特征值的互异性分情况分析:

    • 当特征值两两互异时,直接根据[a0(t)a1(t)an1(t)]=[1λ1λ1n11λ2λ2n11λnλnn1]1[eλ1teλ2teλnt]\begin{bmatrix}a_0(t)\\a_1(t)\\\vdots\\a_{n-1}(t)\end{bmatrix}=\begin{bmatrix}1&\lambda_1&\cdots&\lambda_1^{n-1}\\1&\lambda_2&\cdots&\lambda_2^{n-1}\\\vdots&\vdots&\ddots&\vdots\\1&\lambda_n&\cdots&\lambda_n^{n-1}\end{bmatrix}^{-1}\begin{bmatrix}e^{\lambda_1t}\\e^{\lambda_2t}\\\vdots\\e^{\lambda_nt}\end{bmatrix} 求解有限项。
    • 当特征值存在重根时,对于单根部分列写方程:

eλit=a0(t)+a1(t)λi++an1(t)λin1e^{\lambda_it}=a_0(t)+a_1(t)\lambda_i+\cdots+a_{n-1}(t)\lambda_i^{n-1}

而对于kk 重根部分在列写方程eλit=a0(t)+a1(t)λi++ak1(t)λik1e^{\lambda_it}=a_0(t)+a_1(t)\lambda_i+\cdots+a_{k-1}(t)\lambda_i^{k-1} 外还需要补充方程:

{teλit=a1(t)+2a2(t)λi++(k1)ak1(t)λik2t2eλit=2a2(t)+6a3(t)λi++(k1)(k2)ak1(t)λik3tk1eλit=(k1)!ak1(t)\left\{ \begin{matrix} te^{\lambda_it}=a_1(t)+2a_2(t)\lambda_i+\cdots+(k-1)a_{k-1}(t)\lambda_i^{k-2}\\t^2e^{\lambda_it}=2a_2(t)+6a_3(t)\lambda_i+\cdots+(k-1)(k-2)a_{k-1}(t)\lambda_i^{k-3} \\\vdots\\t^{k-1}e^{\lambda_it}=(k-1)!a_{k-1}(t) \\\end{matrix}\right.

联立nn 条方程求解有限项

  1. 代入求解矩阵指数:

eAt=a0(t)I+a1(t)A++an1An1e^{\bold{A}t}=a_0(t)\bold{I}+a_1(t)\bold{A}+\cdots+a_{n-1}\bold{A}^{n-1}

例题

试求矩阵A=[010001230]\bold{A}=\begin{bmatrix}0&1&0\\0&0&1\\2&3&0\end{bmatrix} 的矩阵指数,利用化为有限项法求解。

求解特征值λIA=λ100λ123λ=(λ+1)2(λ2)|\lambda\bold{I}-\bold{A}|=\begin{vmatrix}\lambda&-1&0\\0&\lambda&-1\\-2&-3&\lambda\end{vmatrix}=(\lambda+1)^2(\lambda-2),得到特征值为λ1,2=1\lambda_{1,2}=-1λ3=2\lambda_3=2

对于单根 λ3=2\lambda_3=2,有 e2t=a0(t)+2a1(t)+4a2(t)e^{2t}=a_0(t)+2a_1(t)+4a_2(t)

对于二重根λ1,2=1\lambda_{1,2}=-1,有 et=a0(t)a1(t)+a2(t)e^{-t}=a_0(t)-a_1(t)+a_2(t),还需要补充方程:

tet=a1(t)2a2(t)te^{-t}=a_1(t)-2a_2(t)

联立三组方程解得:

{a0(t)=19(e2t+8et+6tet)a1(t)=19(2e2t2et+3tet)a3(t)=19(e2tet3tet)\left\{ \begin{matrix} a_0(t)=\frac{1}{9}(e^{2t}+8e^{-t}+6te^{-t})\\ a_1(t)=\frac{1}{9}(2e^{2t}-2e^{-t}+3te^{-t}) \\ a_3(t)=\frac{1}{9}(e^{2t}-e^{-t}-3te^{-t}) \\\end{matrix}\right.

eAt=a0(t)I+a1(t)A++an1An1=19[e2t+(8+6t)ete2t(23t)ete2t(1+3t)et2e2t(2+6t)et4e2t+(53t)et2e2t(23t)et4e2t+(64t)et8e2t+(38t)et4e2t+(53t)et]\begin{aligned}e^{\bold{A}t}&=a_0(t)\bold{I}+a_1(t)\bold{A}+\cdots+a_{n-1}\bold{A}^{n-1}\\&=\frac{1}{9}\begin{bmatrix} e^{2t}+(8+6t)e^{-t}&e^{2t}-(2-3t)e^{-t}&e^{2t}-(1+3t)e^{-t}\\ 2e^{2t}-(2+6t)e^{-t}&4e^{2t}+(5-3t)e^{-t}&2e^{2t}-(2-3t)e^{-t}\\4e^{2t}+(6-4t)e^{-t}&8e^{2t}+(3-8t)e^{-t}&4e^{2t}+(5-3t)e^{-t}\end{bmatrix}\end{aligned}

# 线性时不变系统非齐次状态方程的解

动态系统在控制的作用下的运动称为受控运动。线性时不变系统非齐次状态方程的解即为线性时不变系统的受控运动。考虑系统 x˙(t)=Ax(t)+Bu(t),x(0),t0\dot{x}(t)=Ax(t)+Bu(t),x(0),t\geq0,其动态响应形式为:

x(t)=eA(tt0)x(t0)+t0teA(tτ)Bu(τ)dτ,t0(8)x(t)=e^{A(t-t_0)}x(t_0)+\int_{t_0}^te^{A(t-\tau)}Bu(\tau)\,d\tau, t\geq0 \tag{8}

可理解为由两部分组成:一部分是由初始状态引起的系统自由运动,即零输入响应;另外一部分是由控制输入所产生的受控运动,即零状态响应。

推导过程

对于系统x˙(t)=Ax(t)+Bu(t),x(0),t0\dot{x}(t)=Ax(t)+Bu(t),x(0),t\geq0,左乘eAte^{-At} 后求导可得:

ddt[eAtx(t)]=eAt[x˙(t)Ax(t)]=eAtBu(t)\frac{d}{dt}[e^{-At}x(t)]=e^{-At}[\dot{x}(t)-Ax(t)]=e^{-At}Bu(t)

两边积分得:

0t{ddt[eAtx(t)]}dτ=0teAtBu(t)dτ\int_0^t\{\frac{d}{dt}[e^{-At}x(t)]\}d\tau=\int_0^te^{-At}Bu(t)d\tau

eAtx(t)x(0)I=0teAtBu(t)dτe^{-At}x(t)-x(0)I=\int_0^te^{-At}Bu(t)d\tau

x(t)=eA(tt0)x(t0)+t0teA(tτ)Bu(τ)dτ,t0x(t)=e^{A(t-t_0)}x(t_0)+\int_{t_0}^te^{A(t-\tau)}Bu(\tau)\,d\tau, t\geq0

# 线性时不变系统的状态转移矩阵

在线性时不变系统解 x(t)=eA(tt0)x(t0)+t0teA(tτ)Bu(τ)dτ,t0x(t)=e^{A(t-t_0)}x(t_0)+\int_{t_0}^te^{A(t-\tau)}Bu(\tau)\,d\tau, t\geq0 中,定义状态转移矩阵Φ(t,t0)=eA(tt0)\Phi(t,t_0)=e^{A(t-t_0)}

  1. 线性时不变系统的状态转移矩阵可记为:Φ(t,t0)=Φ(tt0)\Phi(t,t_0)=\Phi{(t-t_0)}
  2. x(t)x(t) 是由初始值引起的零输入解和控制产生的零状态解的叠加。
  3. 解的结构显示了从x(t0)x(t_0)x(t)x(t) 的一种变换关系。
线性连续系统的状态转移矩阵
  1. 定义

对于线性连续系统的状态方程:x˙(t)=A(t)x(t)+B(t)u(t),x(t0)=x0,A(t)Rn×n\dot{x}(t)=A(t)x(t)+B(t)u(t),x(t_0)=x_0,A(t)\in{R^{n\times{n}}},那么称满足以下矩阵方程的解Φ(t,t0)\Phi(t,t_0) 为系统的状态转移矩阵。

Φ˙(t,t0)=A(t)Φ(t,t0),Φ(t0,t0)=I,tt0(9)\dot{\Phi}(t,t_0)=A(t)\Phi(t,t_0),\Phi(t_0,t_0)=I,t\geq{t_0} \tag{9}

  1. 状态转移矩阵的性质
  • dΦ(t,t0)dt=A(t)Φ(t,t0),Φ(t0,t0)=I\frac{d\Phi(t,t_0)}{dt}=A(t)\Phi(t,t_0),\Phi(t_0,t_0)=I
  • Φ(t2,t1)Φ(t1,t0)=Φ(t2,t0)\Phi(t_2,t_1)\Phi(t_1,t_0)=\Phi(t_2,t_0)
  • Φ(mt)=Φ(t+t++t)=[Φ(t)]m\Phi(mt)=\Phi(t+t+\cdots+t)=[\Phi(t)]^m

# 线性时变系统状态方程的解 *

# 线性时变系统齐次状态方程的解

# 线性时变系统的状态转移矩阵

# 线性时变系统非齐次状态方程的解

# 线性连续系统的时间离散化

线性连续系统的时间离散化问题本质上就是在一定的采样方式和保持方式下,由系统的连续时间状态空间描述来得到对应的离散时间状态空间描述,并建立两者的系数矩阵间的关系式。

# 近似离散化

考虑以下线性时变系统:x˙(t)=A(t)x(t)+B(t)u(t)\dot{x}(t)=A(t)x(t)+B(t)u(t),当采样周期TT 较小且精度要求不高时,可将其离散化为:

x˙(kT)1T[x((k+1)T)x(kT)](10)\dot{x}(kT)\approx \frac{1}{T}[x((k+1)T)-x(kT)] \tag{10}

t=kTt=kT,有

1T[x((k+1)T)x(kT)]=A(kT)x(kT)+B(kT)u(kT)\frac{1}{T}[x((k+1)T)-x(kT)]=A(kT)x(kT)+B(kT)u(kT)

x[(k+1)T]=[I+TA(kT)]x(kT)+TB(kT)u(kT)=G(kT)x(kT)+H(kT)u(kT)\begin{aligned}x[(k+1)T]&=[I+TA(kT)]x(kT)+TB(kT)u(kT)\\&=G(kT)x(kT)+H(kT)u(kT)\end{aligned}

其中,G(kT)=I+TA(kT)G(kT)=I+TA(kT)H(kT)=TB(kT)H(kT)=TB(kT)

注:一般而言,当采样周期为系统最小时间系数的110\frac{1}{10} 左右,近似度已经足够。

例题

系统的状态方程为x˙(t)=A(t)x(t)+B(t)u(t)\dot{x}(t)=A(t)x(t)+B(t)u(t),其中A(t)=[05(1e5t)05(e5t1)]A(t)=\begin{bmatrix}0&5(1-e^{-5t})\\0&5(e^{-5t}-1)\end{bmatrix}B(t)=[55e5t05(1e5t)]B(t)=\begin{bmatrix}5&5e^{-5t}\\0&5(1-e^{-5t})\end{bmatrix}。试求采样周期为T=0.2sT=0.2s 时的离散状态方程。

直接代入公式有:

G(kT)=I+TA(kT)=[11ek0ek]G(kT)=I+TA(kT)=\begin{bmatrix}1&1-e^{-k}\\0&e^{-k}\end{bmatrix}

H(kT)=TB(kT)=[1ek01ek]H(kT)=TB(kT)=\begin{bmatrix}1&e^{-k}\\0&1-e^{-k}\end{bmatrix}

那么,离散状态方程为:x[(k+1)T]=G(kT)x(kT)+H(kT)u(kT)x[(k+1)T]=G(kT)x(kT)+H(kT)u(kT)

将状态方程x˙=[0123]x+[01]u\dot{x}=\begin{bmatrix}0&1\\-2&-3\end{bmatrix}x+\begin{bmatrix}0\\1\end{bmatrix}u 近似离散化,T=0.2sT=0.2s

由题:G=I+TA=[1001]+0.2[0123]=[10.20.40.4]G=I+TA=\begin{bmatrix}1&0\\0&1\end{bmatrix}+0.2\begin{bmatrix}0&1\\-2&-3\end{bmatrix}=\begin{bmatrix}1&0.2\\-0.4&0.4\end{bmatrix}H=0.2[01]=[00.2]H=0.2\begin{bmatrix}0\\1\end{bmatrix}=\begin{bmatrix}0\\0.2\end{bmatrix}

故离散状态方程为:

$x[0.2(k+1)]=\begin{bmatrix}1&0.2\\-0.4&0.4\end{bmatrix}x(0.2k)+\begin{bmatrix}0\\0.2\end{bmatrix}u(0.2k)$

# 线性时不变系统状态方程的离散化

在线性时不变系统中,x˙(t)=A(x)+B(u)\dot{x}(t)=A(x)+B(u),其时间离散化状态方程为:

x[(k+1)T]=Gx(kT)+Hu(kT)(11)x[(k+1)T]=Gx(kT)+Hu(kT) \tag{11}

其中G=eATG=e^{AT}H=(0TeATdt)BH=(\int_0^Te^{AT}dt)B。假设条件为:(1) 等采样周期TT;(2)u(t)u(kT),kTt(k+1)Tu(t)\equiv u(kT),kT\leq t\leq (k+1)T

推导证明

对于线性时不变系统 x˙(t)=A(x)+B(u)\dot{x}(t)=A(x)+B(u),其状态方程的解为:

x(t)=eA(tt0)x(t0)+eA(tτ)Bu(τ)dτ(12)x(t)=e^{A(t-t_0)}x(t_0)+\int e^{A(t-\tau)}Bu(\tau)d\tau \tag{12}

假设:(1) 等采样周期TT;(2)x(k)=[x(t)]t=kTx(k)=[x(t)]_{t=kT},u(k)=[u(t)]_

那么令 t=(k+1)Tt=(k+1)Tt0=kTt_0=kT,有:

x[(k+1)T]=eATx(kT)+kT(k+1)TeA[(k+1)Tτ]Bu(τ)dτ=eATx(kT)+kT(k+1)TeA[(k+1)Tτ]Bdτu(kT)\begin{aligned}x[(k+1)T]&=e^{AT}x(kT)+\int_{kT}^{(k+1)T}e^{A[(k+1)T-\tau]}Bu(\tau)d\tau\\&=e^{AT}x(kT)+\int_{kT}^{(k+1)T}e^{A[(k+1)T-\tau]}Bd\tau \cdot u(kT)\end{aligned}

t=(k+1)Tτt=(k+1)T-\taudτ=dtd\tau =-dt,有:

x[(k+1)T]=eATx(kT)+0τeA(t)Bdtu(kT)=eATx(kT)+0τeA(t)dtBu(kT)\begin{aligned}x[(k+1)T]&=e^{AT}x(kT)+\int_{0}^{\tau}e^{A(t)}Bdt\cdot u(kT)\\&=e^{AT}x(kT)+\int_{0}^{\tau}e^{A(t)}dt\cdot Bu(kT)\end{aligned}

G=eATG=e^{AT}H=(0TeATdt)BH=(\int_0^Te^{AT}dt)B,有线性时不变系统的离散状态方程为:

x[(k+1)T]=Gx(kT)+Hu(kT)x[(k+1)T]=Gx(kT)+Hu(kT)

解题步骤
  1. 求解矩阵指数,方法见矩阵指数的计算
  2. 求解系数矩阵:G=eATG=e^{AT}H=(0TeATdt)BH=(\int_0^Te^{AT}dt)B
  3. 列写时间离散化状态方程:x[(k+1)T]=Gx(kT)+Hu(kT)x[(k+1)T]=Gx(kT)+Hu(kT)
例题

将状态方程x˙=[0102]x+[01]u\dot{x}=\begin{bmatrix}0&1\\0&-2\end{bmatrix}x+\begin{bmatrix}0\\1\end{bmatrix}u 离散化,T=0.1sT=0.1s

利用拉氏变换法求解矩阵指数函数。取拉氏变换有:

[sIA]1=[s10s+2]1=[1s1s(s+2)01s+2][sI-A]^{-1}=\begin{bmatrix}s&-1\\0&s+2\end{bmatrix}^{-1}=\begin{bmatrix}\frac{1}{s}&\frac{1}{s(s+2)}\\0&\frac{1}{s+2}\end{bmatrix}

取拉氏逆变换得到矩阵指数函数:

eAt=L1[sIA]1=[10.5(1e2T)0e2T]e^{At}=L^{-1}[sI-A]^{-1}=\begin{bmatrix}1&0.5(1-e^{-2T})\\0&e^{-2T}\end{bmatrix}

进而求解系数矩阵:

G=eAT=[10.5(1e2T)0e2T]=[10.09100.819]G=e^{AT}=\begin{bmatrix}1&0.5(1-e^{-2T})\\0&e^{-2T}\end{bmatrix}=\begin{bmatrix}1&0.091\\0&0.819\end{bmatrix}

H=(0TeATdt)B=[0T[10.5(1e2T)0e2T]dt][01]=[T0.5T+0.25e2T0.2500.5e2T+0.5][01]=[0.0050.091]\begin{aligned}H&=(\int_0^Te^{AT}dt)B=\Bigg[\int_0^T\begin{bmatrix}1&0.5(1-e^{-2T})\\0&e^{-2T}\end{bmatrix}dt\Bigg]\begin{bmatrix}0\\1\end{bmatrix}\\&=\begin{bmatrix}T&0.5T+0.25e^{-2T}-0.25\\0&-0.5e^{-2T}+0.5\end{bmatrix}\begin{bmatrix}0\\1\end{bmatrix}=\begin{bmatrix}0.005\\0.091\end{bmatrix}\end{aligned}

故时间离散化状态方程为:

x[0.1(k+1)]=[10.09100.819]x(0.1k)+[0.0050.091]u(0.1k)x[0.1(k+1)]=\begin{bmatrix}1&0.091\\0&0.819\end{bmatrix}x(0.1k)+\begin{bmatrix}0.005\\0.091\end{bmatrix}u(0.1k)

# 线性离散系统状态方程的解

离散系统的差分方程形状态方程有两种解法:递推法和 z 变换法。其中递推法在时变系统和时不变系统中都适用,而 z 变换法只适用于时不变系统。

# 递推法

  1. 在线性时变系统中,x(k+1)=G(k)x(k)+H(k)u(k)x(k+1)=G(k)x(k)+H(k)u(k),有:

{x(1)=G(0)x(0)+H(0)u(0)x(2)=G(1)x(1)+H(1)u(1)x(3)=G(2)x(2)+H(2)u(2)\left\{ \begin{matrix} x(1)=G(0)x(0)+H(0)u(0)\\ x(2)=G(1)x(1)+H(1)u(1) \\ x(3)=G(2)x(2)+H(2)u(2) \\\vdots \end{matrix}\right.

给定初始条件x(0)x(0) 和输入序列u(0),u(1),u(0),u(1),\cdots 后即可求解x(k)x(k)

  1. 在线性时不变系统中,x(k+1)=Gx(k)+Hu(k)x(k+1)=Gx(k)+Hu(k),其中G,HG,H 均为常数矩阵,因此:

x(k)=Gkx(0)+i=0k1Gk1iHu(i)(13)x(k)=G^kx(0)+\sum_{i=0}^{k-1}G^{k-1-i}Hu(i) \tag{13}

上式称为线性时不变离散系统的状态转移方程,其中Φ(k)=Gk\Phi(k)=G^k 称为线性时不变离散系统的状态转移矩阵。

状态转移矩阵的性质:

  1. Φ(k+1)=GΦk,Φ(0)=I\Phi(k+1)=G\Phi{k},\Phi(0)=I
  2. Φ(k2k0)=Φ(k2k1)Φ(k1k0)\Phi(k_2-k_0)=\Phi(k_2-k_1)\Phi(k_1-k_0)
  3. Φ1(k)=Φ(k)\Phi^{-1}(k)=\Phi(-k)

# z 变换法

考虑时不变离散系统:x(k+1)=Gx(k)+Hu(k)x(k+1)=Gx(k)+Hu(k),取 z 变换有:

zx(z)zx(0)=Gx(z)+Hu(z)zx(z)-zx(0)=Gx(z)+Hu(z)

z(z)=(zIG)1zx(0)+(zIG)1Hu(z)(14)z(z)=(zI-G)^{-1}zx(0)+(zI-G)^{-1}Hu(z) \tag{14}

取 z 逆变换有:

x(k)=z1[(zIG)1z]x(0)+z1[(zIG)1Hu(z)](15)x(k)=z^{-1}\Big[(zI-G)^{-1}z\Big]x(0)+z^{-1}\Big[(zI-G)^{-1}Hu(z)\Big] \tag{15}

对比公式(13)和公式(15),由解的唯一性可知,

z1[(zIG)1z]=Gk(16)z^{-1}\Big[(zI-G)^{-1}z\Big]=G^k \tag{16}

z1[(zIG)1Hu(z)]=i=0k1Gk1iHu(i)(17)z^{-1}\Big[(zI-G)^{-1}Hu(z)\Big]=\sum_{i=0}^{k-1}G^{k-1-i}Hu(i) \tag{17}

例题

考虑离散系统:x(k+1)=Gx(k)+Hu(k)x(k+1)=Gx(k)+Hu(k),其中G=[010.161]G=\begin{bmatrix}0&1\\-0.16&-1\end{bmatrix}H=[11]H=\begin{bmatrix}1\\1\end{bmatrix},初始条件为x(0)=[11]x(0)=\begin{bmatrix}1\\-1\end{bmatrix},试求当u(k)=1u(k)=1 时状态方程的解。

用 z 变换法求解,先计算(zIG)1(zI-G)^{-1},有

(zIG)1=[z10.16z+1]1=1(z+0.2)(z+0.8)[z+110.16z]=[43×1z+0.213×1z+0.853×1z+0.253×1z+0.80.83×1z+0.2+0.83×1z+0.813×1z+0.2+43×1z+0.8]\begin{aligned}(zI-G)^{-1}&=\begin{bmatrix}z&-1\\0.16&z+1\end{bmatrix}^{-1}=\frac{1}{(z+0.2)(z+0.8)}\begin{bmatrix}z+1&1\\-0.16&z\end{bmatrix}\\&=\begin{bmatrix}\frac{4}{3}\times \frac{1}{z+0.2}-\frac{1}{3}\times \frac{1}{z+0.8}&\frac{5}{3}\times \frac{1}{z+0.2}-\frac{5}{3}\times \frac{1}{z+0.8}\\-\frac{0.8}{3}\times \frac{1}{z+0.2}+\frac{0.8}{3}\times \frac{1}{z+0.8}&-\frac{1}{3}\times \frac{1}{z+0.2}+\frac{4}{3}\times \frac{1}{z+0.8}\end{bmatrix}\end{aligned}

由于u(k)=1u(k)=1,则u(z)=zz1u(z)=\frac{z}{z-1},故zx(0)+Hu(z)=[zz]+[zz1zz1]=[z2z1z2+2zz1]zx(0)+Hu(z)=\begin{bmatrix}z\\-z\end{bmatrix}+\begin{bmatrix}\frac{z}{z-1}\\\frac{z}{z-1}\end{bmatrix}=\begin{bmatrix}\frac{z^2}{z-1}\\\frac{-z^2+2z}{z-1}\end{bmatrix}

那么代入公式(15)有:

x(z)=(zIG)1[zx(0)+Hu(z)]=[176×zz+0.2+229×zz+0.8+2518×zz13.46×zz+0.217.69×zz+0.8+718×zz1]\begin{aligned}x(z)&=(zI-G)^{-1}[zx(0)+Hu(z)]\\&=\begin{bmatrix}-\frac{17}{6}\times \frac{z}{z+0.2}+\frac{22}{9}\times \frac{z}{z+0.8}+\frac{25}{18}\times \frac{z}{z-1}\\\frac{3.4}{6}\times \frac{z}{z+0.2}-\frac{17.6}{9}\times \frac{z}{z+0.8}+\frac{7}{18}\times \frac{z}{z-1}\end{bmatrix}\end{aligned}

求 z 逆变换有:

x(k)=[176(0.2)k+229(0.8)k+25183.46(0.2)k17.69(0.2)k+718]x(k)=\begin{bmatrix}-\frac{17}{6}(-0.2)^k+\frac{22}{9}(-0.8)^k+\frac{25}{18}\\\frac{3.4}{6}(-0.2)^k-\frac{17.6}{9}(-0.2)^k+\frac{7}{18}\end{bmatrix}

更新于 阅读次数