系统仿真实验报告5-6

控制系统仿真技术实验

实验报告

实验题目1.MATLAB 的符号计算与SIMULINK的使用

2.连续系统数字仿真及模型变换方法

              


一.MATLAB 符号计算

1)求矩阵对应的行列式和特征根

a=sym('[a11 a12;a21 a22]');

da=det(a)

da =

a11*a22-a12*a21

>> ea=eig(a)

ea =

 1/2*a11+1/2*a22+1/2*(a11^2-2*a11*a22+a22^2+4*a12*a21)^(1/2)

 1/2*a11+1/2*a22-1/2*(a11^2-2*a11*a22+a22^2+4*a12*a21)^(1/2)

>> a=sym('[1 2;2 1]');

a =

[ 1, 2]

[ 2, 1]

>> det(a)

ans =

-3

>> ea(a)

??? Error using ==> subsindex

Function 'subsindex' is not defined for values of class 'sym'.

>> ea=eig(a)

ea =

  3

 -1

(2)求方程的解(包括精确解和一定精度的解)

>> r1=solve('x^2-x-1')

rv=vpa(r1)

rv4=vpa(r1,4)

rv20=vpa(r1,20)

r1 =

 1/2*5^(1/2)+1/2

 1/2-1/2*5^(1/2)

rv =

  1.6180339887498948482045868343656

 -.61803398874989484820458683436560

 rv4 =

  1.618

 -.6180

rv20 =

  1.6180339887498948482

 -.61803398874989484820

(3)

>> a=sym('a');b=sym('b');c=sym('c');d=sym('d'); %定义4个符号变量

>> w=10;x=5;y=-8;z=11; %定义4个数值变量

>> A=[a,b;c,d] %建立符号矩阵A

>> B=[w,x;y,z] %建立数值矩阵B

>> det(A) %计算符号矩阵A的行列式

>> det(B) %计算数值矩阵B的行列式

A =

[ a, b]

[ c, d]

B =

    10     5

    -8    11

ans =

a*d-b*c

ans =

   150

(4)

>> syms x y;

>> s=(-7*x^2-8*y^2)*(-x^2+3*y^2);

>> expand(s) %对s展开

ans =

7*x^4-13*x^2*y^2-24*y^4

>> collect(s,x) %对s按变量x合并同类项(无同类项)

ans =

7*x^4-13*x^2*y^2-24*y^4

>>factor(ans) % 对ans分解因式

ans =

(8*y^2+7*x^2)*(x^2-3*y^2)

>> m=x^4+x^2+3*x^4;

>> collect(m,x)

ans =

4*x^4+x^2

>> factor(m)

ans =

x^2*(4*x^2+1)

>> 

(5)>> A=[34,8,4;3,34,3;3,6,8];

b=[4;6;2];

X=linsolve(A,b) %调用linsolve函数求解

X =

    0.0675

    0.1614

    0.1037

A\b %用另一种方法求解

ans =

    0.0675

    0.1614

    0.1037

>> 

(6)

>> syms a11 a12 a13 a21 a22 a23 a31 a32 a33 b1 b2 b3;

>> A=[a11,a12,a13;a21,a22,a23;a31,a32,a33];

>> b=[b1;b2;b3];

>> XX=A\b %用左除运算求解

