线性定常连续系统状态方程的解

Ch.3线性系统时域分析

第1页本章介绍本章简介本章讨论线性系统运动分析。主要介绍连续系统与离散系统状态空间模型求解、状态转移矩阵性质和计算以及连续系统状态方程离散化。本章最终介绍基于Matlab状态空间模型求解与控制系统运动仿真问题程序设计与仿真计算。第2页概述概述建立了系统数学描述之后,接着而来是对系统作定量和定性分析。定量分析主要包含研究系统对给定输入信号响应问题,也就是对描述系统状态方程和输出方程求解问题。定性分析主要包含研究系统结构性质,如能控性、能观性、稳定性等。第3页概述本章先讨论用状态空间模型描述线性系统定量分析问题,即状态空间模型--状态方程和输出方程求解问题。依据常微分方程理论求解一个一阶定常线性微分方程组,通常是很轻易。可是求解一个时变一阶线性微分方程组却非易事。状态转移矩阵引入,从而使得定常系统和时变系统求解公式含有一个统一形式。为此,本章将重点讨论状态转移矩阵定义、性质和计算方法,并在此基础上导出状态方程求解公式。第4页概述本章讨论另一个中心问题是连续系统状态方程离散化,即建立连续系统离散系统状态方程。伴随计算机在控制系统分析、设计和实时控制中广泛应用,这个问题显得越来越主要。在离散系统状态方程建立基础上,本章也将讨论对应状态方程求解问题,并将导出在形式上与连续系统状态方程解一致离散系统状态方程解。第5页概述(4/4)本章需处理问题:线性定常连续系统状态方程解理论基本概念:状态转移矩阵状态转移矩阵和矩阵指数函数eAt性质和计算怎样将线性定常连续系统离散化线性定常离散系统状态方程解理论第6页线性定常连续系统状态方程解3.1线性定常连续系统状态方程解求解状态方程是进行动态系统分析与综合基础,是进行定量分析主要方法。本节讲授状态方程求解理论是建立在状态空间上,以矩阵代数运算来描述定系数常微分方程解理论。下面基于矩阵代数运算状态方程解理论中,引入了状态转移矩阵这一基本概念。该概念对我们深刻了解系统动态特征、状态变迁(动态演变)等都是非常有帮助,对该概念必须准确掌握和深入了解。第7页线性定常连续系统状态方程解本节需处理主要问题状态转移矩阵?矩阵指数函数？状态转移矩阵和矩阵指数函数性质齐次状态方程求解?非齐次状态方程求解?非齐次状态方程解各部分意义？输出方程解？重点喔!重点与难点喔!要了解喔!第8页线性定常连续系统状态方程解在讨论普通线性定常连续系统状态方程解之前,先讨论线性定常齐次状态方程解,以引入矩阵指数函数和状态转移矩阵等概念。所谓齐次状态方程就是指状态方程中不考虑输入项(u(t)=0)作用,满足方程解齐次性。研究齐次状态方程解就是研究系统本身在无外力作用下自由(自治)运动。所谓非齐次状态方程就是指状态方程中输入项作用,状态方程解对输入含有非齐次性。研究非齐次状态方程解就是研究系统在外力作用下强迫运动。第9页线性定常连续系统状态方程解下面,将依次分别讨论:齐次状态方程解线性定常连续系统状态转移矩阵线性定常连续系统非齐次状态方程解系统脉冲响应第10页线性定常齐次状态方程解3.1.1线性定常齐次状态方程解什么是微分方程齐次方程?齐次方程就是指满足解齐次性方程，即若x是方程解，则对任意非零实数a,ax亦是该方程解。所谓齐次状态方程，即为以下不考虑输入自治方程x’=Ax齐次状态方程满足初始状态解,也就是由初始时刻t0初始状态x(t0)所引发无输入强迫项(无外力)时自由运动。第11页线性定常齐次状态方程解对上述齐次状态方程,惯用常微分方程求解方法有级数展开法和拉氏变换法

2种。第12页级数展开法1.级数展开法在求解齐次状态方程式之前,首先观察标量常微分方程在初始时刻t0=0解。该方程中x(t)为标量变量,a为常数。由常微分方程理论知,该方程解连续可微。所以,该解经泰勒展开可表征为无穷级数,即有式中,qk(k=1,2,...)为待定级数展开系数。第13页级数展开法(2/12)将所设解代入该微分方程,可得假如所设解是方程真实解,则对任意t,上式均成立。所以,使t有相同幂次项各项系数相等,即可求得令x(t)解表示式中t=0,可确定q0=x(0)所以,x(t)解表示式可写为第14页级数展开法(3/12)上述求解标量微分方程级数展开法,可推广至求解向量状态方程解。为此,设其解为t向量幂级数,即x(t)=q0+q1t+q2t2+…+qktk+…

式中,qk(k=1,2,...)为待定级数展开系数向量。将所设解代入该向量状态方程x’=Ax,可得q1+2q2t+3q3t2+…+kqktk-1+…=A(q0+q1t+q2t2+…+qktk+…)假如所设解是方程真实解,则对任意t,上式均成立.所以,使t有相同幂次项各项系数相等,即可求得第15页级数展开法(4/12)若初始时刻t0=0,初始状态x(0)=x0,则可确定q0=x(0)=x0所以,

状态x(t)解可写为该方程右边括号里展开式是n×n维矩阵函数。因为它类似于标量指数函数无穷级数展开式,所以称为矩阵指数函数,且记为第16页级数展开法(5/12)利用矩阵指数函数符号,齐次状态方程解可写为：x(t)=eAtx0第17页拉氏变换法(1/12)2．拉氏变换法若将对标量函数拉氏变换定义扩展到向量函数和矩阵函数,定义对向量函数和矩阵函数拉氏变换为分别对该向量函数和矩阵函数各个元素求对应拉氏变换,那么可利用拉氏变换及拉氏反变换方法求解齐次状态方程解。对该齐次状态方程x’=Ax,设初始时刻t0=0且初始状态x(t)=x0,对方程两边取拉氏变换,可得sX(s)-x0=AX(s)于是可求得该齐次状态方程解x(t)拉氏变换为X(s)=(sI-A)-1x0第18页拉氏变换法(2/12)对上式取拉氏反变换,即得齐次状态方程解为x(t)=L-1[(sI-A)-1]x0下面讨论怎样求解拉氏反变换L-1[(sI-A)-1]。主要思想为将标量函数拉氏变换与反变换平行推广至矩阵函数中。对标量函数,我们有第19页拉氏变换法(3/12)将上述关系式推广到矩阵函数则有其中eAt称为时间t矩阵指数函数，并有第20页拉氏变换法所以,基于上述(sI-A)-1拉氏反变换,该齐次方程解为x(t)=L-1[(sI-A)-1]x0=eAt

