根据基尔霍夫定律[RLC]有:
图1-3
和
令
得
以上方程式在M-文件electsys.m中定义如下:
functionxdot=electsys(t,x);%状态导数
V=1;%阶跃输入
R=1.4;L=2;C=0.32;
xdot=[x
(2)/C;1/L*(V-x
(1)-R*x
(2))];
下面的M-文件使用函数ode23对系统进行仿真。
t0=0;tfinal=15;﹡时间间隔0~15秒
x0=[0.5,0];%初值条件x1=0.5,x2=0
tol=0.001;%精度
trace=0;%若非零值则打印出每一步的计算值
[t,x]=ode23('electsys',t0,tfinal,x0,tol,trace);
subplot(211),plot(t,x);
title('TimeresponseofaRLCseriescircuit')
xlabel('Time-sec')
text(8,1.05,'Capacitorvoltage'),text(8,0.05,'current')
vc=x(:
1),i=x(:
2);
subplot(212),plot(vc,i);
title('currentversuscapacitorvoltage')
xlabel('CapacitorVoltage')
ylabel('current')
subplot(111)
仿真结果见图1-4。
图1-4
1.1.3非线性系统
在变量的一定范围内大多数物理系统是线性的。
然而,随着变量范围的延伸,所有的系
统最终都是非线性的。
对于非线性系统,叠加原理不再成立。
Ode23和ode45可以求解非线性微分方程,见例1-3。
例1.3
图1-5的单摆,一长L米的无重量绳悬挂重为W=mg千克
的物体,阻尼系数为B千克/米/秒。
通常,该系统用线性微分方程近似描
述,实际上它是非线性的。
如果θ为绳的倾角,物体的速度将为Lθ,那么离心力为
根据牛顿定理有
联结上两式,得
令
得
上面方程式在M-文件pendulum.m中定义如下:
functionxdot=pendulum(t,x);%状态导数
W=2;L=.6;B=0.02;g=9.81;m=W/g;
xdot=[x
(2);-B/m*x
(2)-W/(m*L)*sin(x
(1))];
下面的M-文件对系统进行仿真。
t0=0;tfinal=5;%时间间隔0~5秒
x0=[1,0];%初始条件x1=1,x2=0
tol=0.0001;trace=0;%精度,不打印每一步的计算值
[t,x]=ode23('pendulum',t0,tfinal,x0,tol,trace);
subplot(211),plot(t,x)
title('Timeresponseofarigidpendulum'),xlabel('Time-sec')
text(3.2,3.5,'Velocity'),text(3.2,1.2,'Angle-Rad')
th=x(:
1);w=x(:
2);
subplot(212),plot(th,w)
title('Phaseplaneplotofpendulum')
xlabel('Position-Rad'),ylabel('Angularvelocity')
subplot(111)
仿真结果见图1-6。
图1-6
二、用MATLAB建立控制系统模型
MATLAB是基于数学模型分析控制系统性能的功能强大的软件工具。
由于微分方程与传递函数的模型参数完全一致,在MATLAB中只需要用传递函数描述,这里介绍用MATLAB及Simulink建立不同形式的传递函数模型,进行不同形式传递函数之间的转换。
1.有理函数模型
线性系统的传递函数模型可一般地表示为:
(1)
将系统的分子和分母多项式的系数按降幂的方式以向量的形式输入给两个变量
和
,就可以轻易地将传递函数模型输入到MATLAB环境中。
命令格式为:
;
(2)
; (3)
在MATLAB控制系统工具箱中,定义了tf()函数,它可由传递函数分子分母给出的变量构造出单个的传递函数对象。
从而使得系统模型的输入和处理更加方便。
该函数的调用格式为:
G=tf(num,den); (4)
例2.1 一个简单的传递函数模型:
可以由下面的命令输入到MATLAB工作空间中去。
>> num=[1,5];
den=[1,2,3,4,5];
G=tf(num,den)
运行结果:
Transferfunction:
s+5
-----------------------------
s^4+2s^3+3s^2+4s+5
这时对象G可以用来描述给定的传递函数模型,作为其它函数调用的变量。
例2.2 一个稍微复杂一些的传递函数模型:
该传递函数模型可以通过下面的语句输入到MATLAB工作空间。
>>num=6*[1,5];
den=conv(conv([1,3,1],[1,3,1]),[1,6]);
tf(num,den)
运行结果
Transferfunction:
6s+30
-----------------------------------------
s^5+12s^4+47s^3+72s^2+37s+6
其中conv()函数(标准的MATLAB函数)用来计算两个向量的卷积,多项式乘法当然也可以用这个函数来计算。
该函数允许任意地多层嵌套,从而表示复杂的计算。
2.零极点模型
线性系统的传递函数还可以写成极点的形式:
(5)将系统增益、零点和极点以向量的形式输入给三个变量
、Z和P,就可以将系统的零极点模型输入到MATLAB工作空间中,命令格式为:
(6)
(7)
(8)
在MATLAB控制工具箱中,定义了zpk()函数,由它可通过以上三个MATLAB变量构造出零极点对象,用于简单地表述零极点模型。
该函数的调用格式为:
G=zpk(Z,P,KGain) (9)
例2.3 某系统的零极点模型为:
该模型可以由下面的语句输入到MATLAB工作空间中。
>>KGain=6;
z=[-1.9294;-0.0353+0.9287j;-0.0353-0.9287j];
p=[-0.9567+1.2272j;-0.9567-1.2272j;0.0433+0.6412j;0.0433-0.6412j];
G=zpk(Z,P,KGain)
运行结果:
Zero/pole/gain:
6(s+1.929)(s^2+0.0706s+0.8637)
----------------------------------------------
(s^2-0.0866s+0.413)(s^2+1.913s+2.421)
注意:
对于单变量系统,其零极点均是用列向量来表示的,故Z、P向量中各项均用分号(;)隔开。
3.反馈系统结构图模型
设反馈系统结构图如图5所示。
图5 反馈系统结构图
控制系统工具箱中提供了feedback()函数,用来求取反馈连接下总的系统模型,该函数调用格式如下:
G=feedback(G1,G2,sign); (10)
其中变量sign用来表示正反馈或负反馈结构,若sign=-1表示负反馈系统的模型,若省略sign变量,则仍将表示负反馈结构。
G1和G2分别表示前向模型和反馈模型的LTI(线性时不变)对象。
例2.4 若反馈系统图5中的两个传递函数分别为:
,
则反馈系统的传递函数可由下列的MATLAB命令得出
>>G1=tf(1,[1,2,1]);
G2=tf(1,[1,1]);
G=feedback(G1,G2)
运行结果:
Transferfunction:
s+1
---------------------
s^3+3s^2+3s+2
若采用正反馈连接结构输入命令
>>G=feedback(G1,G2,1)
则得出如下结果:
Transferfunction:
s+1
-----------------
s^3+3s^2+3s
例2.5 若反馈系统为更复杂的结构如图6所示。
其中
,
,
则闭环系统的传递函数可以由下面的MATLAB命令得出:
>> G1=tf([1,7,24,24],[1,10,35,50,24]);
G2=tf([10,5],[1,0]);
H=tf([1],[0.01,1]);
G_a=feedback(G1*G2,H)
得到结果:
Transferfunction:
0.1s^5+10.75s^4+77.75s^3+278.6s^2+361.2s+120
--------------------------------------------------------------------
0.01s^6+1.1s^5+20.35s^4+110.5s^3+325.2s^2+384s+120
图6 复杂反馈系统
4.有理分式模型与零极点模型的转换
有了传递函数的有理分式模型之后,求取零极点模型就不是一件困难的事情了。
在控制系统工具箱中,可以由zpk()函数立即将给定的LTI对象G转换成等效的零极点对象G1。
该函数的调用格式为:
G1=zpk(G)(11)
例2.6给定系统传递函数为:
对应的零极点格式可由下面的命令得出
>> num=[6.8,61.2,95.2];
den=[1,7.5,22,19.5,0];
G=tf(num,den);
G1=zpk(G)
显示结果:
Zero/pole/gain:
6.8(s+7)(s+2)
-------------------------
s(s+1.5)(s^2+6s+13)
可见,在系统的零极点模型中若出现复数值,则在显示时将以二阶因子的形式表示相应的共轭复数对。
同样,对于给定的零极点模型,也可以直接由MATLAB语句立即得出等效传递函数模型。
调用格式为:
G1=tf(G)(12)
例2.7 给定零极点模型:
可以用下面的MATLAB命令立即得出其等效的传递函数模型。
输入程序的过程中要注意大小写。
>> Z=[-2,-7];
P=[0,-3-2j,-3+2j,-1.5];
K=6.8;
G=zpk(Z,P,K);
G1=tf(G)
结果显示:
Transferfunction:
6.8s^2+61.2s+95.2
-------------------------------
s^4+7.5s^3+22s^2+19.5s
5.Simulink建模方法
在一些实际应用中,如果系统的结构过于复杂,不适合用前面介绍的方法建模。
在这种情况下,功能完善的Simulink程序可以用来建立新的数学模型。
Simulink是由MathWorks软件公司1990年为MATLAB提供的新的控制系统模型图形输入仿真工具。
它具有两个显著的功能:
Simul(仿真)与Link(连接),亦即可以利用鼠标在模型窗口上“画”出所需的控制系统模型。
然后利用SIMULINK提供的功能来对系统进行仿真或线性化分析。
与MATLAB中逐行输入命令相比,这样输入更容易,分析更直观。
下面简单介绍SIMULINK建立系统模型的基本步骤:
(1)SIMULINK的启动:
在MATLAB命令窗口的工具栏中单击按钮
或者在命令提示符>>下键入simulink命令,回车后即可启动Simulink程序。
启动后软件自动打开Simullink模型库窗口,如图7所示。
这一模型库中含有许多子模型库,如Sources(输入源模块库)、Sinks(输出显示模块库)、Nonlinear(非线性环节)等。
若想建立一个控制系统结构框图,则应该选择File|New菜单中的Model选项,或选择工具栏上newModel
按钮,打开一个空白的模型编辑窗口如图8所示。
图7 simulink模型库
图8 模型编辑窗口
(2)画出系统的各个模块:
打开相应的子模块库,选择所需要的元素,用鼠标左键点中后拖到模型编辑窗口的合适位置。
(3)给出各个模块参数:
由于选中的各个模块只包含默认的模型参数,如默认的传递函数模型为1/(s+1)的简单格式,必须通过修改得到实际的模块参数。
要修改模块的参数,可以用鼠标双击该模块图标,则会出现一个相应对话框,提示用户修改模块参数。
(4)画出连接线:
当所有的模块都画出来之后,可以再画出模块间所需要的连线,构成完整的系统。
模块间连线的画法很简单,只需要用鼠标点按起始模块的输出端(三角符号),再拖动鼠标,到终止模块的输入端释放鼠标键,系统会自动地在两个模块间画出带箭头的连线。
若需要从连线中引出节点,可在鼠标点击起始节点时按住Ctrl键,再将鼠标拖动到目的模块。
(5)指定输入和输出端子:
在Simulink下允许有两类输入输出信号,第一类是仿真信号,可从source(输入源模块库)图标中取出相应的输入信号端子,从Sink(输出显示模块库)图标中取出相应输出端子即可。
第二类是要提取系统线性模型,则需打开Connection(连接模块库)图标,从中选取相应的输入输出端子。
例2.8 典型二阶系统的结构图如图9所示。
用SIMULINK对系统进行仿真分析。
图9 典型二阶系统结构图
按前面步骤,启动simulink并打开一个空白的模型编辑窗口。
(1)画出所需模块,并给出正确的参数:
●在sources子模块库中选中阶跃输入(step)图标,将其拖入编辑窗口,并用鼠标左键双击该图标,打开参数设定的对话框,将参数steptime(阶跃时刻)设为0。
●在Math(数学)子模块库中选中加法器(sum)图标,拖到编辑窗口中,并双击该图标将参数Listofsigns(符号列表)设为|+-(表示输入为正,反馈为负)。
●在continuous(连续)子模块库中、选积分器(Integrator)和传递函数(TransferFcn)图标拖到编辑窗口中,并将传递函数分子(Numerator)改为〔900〕,分母(Denominator)改为〔1,9〕。
●在sinks(输出)子模块库中选择scope(示波器)和Out1(输出端口模块)图标并将之拖到编辑窗口中。
(3)将画出的所有模块按图9用鼠标连接起来,构成一个原系统的框图描述如图10所示。
(4)选择仿真算法和仿真控制参数,启动仿真过程。
●在编辑窗口中点击Simulation|Simulationparameters菜单,会出现一个参数对话框,在solver模板中设置响应的仿真范围StartTime(开始时间)和StopTime(终止时间),仿真步长范围Maxinumstepsize(最大步长)和Mininumstepsize(最小步长)。
对于本例,StopTime可设置为2。
最后点击Simulation|Start菜单或点击相应的热键启动仿真。
双击示波器,在弹出的图形上会“实时地”显示出仿真结果。
输出结果如图11所示。
图10二阶系统的simulink实现
在命令窗口中键入whos命令,会发现工作空间中增加了两个变量――tout和yout,这是因为Simulink中的Out1模块自动将结果写到了MATLAB的工作空间中。
利用MATLAB命令plot(tout,yout),可将结果绘制出来,如图12所示。
比较11和12,可以发现这两种输出结果是完全一致的。
图11仿真结果示波器显示
图12MATLAB命令得出的系统响应曲线
三、MATLAB在控制系统模型转换中的应用
3.1传递函数向状态空间描述的转换
下面一组函数允许控制系统模型的不同表示形式之间可以互相转换。
(1)状态空间模型到传递函数模型的转换。
[num,den]=ss2tf(a,b,c,d,u)
(2)状态空间模型到零-极点增益模型的转换。
[z,p,k]=sstzp(a,b,c,d,u)
(3)传递函数模型到状态空间模型的转换。
[a,b,c,d]=tf2ss(num,den)
(4)传递函数模型到零-极点增益模型的转换。
[z,p,k]=tf2zp(num,den)
(5)零一极点增益模型到状态空间模型的转换。
[a,b,c,d]=zp2ss(z,p,k)
(6)零-极点增益模型到传递函数模型的转换。
[num,den]=zp2tf(z,p,k)(7)传递函数模型到部分分式模型的转换。
[z,p,k]=residue(num,den)(8)部分分式模型到传递函数模型的转换。
[num,den]=residue(z,p,k)(9)连续系统模型到离散系统模型的转换。
[ad,bd]=c2d(a,b,Ts)
(10)离散系统模型到连续系统模型的转换。
[a,b]=d2c(ab,bd,Ts)
控制系统工具箱包含一组模型转换的函数。
[A,B,C,D]=tf2ss(num,den)将传递函数转换为状态空间描述。
例3.1求下面传递函数的状态空间描述
>>num=[172];den=[192624];
>>[A,B,C,D]=tf2ss(num,den)
状态方程各矩阵如下:
D=0
3.2状态空间描述向传递函数的转换
已知状态方程和输出方程
y=Cx+Du
采用拉普拉斯变换有
Y(s)=C(sI-A)-1Bu(s)+Du(s)
则
函数ss2tf(A,B,C,D,i)是将状态空间描述转换为对第一个输入的传递函数。
[num,den]=ss2tf(A,B,C,D,i)是将状态空间描述化为分子、分母多项式形式的传递函数。
[z,p,k]=ss2zp(A,B,C,D,i)将状态空间描述化为零极点形式表示的传递函数。
例3.3一个系统的状态空间描述如下
y=[100]x
求传递函数G(s)=Y(s)/U(s)
>>A=[010;0001;-1-2-3];B=[10;0;0];
>>C=[100];D=[0];
>>[num,den]=ss2tf(A,B,C,D,1)
>>[z,p,k]=ss2zp(A,B,C,D,1)
其中,ss2tf(A,B,C,D,1)中“1”表示对第一个输入。
传递函数的分子、分母多项式系数如下:
num=
010.000030.000020.0000
den=
1.00003.00002.00001.0000
传递函数的零、极点如下:
z=
-1
-2
p=
-0.3376+0.5623i
-0.3376-0.5623i
-2.3247
k=10
因而传递函数为
3.3由方框图求状态空间描述和传递函数
控制系统工具箱中提供了函数[A,B,C,D]=connect(a,b,c,q,iu,iy)。
将方框图描述转换成状态空间描述和传递函数。
其中q