2022年Matlab程序解现代控制理论与工程中的状态方程.docx
精品学习资源用矩阵指数法解状态方程的 MATLAB函数 vslove1：函数 vslove1: 求解线性定常连续系统状态方程的解function Phit,PhitBu=vsolves1A,B,ut%vsolves1 求线性连续系统状态方程X =AX+Bu 的解%Phit,phitBu=vsolves1A,B,ut%A,B 系数矩阵%ut 把握输入，必需为时域信号的符号表达式，符号变量为t%Phit 输出 Phit%PhitBu 输出 phit-tao*B*utao在区间（ 0， t）的积分syms t tao% 定义符号变量 t,taoPhit=expmA*t ；% 求矩阵指数 expAt if B=0B=zerossizeA,l,l ；% 重构系数矩阵B endphi=subPhit, t,t-tao；% 求 expAt-taoPhitBu=intphi*B*ut,tao,0,t； %求 expAt-tao*B*utao在 0t 区间的积分用拉氏变换法解状态方程的 MATLAB函数 vslove2：函数 vslove2: 求解线性定常连续系统状态方程的解function sl_A,sl_ABu=vsolves1A,B,us%vsolves2 求线性连续系统状态方程X =AX+Bu 的解%sl_A,sl_ABu=vsolves1A,B,ut%A,B 系数矩阵%us 把握输入，必需为拉氏变换后的符号表达式，符号变量为s%sl_A 输出矩阵 sl-A-1 拉式反变换的结果%sl_ABu 输出 sl-A-1*B*us拉式反变换后的结果syms s%定义符号变量t,tao AA=s*eyesizeA-A；% 求 sI-AinvAA=invAA；% 求（ sI-A ）矩阵的逆 intAA tAA=ilaplaceintAA；% 求 intAA 的拉氏反变换sI_A=simplify ；% 简化拉式反变换的结果if B=0B=zerossizeA,l,l ；% 重构系数矩阵B endtAB=ilaplaceintAA*B*us；%求 intAA*B*us的拉氏反变换sI_ABu=simplifytAB；%化简拉式反变换的结果求解时变系统状态方程的 MATLAB函数 tslove：函数 tslove: 求解线性时变连续系统状态方程的解function Phi,PhiBu=tsolvesA,B,u,x,a,n%tsolves 求时变系统状态方程%Phi,phiBu=vsolves1A,B,u,x,a,n%A,B 时变系数矩阵欢迎下载精品学习资源%Phi 状态转移矩阵运算结果%PhiBu 受控解重量%u把握输入向量，时域形式%x 符号变量，指明矩阵A 中的时变参数，通常为时间t%a积分下限%n时变状态转移矩阵中运算重积分的最大项数，n=0 时无重积分项% n=1 时包含二重积分项， .Phi=transmtxA,x,a,n ；% 运算状态转移矩阵Phitao=subsPhi,x, tao；% 求 Phitao if B=0Btao=zerossizeA,l,l ；%求 Btao endutao=subsu,x,tao；%求 utaoPhiBu=simpleintPhitao*Btao*utao,tao,a,x ； %运算受控重量求解时变系统转移矩阵的 MATLAB函数 transmtx ：函数 transmtx: 求解线性时变系统状态转移矩阵function Phi=transmtxA,x,a,n%transmtx 运算时变系统状态转移矩阵%Phi=transmtxA,x,a,n%Phi 状态转移矩阵运算结果%A时变系数矩阵%x 符号变量，指明矩阵A 中的时变参数，通常为时间t%a积分下限%n时变状态转移矩阵中运算重积分的最大项数，n=0 时无重积分项% n=1 时包含二重积分项， .Phi=eyesizeA ；% 初始化 Phi=I for lop=0:nAA=A ；for i=1:lopif AA=0break；endAtemp=subsAA,x, taoi ； AA=simplifyA*intAtemp,tao,a,x； endif AA=0break；endAtemp=subsAA,x, taoi ；AA=simplifyA*intAtemp,tao,a,x；% 运算重积分 Phi=simplifyPhi+AA；% 修正 Phiend求解线性定常离散系统状态方程的MATLAB函数 disolve：欢迎下载精品学习资源函数 disolve: 求解线性定常离散系统状态方程的解function Ak,AkBu=disolveA,B,uz%disolve 求线性离散系统状态方程xk+1=Axk+Buk的解%Ak,AkBu=disolveA,B,uz%A,B 系数矩阵%uz 把握输入，必需为Z 变换后的符号表达式，符号变量为z%Ak 输出矩阵 zI-A-1zZ反变换后的结果%AkBu 输出矩阵 zI-A-1*B*uzZ反变换后的结果syms z% 定义符号变量z AA=z*eyesizeA-A；% 求 zI-AinvAA=invAA ； % 求（ zI-A ）矩阵的逆 intAA tAA=iztransintAA*z ； % 求 intAA*z 的 Z 反变换Ak=simpletAA ； % 简化 Z 反变换的结果if B=0B=zerossizeA,l,l ； % 重构系数矩阵 B endtAB=iztransintAA*B*uz ； % 求 intAA*B*uz 的 Z 反变换AkBu=simpletAB ； %化简 Z 反变换的结果求解线性时变离散系统状态方程的MATLAB函数 tdsolve：函数 tdsolve: 求解线性时变离散系统状态方程的解function xk=tsolveAk,Bk,uk,x0,kstart,kend%tdsolve 求线性时变离散系统状态方程xk+1=Akxk+Bkuk的解%xk=tsolveAk,Bk,uk,x0,kstart,kend%Ak,Bk系数矩阵%uk 把握输入，必需为时域符号表达式，符号变量为k%x0初始状态%kstart 初始时刻%kend终止时刻%xk 输出结果，矩阵每一列分别对应xk0+1,xk0+2.syms k% 定义符号变量 k if Bk=0Bk=zerossizeA,l,l ；% 重构系数矩阵B endxk= ；for kk=kstart+1: kend AA=eyesizek ；for i=kstart : kk-1% 运算 Ak-1Ak-2.Ak0+1Ak0A=subsAk, k ,i ； AA=A*AA；end AAB=eyesizeAk；BB=zerossizeBk ；for i=kk-1 :-1:kk+1%运算 Ak-1Ak-2.Aj+1Bjuj的累加和欢迎下载精品学习资源A=subsAk, k ,i ； AAB=AAB*A；B=subsBk, k,kk-1+i+kstart ； u=subsuk,k,kk-1+i+kstart ； BB=BB+AAB*B*u；endB=subsBk, k,kk-1 ；u=subsuk,k,kk-1 ； BB=BB+B*u ；xk=xk AA*x0+BB；%运算 xk end连续系统状态方程离散化的 MATLAB符号函数 sc2d：函数 sc2d: 线性连续系统状态方程的离散化function Ak,Bk=sc2dA,B%sc2d 离散化线性连续系统状态方程X =AX+Bu%sysd=sc2dA,B%A,B 连续系统的系数矩阵%Ak,Bk离散系统系数符号矩阵%离散状态方程为：xk+1=Ak*xk+Bk*uk%Ak,Bk 中变量 T 为采样周期syms t T% 定义符号变量 t T Phit=expmA*t ；% 求矩阵指数 expAt if B=0B=zerossizeA,l,l ；% 重构系数矩阵B endPhitB=intPhit*B,t,0,T ；%求 expAt*B在 0T 区间的积分Ak=simplesubsPhit, t,T ； Bk=simplePhitB ；线性时变系统离散化的 MATLAB函数 tc2d：函数 tc2d: 线性时变系统的离散化function Ak,Bk=tc2dA,B,x,n%tc2d 线性时变系统的离散化%Ak,Bk=tc2dA,B,x,n%A,B 连续系统的系数矩阵%Ak离散化后的系数矩阵Akt%Bk离散化后的系数矩阵Bkt%x符号变量，指明矩阵AB 中的时变参数，通常为时间t%n时变状态转移矩阵中运算重积分的最大项数，n=0 时无重积分项，%n=1 时包含二重积分项， .syms t TPhit=transmtxA,x,k*T,n；% 运算时变系统的状态转移矩阵Ak=simplifysubsPhi,x,k+1*T；% 运算离散化后的系数矩阵AkT欢迎下载精品学习资源Phitao=subsPhi,x, tao；% 求 Phitao if B=0Btao=zerossizeA,l,l ；elseBtao=subsB,x, tao；% 求 Btao endPhitB=simpleintPhitao*Btao,tao,k*T,x ；%运算受控重量 Bk=simplifysubsPhiB,x,k+1*T；% 运算离散化后的系数矩阵BkT定常系统可控规范 I 型变换函数 ccanonl：函数 ccanonl: 求线性定常系统的可控规范I 型形式function Abar,Bbar,Cbar,T=ccanonlA,B,C%ccanonl 求系统 X =AX+Bu ， y=Cx 的可控规范 I 型系数矩阵%Abar,Bbar,Cbar, 变换后的可控规范I 型系数矩阵%T相像变换矩阵n=lengthA ；Co=ctrbA,B ；if rankCo=n,% 判定系统可控性error 系统不行控！ ； endRs=sympolymtxA；% 运算 R 矩阵并转变为符号矩阵Rs 形式As=symA ；% 讲矩阵 A 转变为符号矩阵 AsBs=symB ；%讲矩阵 B 转变为符号矩阵 BsTs=Bs；for i=1:n-1,Ts=Asi*Bs Ts ；endTs=Ts*Rs；% 运算相像变换符号矩阵TsAbar=numericinvTs*As*Ts；%实现矩阵 A 的相像变换并转变为数值形式Bbar=numericinvTs*Bs；%实现矩阵 B 的相像变换并转变为数值形式Cbar=numericinvC*Ts ；%实现矩阵 C 的相像变换并转变为数值形式T=numericTc ；% 相像变换矩阵 T 转变为数值形式求系统相像变换的 R 矩阵函数 polymtx ：函数 polymtx: 求系统相像矩阵变换的R 矩阵function R=polymtxA%R=polymtxA A 须为方针%1000%an-1100% R= an-2an-100%. .%a2a310%a1a2an-1 1%其中 aii=1,.n-1 为矩阵 A 的特点多项式欢迎下载精品学习资源%|s*I-A|=sn+an-1*sn-1+.+a1*s+a0%的各项系数d=sizeA ；if lengthd=2,% 判定系统可控性error 错误：非二维矩阵！ ；endif d1=d2,% 判定系统可控性error 错误： A 非方针！ ；end n=d1 ；As=symA ；p=sym2polypolyAs；R= ；for i=1:nR=R p1:n ；p=0 p1:n ； end R=numericR ；线性定常系统可控分解函数 cdecomp：函数 cdecomp: 线性定常系统的可控性分解function Abar,Bbar,Cbar,P=cdecompA,B,C%cdecomp 可控性分解%Abar,Bbar,Cbar,P=cdecompA,B,C%如系统不完全可控，就存在相像变换矩阵P，使得%-1-1%Abar=P*A*P, Bbar=P*B, Cbar=C*P%其中%Abar=|Ac A12| ,Bbar=|Bc |,Cbar=|Cc Cuc|0 Auc|0|%Ac,Bc 构成系统的可控子空间As=symA ；% 转变为符号矩阵求解Bs=symB ；Cs=symC ；nA=sizeAs,l ；Ms=symctrbAs,Bs ；% 求可控性矩阵 M n=numericrankMs ；if n<nA,P= ；% 系统不行控，运算变换矩阵P i=1 ；while numeric rankP<n,%取 M 的 n 个线性无关列矢量至P 矩阵P=P Ms:,i ；if numericrankP<sizeP,2,P :,sizeP,2= ；欢迎下载精品学习资源end i=i+1 ；endE=symeyesizeA ；i=1 ；while numeric rankP<nA,%取某一单位向量至P 阵使 P 满秩P=P E:,i ；if numericrankP<sizeP,2,P :,sizeP,2= ；end i=i+1 ；endelseP=eyesizeA ；% 如系统可控，取 P=1 endAbar=numericinvP*A*P；% 转变为数值矩阵输出Bbar=numericinvP*B；Cbar=numericinvC*P；P=numericP ；线性定常系统可观分解函数 odecomp：函数 odecomp: 线性定常系统的可观性分解function Abar,Bbar,Cbar,P=odecompA,B,C%odecomp 可控性分解%Abar,Bbar,Cbar,P=odecompA,B,C%如系统不完全可观，就存在相像变换矩阵P，使得%-1-1%Abar=P*A*P, Bbar=P*B, Cbar=C*P%其中%Abar=|Ao0| ,Bbar=|Bo |,Cbar=|Co 0|A21 Ano|Bno|%Ao,Bo 构成系统的可观子空间%依据对偶原理，应用可控性分解函数cdecomp 实现可观性分解abar,bbar,cbar,P=cdecompA,B,C；Abar=abar ； Bbar=cbar ； Cbar=bbar ；欢迎下载