x0上述拉氏反变换法求解结果与前面级数展开法求解结果一致。若初始时刻t0

0,对上述齐次状态方程解作坐标变换,则可得解另一个表述形式:状态方程解表示式说明了齐次状态方程解实质上是初始状态x(t0)从初始时刻t0到时刻t系统运动状态转移,其转移特征和时刻t状态完全由矩阵指数函数和初始状态x(t0)所决定。第21页拉氏变换法为讨论方便,引入能描述系统状态转移特征线性定常连续系统状态转移矩阵以下:(t)=eAt所以,有以下关系式x(t)=(t)x0=(t-t0)x(t0)由上述状态转移矩阵定义和齐次状态方程解,系统状态转移矩阵有以下关系(t)=L-1[(sI-A)-1]第22页拉氏变换法齐次状态方程解描述了线性定常连续系统自由运动。由解表示式能够看出,系统自由运动轨线是由从初始时刻初始状态到t时刻状态转移刻划,如图3-1所表示。图3-1状态转移特征第23页拉氏变换法当初始状态给定以后,系统状态转移特征就完全由状态转移矩阵所决定。所以,状态转移矩阵包含了系统自由运动全部信息。可见,状态转移矩阵计算是齐次状态方程求解关键。第24页拉氏变换法—例3-1例3-1

试求以下状态方程在初始状态x0下解解(1)

首先求出矩阵指数函数eAt,其计算过程为第25页拉氏变换法(9/12)—例3-1(3)

状态方程解为(2)

计算矩阵指数函数eAt。第26页线性定常连续系统状态转移矩阵3.1.2线性定常连续系统状态转移矩阵下面深入讨论前面引入状态转移矩阵,主要内容为:基本定义矩阵指数函数和状态转移矩阵性质第27页基本定义:状态转移矩阵定义1.基本定义定义3-1

对于线性定常连续系统x’=Ax,当初始时刻t0=0时,满足以下矩阵微分方程和初始条件:

’(t)=A

(t),

(t)|t=0=I

解

(t)为线性定常连续系统x’=Ax状态转移矩阵。这里定义状态转移矩阵与前面定义是一致。引入上述状态转移矩阵新定义,主要是为了使状态转移矩阵概念易于推广到时变系统、离散系统等,使得有可能对各种类型系统状态方程解作统一描述，更加好地刻划系统状态运动改变规律。第28页几类特殊形式状态转移矩阵当系统矩阵A为n×n维方阵时,状态转移矩阵Φ(t)亦为n×n维方阵,且其元素为时间t函数。下面讨论几个特殊形式系统矩阵A状态转移矩阵1)对角线矩阵。当A为以下对角线矩阵:A=diag{

1

2…

n}则状态转移矩阵为式中，diag{…}表示由括号内元素组成对角线矩阵。第29页几类特殊形式状态转移矩阵(2)块对角矩阵。当A为以下块对角矩阵:A=block-diag{A1

A2…Al}其中Ai为mi

mi维分块矩阵,则状态转移矩阵为式中，block-diag{…}表示由括号内各方块矩阵组成块对角矩阵。第30页几类特殊形式状态转移矩阵(3)约旦块矩阵。当Ai为特征值为

imi

mi维约旦块,则分块矩阵矩阵指数函数为对上述三种特殊形式矩阵状态转移矩阵和矩阵指数函数,可利用矩阵指数函数展开式证实。第31页状态转移矩阵性质2.矩阵指数函数和状态转移矩阵性质由矩阵指数函数展开式和状态转移矩阵定义,可证实矩阵指数函数和状态转移矩阵含有以下性质(Φ(t)为方阵A状态转移矩阵)1)

Φ(0)=eA0=I第32页状态转移矩阵性质2)

eA(t+s)=eAteAs，

Φ(t+s)=Φ(t)Φ(s)式中t和s为两个独立标量自变量证实由指数矩阵函数展开式,有3)[Φ(t2-t1)]-1=Φ(t1-t2)第33页状态转移矩阵性质4)对于n

n阶方阵A和B,下式仅当AB=BA时才成立e(A+B)t=eAteBt5)

6)[Φ(t)]n=Φ(nt)7)

Φ(t2-t1)Φ(t1-t0)=Φ(t2-t0)第34页状态转移矩阵性质由状态转移矩阵意义,有x(t2)=Φ(t2-t1)x(t1)=Φ(t2-t1)[Φ(t1-t0)x(t0)]=[Φ(t2-t1)Φ(t1-t0)]x(t0)而x(t2)=Φ(t2-t0)x(t0)所以，性质(7)表明,在系统状态转移过程中,既能够将系统一步状态转移分解成多步状态转移,也能够将系统多步状态转移等效为一步状态转移,如图3-2所表示。图3-2系统状态转移第35页状态转移矩阵性质例3-2

求以下系统状态转移矩阵逆矩阵。解：对于该系统,在例3-1已求得状态转移矩阵为因为Φ-1(-t)=Φ(t),所以求得状态转移矩阵逆矩阵为第36页非齐次状态方程解3.1.3非齐次状态方程解当线性定常连续系统含有输入作用时,其状态方程为以下非齐次状态方程:x’=Ax+Bu该状态方程在初始状态下解,也就是由初始状态x(t0)和输入作用u(t)所引发系统状态运动轨迹。第37页非齐次状态方程解下面用两种求解常微分方程方法直接求解法拉氏变换法讨论非齐次状态方程解,以及解表示式意义输出方程解第38页直接求解法1.直接求解法将状态方程x’=Ax+Bu移项,可得x’-Ax=Bu将上式两边左乘以e-At,则有e-At[x’-Ax]=e-AtBu即d(e-Atx)/dt=e-AtBu在区间[t0,t]内对上式积分,则有第39页直接求解法即上式便是非齐次状态方程解。当t0=0时,解x(t)又可记为所以第40页直接求解法若用状态转移矩阵来表示,上述非齐次状态方程解又可分别记为第41页拉氏变换法2.拉氏变换法将该非齐次状态方程两边取拉氏变换,可得sX(s)-x0=AX(s)+BU(s)即X(s)=(sI-A)-1[x0+BU(s)]其中X(s)和U(s)分别为x(t)和u(t)拉氏变换。对上式两边取拉氏反变换,并利用卷积分公式,则有上述求解关键为等式右边第二项。第42页拉氏变换法下面先回顾卷积积分拉氏变换法则。设W1(s)和W2(s)分别为原函数f1(t)和f2(t)拉氏变换,则f1(t)和f2(t)卷积拉氏变换为结果与直接求解法完全相同。对上述状态方程求解式利用卷积分公式,则有第43页状态方程解意义3.状态方程解意义由前面讨论非齐次状态方程解知,线性定常连续系统状态方程解由两个部分相加组成。第一个部分是由初始状态所引发自由运动,它是系统初始状态对系统状态转移影响,与初始时刻后输入无关,称为状态零输入响应。第二个部分是由输入所引发系统强迫运动,其值为输入函数与矩阵指数函数卷积。所以,它与输入相关,与系统初始状态无关,称为状态零状态响应。记得电路理论中暂态电路解否?第44页状态方程解意义状态方程解表明,系统在任意时刻状态取决于系统初始状态x(t0)和从初始时刻t0以来输入。假如人为地选择输入信号(施以控制),就能够使系统状态在状态空间中取得所期望状态轨线。第45页输出方程解4.输出方程解由非齐次状态方程解x(t),可得输出方程y=Cx+Du输出响应为或或第46页输出方程解线性定常连续系统输出解由3个部分相加组成。第一个部分是由初始状态所引发自由运动第二个部分是由输入所引发系统强迫运动。第三个部分是由直联项引发前馈响应。或第47页输出方程解--例3-3例3-3