XX =

 -(a12*a33*b2-a12*b3*a23+a13*a22*b3-a13*b2*a32+b1*a23*a32-b1*a22*a33)/(-a11*a23*a32+a23*a31*a12+a21*a13*a32-a12*a21*a33-a22*a31*a13+a11*a22*a33)

  (a11*a33*b2-a11*b3*a23-a31*a13*b2+b3*a21*a13-a33*a21*b1+a31*b1*a23)/(-a11*a23*a32+a23*a31*a12+a21*a13*a32-a12*a21*a33-a22*a31*a13+a11*a22*a33)

  (a21*b1*a32+a11*a22*b3-a11*b2*a32+b2*a31*a12-a12*a21*b3-a22*

>> A=[1,2,3,4;4,5,6,7]

A =

     1     2     3     4

     4     5     6     7

>> b=[1;2]

b =

     1

     2

>> x=linsolve(A,b)

x =

    0.1111

         0

         0

    0.2222

>> 

(7)

>> syms a b t x y z;

>> f=sqrt(1+exp(x));

>> diff(f) %未指定求导变量和阶数,按缺省规则处理

ans =

1/2/(1+exp(x))^(1/2)*exp(x)

>> f=x*cos(x);

>> diff(f,x,2) %求f对x的二阶导数

ans =

-2*sin(x)-x*cos(x)

>> diff(f,x,3) %求f对x的三阶导数

ans =

-3*cos(x)+x*sin(x)

>> f1=a*cos(t);f2=b*sin(t);

>> diff(f2)/diff(f1) %按参数方程求导公式求 y 对 x 的导数

ans =

-b*cos(t)/a/sin(t)

 

二. 控制系统时间响应分析

1.  step 2. impulse 3. initial 4. lsim 5. rlocfind

(1)

>> num=[4] ; den=[1 , 1 , 4] ;

>> step(num , den)

>> [y , x , t]=step(num , den) ;

>> tp=spline(y , t , max(y)) %计算峰值时间

>> max(y) %计算峰值

tp =

    1.6062

ans =

    1.4441

>> 

(2)>> a=[0,1;-6,-5];b=[0;1];c=[1,0];d=0;

     >> [y,x]=step(a,b,c,d);

     >> plot(y)

>>

(3)>> num=[4] ; den=[1 , 1 ,4] ;

     >> impulse(num,den)

>>

(4)

>> a=[0 , 1 ; -10 , -2] ; b=[0 ; 1] ;

>> c=[1 , 0] ; d=[0] ;

>> x0=[1 ,0];

>> subplot(1 , 2 , 1) ; initial(a , b , c ,d,x0)

>> subplot(1 , 2 , 2) ; impulse(a , b , c , d)

(5)

>> num=[1] ; den=[1 , 1] ;

>> t=0 : 0.01 : 10 ;

>> u=sin(2*t) ;

>> hold on

>> plot(t , u , 'r')

>> lsim(num , den , u , t)

 

第二篇:信号与系统仿真实验报告6

上机实验6:连续系统的复频域分析

验证试验

1

clear all;clc;

b=[1 0 -1];a=[1 2 3 2];

zs=roots(b);ps=roots(a);

plot(real(zs),imag(zs),'go',real(ps),imag(ps),'mx','markersize',12);

brighten(0.1);

grid;

legend('零点','极点');

brighten(0.1);

2

clear all;clc;

b=[1 0 -1];

a=[1 2 3 2];

zplane(b,a);

legend('零点','极点');

3

%稳态滤波法

clear all;clc;

w=8000;s=j*w;

num=[0 1e4 6e7];den=[1 875 88e6];

H=polyval(num,s)/polyval(den,s);

mag=abs(H);

phase=angle(H)/pi*180;

t=2:1e-6:2.002

vg=12.5*cos(w*t);

vo=12.5*mag*cos(w*t+phase*pi/180);

plot(t,vg,t,vo);grid;

text(0.25,0.85,'Output Voltage','sc');

text(0.07,0.35,'Input Voltage','sc');

title('稳态滤波输出'); ylabel('电压(v)');xlabel('时间(s)');

4

%拉氏变换法

syms s t;

Hs=sym('(10^4*(s+6000))/(s^2+875*s+88*10^6)');

Vs=laplace(12.5*cos(8000*t));Vos=Hs*Vs;

Vo=ilaplace(Vos);Vo=vpa(Vo,4);

ezplot(Vo,[1,1+5e-3]);hold on;

ezplot('12.5*cos(8000*t)',[1,1+5e-3]);axis([1 1+2e-3 -50 50]);

5

clear all;clc;

num=[1e11];den=[1 2.5e6 1e12 0];

[r,p,k]=residue(num,den);

%%%%%%%%%%%%以下是运行结果

% >> r

%

% r =

%

%     0.0333

%    -0.1333

%     0.1000

%

% >> p

%

% p =

%

%     -2000000

%      -500000

%            0

%

% >> k

%

% k =

%

%      []

%

% >>

设计实验

1_1

clear all;clc;

syms s t;

Hs=sym('(s+2)/(s^2+4*s+3)');

f1=laplace(diff(Heaviside(t)));f1s=Hs*f1;f1ss=ilaplace(f1s);f1ss=vpa(f1ss,4);

subplot(2,2,1);ezplot(f1ss,[1,1+5e-3]);hold on;

ezplot('(s+2)/(s^2+4*s+3)',[1,1+5e-3]);axis([1 1+2e-3 -0.25 0.25]);

f2=laplace(Heaviside(t));f2s=Hs*f2;f2ss=ilaplace(f2s);f2ss=vpa(f2ss,4);

subplot(2,2,2);ezplot(f2ss,[1,1+5e-3]);hold on;

ezplot('(s+2)/(s^2+4*s+3)',[1,1+5e-3]);axis([1 1+2e-3 -0.25 0.25]);

f3=laplace(cos(20*t)*Heaviside(t));f3s=Hs*f3;f3ss=ilaplace(f3s);f3ss=vpa(f3ss,4);

subplot(2,2,3);ezplot(f3ss,[1,1+5e-3]);hold on;

ezplot('(s+2)/(s^2+4*s+3)',[1,1+5e-3]);axis([1 1+2e-3 -0.25 0.25]);

f4=laplace(exp(-t)*Heaviside(t));f4s=Hs*f4;f4ss=ilaplace(f4s);f4ss=vpa(f4ss,4);

subplot(2,2,4);ezplot(f4ss,[1,1+5e-3]);hold on;

ezplot('(s+2)/(s^2+4*s+3)',[1,1+5e-3]);axis([1 1+2e-3 -0.25 0.25]);

1_3

clear all;clc;

b=[1.65 -0.331 -576 90.6 19080];a=[1 0.996 463 97.8 12131 8.11 0];

zs=roots(b);ps=roots(a);

subplot(2,1,1);

plot(real(zs),imag(zs),'go',real(ps),imag(ps),'mx','markersize',12);

brighten(0.1);

grid;

legend('零点','极点');

brighten(0.1);

syms s t;

Hs=sym('(1.65*s^4-0.331*s^3-576*s^2+90.6*s+19080)/(s^6+0.996*s^5+463*s^4+97.8*s^3+12131*s^2+8.11*s)');

f1=laplace(diff(Heaviside(t)));f1s=Hs*f1;f1ss=ilaplace(f1s);f1ss=vpa(f1ss,4);

subplot(2,1,2);ezplot(f1ss,[1,1+5e-3]);hold on;

ezplot('(1.65*s^4-0.331*s^3-576*s^2+90.6*s+19080)/(s^6+0.996*s^5+463*s^4+97.8*s^3+12131*s^2+8.11*s)',[1,1+5e-3]);

axis([1 1+2e-3 -1 1]);

   

相关推荐