已知线性定常系统为试求系统在单位阶跃输入作用下,状态方程解。解在例3-1中已求出状态转移矩阵

(t)为于是,系统状态方程在阶跃输入u(t)=1(t)下解为第48页输出方程解—例3-3第49页系统脉冲响应3.1.4系统脉冲响应当系统输入为单位脉冲函数时,系统在零初始状态时输出响应称为脉冲响应。单位脉冲函数

(t)可用下式来定义:下面讨论线性定常系统脉冲响应。第50页系统脉冲响应由线性定常连续系统输出y(t)表示式,可得系统脉冲响应H(t),即为由卷积分性质可得,上式积分结果为H(t)=CeAtB=L-1[C(sI-A)-1B]所以,脉冲响应也反应了系统输入与输出间动态传递关系。第51页状态转移矩阵计算3.2状态转移矩阵计算

在状态方程求解中,关键是状态转移矩阵

(t)计算。对于线性定常连续系统,该问题又归结为矩阵指数函数eAt计算。上一节已经介绍了基于拉氏反变换技术矩阵指数函数eAt计算方法,下面讲述计算矩阵指数函数下述其它3种惯用方法。级数求和法约旦规范形法

化eAt为A有限多项式矩阵函数法第52页级数求和法3.2.1级数求和法

由上一节对矩阵指数函数定义过程中可知:矩阵指数函数eAt计算可由上述定义式直接计算。因为上述定义式是一个无穷级数,故在用此方法计算eAt时必须考虑级数收敛性条件和计算收敛速度问题。类似于标量指数函数eat,对全部有限常数矩阵A和有限时间t来说,矩阵指数函数eAt这个无穷级数表示收敛。第53页级数求和法显然,用此方法计算eAt普通不能写成封闭、简练解析形式,只能得到数值计算近似计算结果。其计算精度取决于矩阵级数收敛性与计算时所取项数多少。假如级数收敛较慢,则需计算级数项数多,人工计算是非常麻烦,普通只适合用于计算机计算。所以,该方法缺点:计算量大精度低非解析方法,难以得到计算结果简练解析表示式。第54页级数求和法—例3-4例3-4

用直接计算法求下述矩阵矩阵指数函数:解按矩阵指数函数展开式计算以下:第55页约旦规范形法

3.2.2约旦规范形法

上节给出了对角线矩阵、块对角矩阵和约旦块三种特殊形式矩阵矩阵指数函数。因为任何矩阵都可经线性变换成为对角线矩阵或约旦矩阵,所以可经过线性变换将普通形式矩阵变换成对角线矩阵或约旦矩阵,再利用上述特殊形式矩阵矩阵指数函数来快速计算矩阵矩阵指数函数。下面讨论之。第56页约旦规范形法下面首先讨论矩阵指数函数一条性质:对矩阵A,经变换矩阵P作线性变换后,有则对应地有以下矩阵指数函数变换关系第57页约旦规范形法依据上述性质,对矩阵A,可经过线性变换方法得到对角线矩阵或约旦矩阵,然后利用该类特殊矩阵矩阵指数函数,由矩阵指数函数变换关系来求原矩阵A矩阵指数函数。该结论可简单证实以下:第58页约旦规范形法—例3-5例3-5

试求以下系统矩阵矩阵指数函数解1.先求A特征值。由特征方程可求得特征值为

1=-1

2=-2

3=-32.求特征值所对应特征向量。由前述方法可求得特征值

1,

2和

3所对应特征向量分别为p1=[101]

p2=[124]

p3=[169]

第59页约旦规范形法—例3-5故将A变换成对角线矩阵变换矩阵P及其逆阵P-1为3.由系统矩阵和矩阵指数函数变换关系,分别有第60页约旦规范形法—例3-6例3-6

试求以下系统矩阵矩阵指数函数第61页约旦规范形法—例3-6解1.先求A特征值。由特征方程可求得特征值为

1=2

2=

3=-12.

因为矩阵A为友矩阵,故将A变换成约旦矩阵变换矩阵P和其逆阵P-1分别为3.

由系统矩阵和矩阵指数函数变换关系,分别有第62页约旦规范形法--例3-6第63页塞尔维斯特内插法3.2.3塞尔维斯特内插法在讨论塞尔维斯特(Sylvester)内插法计算矩阵指数函数eAt时,需要用到关于矩阵特征多项式凯莱-哈密顿(Cayley-Hamilton)定理以及最小多项式概念。所以,首先给出凯莱-哈密顿定理及最小多项式概念,再讨论塞尔维斯特内插法。下面依次介绍:凯莱-哈密顿定理最小多项式塞尔维斯特内插法计算矩阵指数函数第64页凯莱-哈密顿定理1.

凯莱-哈密顿定理凯莱-哈密顿定理是矩阵方程分析和求解中非常主要定理,其表述和证实以下。定理3-1(凯莱-哈密顿定理)

设n

n矩阵A特征多项式为f(

)=|

I-A|=

n+a1

n-1+…+an-1+an则矩阵A必使由上述特征多项式决定矩阵多项式函数f(A)=An+a1An-1+…+an-1A+anI=0上述特征多项式亦称为矩阵A零化特征多项式。□第65页凯莱-哈密顿定理证实因为I=(

I-A)-1(

I-A)=[adj(

I-A)/|

I-A|](

I-A)故|

I-A|I=adj(

I-A)(

I-A)由伴随矩阵定义可知,伴随矩阵adj(

I-A)可表示为以下多项式矩阵函数:adj(

I-A)=

n-1I+n-2B2+…+Bn-1+Bn其中矩阵B2,B3,…,Bn为n

n维常数矩阵。第66页凯莱-哈密顿定理所以由前面两式,有(

n+a1

n-1+…+an-1+an)I=(

n-1I+n-2B2+…+Bn-1+Bn)(

I-A)整理得(

n+a1

n-1+…+an-1+an)I

=

nI+(B2-A)

n-1+…+(Bn-Bn-1A)-BnA第67页凯莱-哈密顿定理上式中,令等号两边

同幂次项系数相等,则有a1I-B2+A=0a2I-B3+AB2=0…an-1I-Bn+ABn-1=0anI+ABn=0所以,将上述各等式从上至下依次右乘以An-1,…,A,I,然后将各等式相加,即得An+a1An-1+…+an-1A+anI=0故矩阵A满足其本身零化特征多项式。第68页最小多项式

2.最小多项式

依据凯莱-哈密尔顿定理,任一n×n维矩阵A满足其本身特征方程,即特征多项式为A一个零化多项式。然而特征多项式不一定是A最小阶次零化多项式。将矩阵A满足最小阶次首一零化多项式称为最小多项式,也就是说,定义n×n维矩阵A最小多项式为满足

(A)=Am+

1Am-1+…+

m-1A+

mI=0,m

n阶次最低首一多项式

(

)=

m+

1

m-1+…+

m-1

+

m第69页最小多项式最小多项式在矩阵多项式分析与计算中起着主要作用。定理3-2给出了特征多项式与最小多项式关系。定理3-2

设首一多项式d(

)是

I-A伴随矩阵adj(

I-A)全部元素最高条约式,则最小多项式为第70页最小多项式证实

由假设知,矩阵adj(

I-A)最高条约式为d(

),故adj(

I-A)=d(

)B(

),式中,B(

)n2个元素(为

函数)最高条约式为1。因为:(

I-A)adj(

I-A)=|

I-A|I可得:d(

)(

I-A)B(

)=|

I-A|I由上式可知,特征多项式|

I-A|可被整除d(

)。所以设d(

)整除|

I-A|得到因式记为

(

),故有

|

I-A|=d(

)

(

),第71页最小多项式因为首一多项式d(

)最高阶次系数为1,所以

(

)最高阶次系数也应为1。所以,综合上两式,可得(

I-A)B(

)=

(

)I因而:

(A)=0即

(

)亦为A零化多项式。设

(

)为A最小多项式,所以零化多项式

(

)可写为

(

)=g(

)

(

)+e(

)其中g(

)和e(

)分别是多项式

(

)除以

(

)商和余项,且e(

)阶次低于

(

)。第72页最小多项式因为

(A)=0和

(A)=0,所以必定有e(A)=0。考虑到

(

)为矩阵A最小多项式,所以不存在比

(

)阶次还低A零化多项式,故e(

)必为零,即有

(

)=g(

)

(

)又因为

(A)=0,所以

(

)可写为

(

)I=(

I-A)H(

)式中,H(

)为

(

)一个因子矩阵,故

(

)I=g(

)

(

)I=g(

)(

I-A)H(

)将上式与(

I-A)B(

)=

(

)I比较,有B(

)=g(

)H(

)第73页最小多项式又因为B(

)n2个元素最高条约式为1,所以g(

)=1于是

(

)=

(

)所以,由前面证实|

I-A|=d(

)

(

)而证实了最小多项式

(

)为第74页最小多项式依据上述定理3-2,n×n维矩阵A最小多项式可按以下步骤求出。1)

依据伴随矩阵adj(

I-A),写出作为

因式分解多项式adj(

I-A)各元素;2)确定作为伴随矩阵adj(

I-A)各元素最高条约式d(

)。选取d(

)最高阶次系数为1。假如不存在条约式,则d(

)=1;3)

最小多项式

(

)可由|

I-A|除以d(

)得到。第75页塞尔维斯特内插法计算矩阵指数函数3.塞尔维斯特内插法计算矩阵指数函数基于最小多项式(或特征多项式),塞尔维斯特内插法能够非常简练、快速地计算出矩阵指数函数,其计算思想与过程可描述以下。若

(

)=

m+

1

m-1+…+

m-1

+

m为矩阵A最小多项式,则由

(A)=0有Am=-

1Am-1-…-

m-1A-

mI即Am可用有限项Am-1,…,A,I线性组合来表示。第76页塞尔维斯特内插法计算矩阵指数函数将上式两边乘以矩阵A,则有即Am+1可用有限项Am-1,…,A,I线性组合来表示。第77页塞尔维斯特内插法计算矩阵指数函数依次类推,则可知,Ai(i>m)可用有限项Am-1,…,A,I线性组合来表示。所以,我们有其中

i(t)(i=0,1,…,m-1)为待定关于时间t函数。即,矩阵指数函数eAt亦能够用有限项Am-1,…,A,I线性函数组合表示。第78页塞尔维斯特内插法计算矩阵指数函数利用上式去计算矩阵指数函数eAt关键是怎样计算待定函数

i(t)。下面分A特征值互异A有重特征值两种情况来讨论怎样计算

i(t)以及eAt。第79页A特征值互异(1)A特征值互异设矩阵An个互异特征值为

1,

2,…,

n,则矩阵A最小多项式

(

)等于特征多项式f(

)=|

I-A|=

n+a1

n-1+…+an-1

+an。因系统全部特征值

i使特征多项式f(

i)=0,故与前面证实过程类似,我们亦有其中待定函数

i(t)(i=0,1,…,n-1)与矩阵指数函数eAt表示式中

i(t)一致。第80页A特征值互异所以,可得以下待定函数

i(t)(i=0,1,…,n-1)线性方程组:求解上述方程得函数

i(t)后,由式(3-49)可计算得矩阵指数函数eAt。第81页A特征值互异--例3-7

例3-7

试求以下系统矩阵矩阵指数函数解

因为矩阵A3个特征值互异,并分别为-1,-2和-3,所以解方程组(3-52)可得第82页A特征值互异则系统状态转移矩阵为第83页A有重特征值(2)A有重特征值因为矩阵A与它约旦矩阵

含有相同最小多项式

(

),所以由前面推导过程可知,约旦矩阵

也满足设A与

特征值

i代数重数为mi,则由上式很轻易证实

i(t)满足求解上述方程,则可求得待定函数

i(t)。第84页A有重特征值为清楚说明问题,设A和

有以下6个特征值:

1,

1,

1,

2,

2,

3。则对应矩阵指数函数计算式(3-49)中待定函数

i(t)(i=0,

1,…,5)计算式为第85页A有重特征值—例3-8值得指出是,上述塞尔维斯特内插法不但对矩阵A最小多项式成立,而且对全部矩阵A零化多项式也成立。所以,在难以求解最小多项式时,上述方法中最小多项式可用矩阵A特征多项式代替,所得结果一致,仅计算量稍大。例3-8

试求以下系统矩阵矩阵指数函数第86页A有重特征值—例3-8解

解矩阵A特征方程,得特征值为1,1和2。因为特征值2为二重特征值,下面按基于最小多项式和特征多项式两种多项式用塞尔维斯特插值法计算矩阵指数函数。第87页A有重特征值—例3-8(1)

基于最小多项式计算。先计算伴随矩阵所以,伴随矩阵adj(

I-A)各元素最高条约式为(

-2),故最小多项式

(

)为第88页A有重特征值—例3-8因为最小多项式阶次为2,则依据塞尔维斯特插值法,矩阵指数函数能够表示为所以,待定函数

i(t)(i=0,

1)计算以下则系统矩阵指数函数为第89页A有重特征值—例3-8(2)

基于特征多项式计算。因为特征多项式阶次为3,则依据塞尔维斯特插值法,矩阵指数函数能够表示为所以,待定函数

i(t)(i=0,

1,2)计算以下第90页A有重特征值—例3-8则系统矩阵指数函数为第91页线性时变连续系统状态方程解3.3线性时变连续系统状态方程解严格说来,实际控制对象都是时变系统,其系统结构或参数随时间改变。如电机温升造成电阻以及系统数学模型改变;电子器件老化使其特征也发生改变;火箭燃料消耗造成其质量以及运动方程参数改变等。不过,因为时变系统数学模型较复杂,且不易于系统分析、优化和控制,所以只要实际工程允许,都可将慢时变系统在一定范围内近似地作为定常系统处理。但对控制目标要求较高高精度控制系统,需作为时变系统处理。第92页线性时变连续系统状态方程解下面将讨论线性时变连续系统状态方程求解问题,依次讨论:线性时变连续系统齐次状态方程解线性时变连续系统状态转移矩阵非齐次状态方程解第93页线性时变连续系统齐次状态方程解3.3.1线性时变连续系统齐次状态方程解当系统没有外部输入作用时,线性时变连续系统状态方程为齐次状态方程,可表示为x’(t)=A(t)x(t)这里讨论其满足初始状态解,也就是由初始时刻t0初始状态x(t0)所引发无输入强迫项(无外力)时自由运动。为确保该齐次状态方程解存在性和唯一性,在系统时间定义域[t0,tf]内,A(t)各元素为时间t分段连续函数。第94页线性时变连续系统齐次状态方程解下面证实时变系统齐次状态方程解为x(t)=(t,t0)x(t0)式中,

(t,t0)为时变系统状态转移矩阵,它定义为以下矩阵微分方程解。第95页线性时变连续系统齐次状态方程解证实对解表示式x(t)=(t,t0)x(t0)求导,则有且:x(t0)=(t0,t0)x(t0)=x(t0)说明式x(t)=(t,t0)x(t0)满足齐次状态方程及其初始条件。依据微分方程解唯一性，所以它是齐次状态方程解。时变系统齐次状态方程解表示了系统自由运动特征,也代表了初始状态x(t0)转移,其转移特征完全由状态转移矩阵Φ(t,t0)决定。第96页线性时变连续系统状态转移矩阵3.3.2线性时变连续系统状态转移矩阵下面深入讨论前面引入状态转移矩阵,主要内容为:状态转移矩阵求解状态转移矩阵性质第97页状态转移矩阵求解1.状态转移矩阵求解对于线性时变连续系统,状态转移矩阵Φ(t,t0)是以下矩阵微分方程和初始条件

’(t)=A(t)

(t),

(t)|t=0=I

解,它是一个n×n维关于时间变量t和t0矩阵函数。为了求得状态转移矩阵Φ(t,t0)表示式,可在时间域内对该矩阵微分方程积分,即有第98页状态转移矩阵求解假如将上式中积分号内Φ(

1,t0)再按上式展开,则有然后按此法继续迭代下去,第99页状态转移矩阵求解于是,可得一个由无穷项之和组成状态转移矩阵

(t,t0),即上式就是线性时变连续系统状态转移矩阵计算公式。在普通情况下,它不能写成封闭解析形式。在实际应用此公式时,可按一定精度要求,用数值积分计算方法去近似计算t1时刻Φ(t1,t0)值。第100页状态转移矩阵求解当初变系统矩阵A(t)满足以下条件时,时变系统状态转移矩阵解能够表示为指数形式。也就是说,只有A(t)与

A(

)d

满足矩阵乘法可交换条件时,上述指数表示形式解才成立。下面对这个条件给予证实。第101页状态转移矩阵求解将该指数表示形式右边展开成级数形式,有假如上式是系统状态转移矩阵,它必须满足状态转移矩阵定义式。于是,将上式两边对时间取导数,依据状态转移矩阵解表示式,状态转移矩阵Φ(t,t0)导数可表示为第102页状态转移矩阵求解比较上述两式可知,只有A(t)和

A(

)d

满足乘法可交换条件时,时变系统状态转移矩阵能够表示为指数形式。所以,线性时变连续系统齐次状态方程解也可表示为指数形式,即第103页状态转移矩阵求解上述A(t)和

A(

)d

可交换条件普通较难以检验是否成立。实际上,依据该可交换条件有上式对于任意时间变量t和t0都成立充分必要条件是:对于任意t1和t2,下式成立A(t1)A(t2)=A(t2)A(t1)所以,实际上较易于检验条件可取代A(t)和

A(

)d

可交换条件,成为时变系统状态转移矩阵解可表示为指数矩阵形式充分必要条件。第104页状态转移矩阵性质2.状态转移矩阵性质时变系统状态转移矩阵性质以下。1)

Φ(t,t)=I2)

传递性

(t2,t1)

(t1,t0)=

(t2,t0)第105页状态转移矩阵性质证实因为x(t2)=

(t2,t0)x(t0)且x(t2)=

(t2,t1)x(t1)=

(t2,t1)

(t1,t0)x(t0)故有

(t2,t0)x(t0)=

(t2,t1)

(t1,t0)x(t0)因为上式对任意初始状态x(t0)都成立,所以有

(t2,t0)=

(t2,t1)

(t1,t0)第106页状态转移矩阵性质3)

可逆性

-1(t,t0)=

(t0,t)证实由性质1)和2),有

(t,t0)

(t0,t)=

(t,t)=I故

-1(t,t0)=

(t0,t)成立。第107页状态转移矩阵性质4)对角线矩阵状态转移矩阵。假如时变系统矩阵A(t)以下表示对角线矩阵。A(t)=diag{a11(t)a22(t)…ann(t)}式中,aii(t)(i=1,2,…,n)为标量函数,则A(t)状态转移矩阵Φ(t,t0)为以下对角线矩阵。

(t,t0)=diag{

11(t,t0)

22(t,t0)…

nn(t,t0)}式中,

ii(t,t0)(i=1,2,…,n)为满足以下标量微分方程状态转移函数即第108页状态转移矩阵性质5)块对角矩阵状态转移矩阵。假如时变系统矩阵A(t)以下表示块对角矩阵。A(t)=block-diag{A1(t)A2(t)…Al(t)}式中,Ai(t)(i=1,2,…,l)为mi×mi维分块矩阵函数,则A(t)状态转移矩阵

(t,t0)为以下块对角矩阵。

(t,t0)=block-diag{

1(t,t0)

2(t,t0)…

l(t,t0)}式中,i(t,t0)(i=1,2,…,l)为满足以下矩阵微分方程状态转移矩阵第109页状态转移矩阵性质—例3-9例3-9

求以下时变系统状态转移矩阵Φ(t,t0)。解首先检验矩阵A(t)和

A(

)d

与是否可交换。为此计算第110页状态转移矩阵性质—例3-9所以A(t1)A(t2)=A(t2)A(t1)

t1，t2即矩阵A(t)和

A(

)d

与满足可交换条件,可由指数展开式方法计算状态转移矩阵,即第111页状态转移矩阵性质—例3-9因为而于是第112页非齐次状态方程解3.3.3非齐次状态方程解当含有外加输入作用时,其状态方程为以下非齐次状态方程:x’(t)=A(t)x(t)+B(t)u(t)该状态方程在初始状态下解,也就是由初始状态x(t0)和输入作用u(t)所引发系统状态运动轨迹。第113页非齐次状态方程解下面将证实当输入u(t)为分段连续时,该非齐次状态方程解为证实先设该非齐次状态方程解为显然,有式中,η(t)为待定函数。第114页非齐次状态方程解将所设解代入该状态方程左边,有将所设解代入该非齐次状态方程右边,有所以有即第115页非齐次状态方程解对上式两端积分,可得故该非齐次状态方程解为当系统状态空间模型中输出方程为y(t)=C(t)x(t)+D(t)u(t)时,系统输出为第116页非齐次状态方程解比较线性定常连续系统与线性时变连续系统状态方程解表示形式:定常系统时变系统初始状态影响初始时刻后输入影响,为脉冲响应函数与输入卷积第117页非齐次状态方程解与线性定常连续系统状态方程和输出方程解比较可知,线性时变连续系统与线性定常连续系统解结构和形式相同,都为状态零输入响应和零状态响应和。线性定常连续系统状态方程和输出方程解可视为线性时变连续系统对应解一个特殊形式。在A(t)为时不变时,时变系统状态转移矩阵Φ(t,t0)即为定常系统状态转移矩阵Φ(t-t0)。由此能够看出引入状态转移矩阵主要性。只有引入状态转移矩阵,才能使时变系统和定常系统求解公式建立统一形式。第118页非齐次状态方程解例3-10

求以下时变系统在阶跃输入时状态变量值。解

由例3-9有第119页非齐次状态方程解由时变系统状态方程解表示式,有第120页线性连续系统状态空间模型离散化3.4线性连续系统状态空间模型离散化离散系统工作状态能够分为以下两种情况。整个系统工作于单一离散状态。对于这种系统,其状态变量、输入变量和输出变量全部是离散量,如现在全数字化设备、计算机集成制造系统等。系统工作在连续和离散两种状态混合状态。对于这种系统,其状态变量、输入变量和输出变量现有连续时间型模拟量,又有离散时间型离散量，如连续被控对象采样控制系统就属于这种情况。第121页线性连续系统状态空间模型离散化对于第2种情况系统,其状态方程现有一阶微分方程组又有一阶差分方程组。为了能对这种系统利用离散系统分析方法和设计方法,要求整个系统统一用离散状态方程来描述。由此,提出了连续系统离散化问题。在计算机仿真、计算机辅助设计中利用数字计算机分析求解连续系统状态方程,或者进行计算机控制时,都会碰到离散化问题。第122页线性连续系统状态空间模型离散化图3-3所表示为连续系统化为离散系统系统框图。图3-3连续系统离散化实现第123页线性连续系统状态空间模型离散化线性连续系统时间离散化问题数学实质,就是在一定采样方式和保持方式下,由系统连续状态空间模型来导出等价离散状态空间模型,并建立起二者各系数矩阵之间关系式。为使连续系统离散化过程是一个等价变换过程,必须满足以下条件和假设。在离散化之后,系统在各采样时刻状态变量、输入变量和输出变量值保持不变。保持器为零阶,即加到系统输入端输入信号u(t)在采样周期内不变,且等于前一采样时刻瞬时值,故有:u(t)=u(kT)kT≤t<(k+1)T

第124页线性连续系统状态空间模型离散化采样周期T选择满足申农(Shannon)采样定理,即采样频率2

/T大于2倍连续信号x(k)上限频率。满足上述条件和假设,即可推导出连续系统离散化状态空间模型。下面分别针对线性定常连续系统和线性时变连续系统讨论离散化问题。第125页线性定常连续系统离散化3.4.1线性定常连续系统离散化本节主要研究线性定常连续系统状态空间模型离散化,即研究怎样基于采样将线性定常连续系统进行离散化,建立对应线性定常离散系统状态空间模型。主要讨论问题为两种离散化方法:准确法和近似法第126页线性定常连续系统离散化线性定常连续系统状态空间模型离散化,实际上是指在采样周期T下,将状态空间模型变换成离散系统以下状态空间模型:因为离散化主要是对描述系统动态特征状态方程而言,输出方程为静态代数方程,其离散化后应保持不变,即C(T)=C

D(T)=D离散化主要针对连续系统状态方程(A,B)怎样经过采样周期T,变换成离散系统状态方程(G,H)。第127页线性定常连续系统离散化在上述条件和假设下,即可推导出连续系统离散化状态空间模型。下面介绍两种离散化方法:准确法、近似法。第128页准确离散化方法1.准确离散化方法所谓线性定常连续系统状态方程准确离散化方法,就是利用状态方程求解公式以确保状态在采样时刻连续状态方程和离散化状态方程有相同解来进行离散化。连续系统状态方程求解公式以下:现在只考虑在采样时刻t=kT和t=(k+1)T时刻之间状态响应,即对于上式,取t0=kT,t=(k+1)T,于是第129页准确离散化方法考虑到u(t)在采样周期内保持不变假定,所以有将上式与线性定常离散系统状态方程x((k+1)T)=(I+AT)x(kT)+BTu(kT)比较,可知两式对任意x(kT)和u(kT)成立条件为G(T)=(T)=eAT对上式作变量代换,令t=(k+1)T-

,则上式可记为上两式即为准确离散化法计算式。第130页准确离散化方法—例3-11例3-11

试用准确离散化方法写出以下连续系统离散化系统状态方程:解首先求出连续系统状态转移矩阵:第131页准确离散化方法—例3-11依据准确法计算式有于是该连续系统离散化状态方程为第132页近似离散化方法2.近似离散化方法所谓线性定常连续系统状态方程近似离散化方法是指在采样周期较小,且对离散化精度要求不高情况下,用状态变量差商代替微商来求得近似差分方程。即,因为x’(kT)=LimT0[x((k+1)T)-x(kT)]/T故当采样周期较小时,有x’(kT)[x((k+1)T)-x(kT)]/T第133页近似离散化方法将上式代入连续系统状态方程,有[x((k+1)T)-x(kT)]/T=Ax(kT)+Bx(kT)即x((k+1)T)=(I+AT)x(kT)+BTu(kT)将上式与线性定常离散系统状态空间模型状态方程比较,则可得以下近似离散化计算公式:G(T)=I+ATH(T)=BT将上述近似离散法和准确离散法比较知,因为I+AT和BT分别是eAT和

eAtdtBTaylor展开式中一次近似,所以近似离散化方法其实是取准确离散化方法对应计算式一次Taylor近似展开式。第134页近似离散化方法—例3-12由上述推导过程可知,普通说来,采样周期T越小,则离散化精度越高。但考虑到实际计算时舍入误差等原因,采样周期T不宜太小。例3-12

试用近似离散化方法写出以下连续系统离散化系统状态方程:第135页近似离散化方法—例3-12解

由近似离散化法计算公式,对本例有于是该连续系统离散化状态方程为第136页近似离散化方法—例3-12对上述近似离散化法精度可检验以下:1.当T=1s时,准确法计算结果为近似法计算结果为2.当T=0.001s时,准确法计算结果为第137页近似离散化方法—例3-12从上述计算结果可知,近似离散法只适合用于较小采样周期。近似法计算结果为第138页线性时变连续系统离散化3.4.2线性时变连续系统离散化线性时变连续系统状态空间模型离散化,实际上是指在指定采样周期T下,将连续系统状态方程变换成线性时变离散系统以下状态方程:第139页线性时变连续系统离散化线性时变连续系统状态方程离散化,就是利用时变系统状态轨迹求解公式来进行离散化。由3.3节可知,连续系统状态方程解可表示为:现在只考虑在采样时刻t=kT和t=(k+1)T时刻之间状态响应,即对于上式,取t0=kT,t=(k+1)T,于是考虑到u(t)在采样周期内保持不变,所以有第140页线性时变连续系统离散化比较下述两式可得线性时变连续系统离散化模型各矩阵以下第141页线性时变连续系统离散化例3-13

试写出以下线性时变连续系统离散化系统状态方程。解由例3-9,该系统转移矩阵函数为第142页线性时变连续系统离散化所以,由上述离散化计算公式，可分别计算第143页线性时变连续系统离散化将上述计算所得G(k)和H(k)代入,则求得离散化状态方程以下第144页线性离散系统状态方程解3.5线性离散系统状态方程解本节研究线性定常离散系统方程解,需处理主要问题:状态转移矩阵状态转移矩阵性质状态方程求解状态方程解各部分意义输出方程解第145页线性离散系统状态方程解线性定常离散时间系统状态方程求解有递推法和Z变换法两种主要方法:Z变换法只能适合用于线性定常离散系统,递推法可推广到时变系统和非线性系统。下面将分别讨论线性定常离散系统线性时变离散系统状态空间模型求解。第146页线性定常离散系统状态方程解3.5.1线性定常离散系统状态方程解下面介绍线性定常离散系统状态方程求解递推法和Z变换法。最终讨论输出方程解第147页递推法1.递推法递推法亦称迭代法。用递推法求解线性定常离散时间系统状态方程x(k+1)=Gx(k)+Hu(k)时,只需在状态方程中依次令k=0,1,2,…,从而有x(1)=Gx(0)+Hu(0)x(2)=Gx(1)+Hu(1)=G2x(0)+GHu(0)+Hu(1)……第148页递推法若给出初始状态x(0),即可递推算出x(1),x(2),x(3),…重复以上步骤,能够得到以下线性离散系统状态方程递推求解公式:上述递推计算公式中第2项为离散卷积,所以有以下另一形式线性离散系统状态方程解表示式第149页递推法若初始时刻k0不为0,则上述状态方程解可表示为:或第150页递推法与连续系统状态方程求解类似,对线性离散系统状态方程求解,亦可引入状态转移矩阵。该状态转移矩阵是以下差分方程初始条件解:

(k+1)=G

(k)

(0)=I用递推法求解上述定义式,可得

(k)=Gk所以,可得线性定常离散系统状态方程另一个解表示形式:第151页递推法比较连续系统与离散系统状态方程解表示形式:连续系统离散系统初始状态影响初始时刻后输入影响,为脉冲响应函数与输入卷积第152页递推法对上述离散系统状态方程求解公式,有以下几点说明:1.与连续系统类似,离散系统状态响应也由两部分组成,一部分为由初始状态引发响应,与初始时刻后输入无关,称为系统状态零输入响应;另一部分是由初始时刻后输入所引发响应,与初始时刻状态值无关,称为系统状态零状态响应。2.引入状态转移矩阵概念和表示之后,线性连续系统和线性离散系统状态方程求解公式在形式上一致,都由零输入响应和零状态响应叠加组成,只是对应零状态响应在形式上略有不一样,一为求积分(卷积),一为求和(离散卷积),但本质是一致。第153页递推法3.在由输入所引发状态响应中,第k个时刻状态只取决于此采样时刻以前输入采样值,而与该时刻输入采样值u(k)无关。这即为计算机控制系统固有一步时滞。第154页递推法下面讨论几个特殊形式系统矩阵G状态转移矩阵(1)对角线矩阵。

当G为以下对角线矩阵:G=diag{

1

2…

n}则状态转移矩阵为(2)块对角矩阵。

当G为以下块对角矩阵:G=block-diag{G1

G2…Gl}其中Gi为mi

mi维分块矩阵,则状态转移矩阵为第155页递推法其中

kj=k!/[(k-j)!j!]为二项式系数。(3)约旦块矩阵。

当Gi为特征值为

imi

mi维约旦块,则分块矩阵矩阵指数函数为第156页递推法(4)

对系统矩阵G,当存在线性变换矩阵P,使得G~=P-1GP则有第157页Z变换法2.Z变换法已知线性定常离散系统状态方程为x(k+1)=Gx(k)+Hu(k)对上式两边求Z变换,可得zX(z)-zx(0)=GX(z)+HU(z)于是(zI-G)X(z)=zx(0)+HU(z)用(zI-G)-1左乘上式两边,有X(z)=(zI-G)-1zx(0)+(zI-G)-1HU(z)对上式进行Z反变换,有x(k)=Z-1[(zI-G)-1zx(0)]+Z-1[(zI-G)-1HU(z)]第158页Z变换法其中W1(z)和W2(z)分别为w1(k)和w2(k)Z变换。将上述公式推广到向量函数和矩阵函数,则可得在Z反变换中对标量函数存在下述公式和性质:第159页Z变换法—例3-14该表示式与前面递推法求解结果一致。例3-14

已知某系统状态方程和初始状态分别为试求系统状态在输入u(k)=1时响应。所以,离散系统状态方程解为:第160页Z变换法—例3-14类似地,可继续递推下去,直到求出所需要时刻解为止。2.用Z变换法求解。先计算(zI-G)-1解1.用递推法求解。分别令k=1,2,3,…,则由状态方程有第161页Z变换法—例3-14所以,有第162页Z变换法—例3-14由Z变换,有u(k)=1U(z)=z/(z-1)所以,有X(z)=(zI-G)-1[zx(0)+HU(z)]第163页Z变换法—例3-14令k=0,1,2,3代入上式,可得第164页输出方程解3.输出方程解将状态方程解代入以下线性定常离散系统输出方程:y(k)=Cx(k)+Du(k)中,可得输出y(k)解为第165页输出方程解或第166页线性时变离散系统状态方程解3.5.2线性时变离散系统状态方程解设线性时变离散系统状态空间模型为式中,初始时刻为k0;初始状态为x(k0)。假定系统状态方程解存在且惟一,则解为式中,

(k,k0)称为线性时变离散系统状态转移矩阵。第167页线性时变离散系统状态方程解线性时变离散系统状态转移矩阵

(k,k0)满足以下矩阵差分方程及初始条件：其解为第168页线性时变离散系统状态方程解与线性定常离散系统类似,线性时变离散系统状态求解公式可用迭代法证实。对线性时变离散系统状态方程，依次令k=k0,k0+1,k0+2,…,从而有第169页线性时变离散系统状态方程解所以有第170页线性时变离散系统状态方程解由上述状态方程解公式可知,线性时变离散系统状态方程解也包含两项。其中,第1项是由初始状态激励,为零输入响应,描述了输入向量为零时系统自由运动。第2项对应初始状态为零时,由输入向量激励响应,称为强迫运动或受控运动。线性时变离散系统运动状态取决于状态转移矩阵

(k,k0),而又是由

(k,k0)唯一决定。第171页线性时变离散系统状态方程解将状态响应代入输出方程,得到系统输出为,可见,系统输出响应也是由零输入响应、零状态响应和直接传输部分3项组成。第172页Matlab问题3.6Matlab问题本章中包括计算问题主要有矩阵指数函数计算、系统运动轨迹计算(即状态空间模型求解)以及连续系统离散化(采样)。基于Matlab基本函数和工具箱,能够进行上述系统运动分析计算和仿真。第173页Matlab问题为更加好地进行动态系统运动分析计算和仿真,编著者设计了一个Matlab符号化和图形化控制系统运动分析软件平台lti_analysis。这里将包括新Matlab程序设计方法,如使用符号计算工具箱进行矩阵指数函数和运动轨迹符号计算、使用图形用户界面(GUI)设计控制系统仿真与试验软件平台。第174页Matlab问题下面分别介绍基于Matlab上述问题程序编制和计算方法，主要有矩阵指数函数计算线性定常连续系统状态空间模型求解

连续系统离散化线性定常离散系统状态空间模型求解

线性定常系统运动分析符号计算和仿真平台

第175页矩阵指数函数计算3.6.1矩阵指数函数计算矩阵指数函数计算问题有两类,一类是数值计算,即给定矩阵A和详细时间t值,计算矩阵指数eAt值;另一类是符号计算,即在给定矩阵A下,计算矩阵指数函数eAt封闭(解析)矩阵函数表示式。数值计算问题可由基本Matlab函数完成,符号计算问题后一类则需要用到Matlab符号工具箱。下面就分别介绍eAt数值计算eAt符号计算第176页eAt数值计算1.eAt数值计算在Matlab中,给定矩阵A和时间t值,计算矩阵指数eAt值能够直接采取基本矩阵函数expm()。Matlabexpm()函数采取帕德(Pade)迫近法计算矩阵指数eAt,精度高,数值稳定性好。expm()函数主要调用格式为Y=expm(X)其中,X为输入需计算矩阵指数矩阵,Y为计算结果。第177页eAt数值计算Matlab问题3-1试在Matlab中计算例3-1中矩阵A在t=0.3时矩阵指数eAt值。Matlab程序m3-1以下。Matlab程序m3-1执行结果以下。A=[01;-2-3];t=0.3;eAt=expm(A*t)eAt=0.93280.1920-0.38400.3568第178页eAt数值计算在Matlab中还有3个计算矩阵指数eAt函数,分别是expmdemo1(),expmdemo2()和expmdemo3()。expmdemo1()就是expm(),采取帕德迫近法计算矩阵指数;而expmdemo2()采取3.2.1节中介绍利用泰勒级数展开法来计算,精度较低;expmdemo3()采取3.2.2节中介绍利用特征值和特征向量来计算对角线矩阵,进而经过对角线矩阵矩阵指数计算原矩阵矩阵指数。第179页eAt数值计算expmdemo3()计算精度取决于特征值、特征向量、指数函数exp()计算精度,因为这3种计算有良好计算方法,所以expmdemo3()计算精度最高。但expmdemo3()只能计算矩阵独立特征向量数等于矩阵维数,即矩阵能变换为对角线矩阵情况,所以,在不能判定矩阵是否能变换为对角线矩阵时,尽可能采取函数expm()。第180页eAt符号计算2.eAt符号计算在Matlab中,对给定矩阵A,可经过符号计算工具箱函数expm()计算变量t矩阵指数函数eAt表示式。在使用Matlab符号计算工具箱计算时,需要定义符号变量,输入符号表示式与符号矩阵。下面介绍使用符号计算工具箱需要基本操作。第181页eAt符号计算1)定义符号变量定义(指定)符号变量语句格式为symstsx...该语句将符号t,s,x,…定义为符号变量。在该语句后,就能够输入和计算符号表示式与符号矩阵。第182页eAt符号计算2)输入符号表示式符号表示式输入可采取赋值语句方式,如赋值语句f1=

sin(x)^2+cos(y)^3-3

为定义符号表示式变量f1为表示式。在Matlab中,符号表示式输入采取符号串形式,其表示式格式与Matlab数值计算格式基本一致。第183页eAt符号计算3)输入符号矩阵Matlab中符号矩阵输入采取函数sym()。sym()调用格式为S=sym(A)该函数功效为将符号串A