第一篇:基于 Matlab 的離散控制系統(tǒng)仿真
2014 / 2015 學(xué)年第 1 學(xué)期
計算機(jī)控制技術(shù) 實
班 級 學(xué) 生 指 導(dǎo) 驗 報 告
學(xué) 號 1108030301 姓 名 蔡 夢 教 師 張 坤 鰲
實驗二 基于 Matlab 的離散控制系統(tǒng)仿真
一、實驗?zāi)康暮鸵螅?/p>
1、學(xué)習(xí)使用 Matlab 的命令對控制系統(tǒng)進(jìn)行仿真的方法
2、學(xué)習(xí)使用 Matlab 中的 Simulink 工具箱進(jìn)行系統(tǒng)仿真的方法
二、實驗環(huán)境
X86系列兼容型計算機(jī),Matlab軟件
三、實驗原理
1、控制系統(tǒng)命令行仿真
1)建立如圖所示一階系統(tǒng)控制模型并進(jìn)行系統(tǒng)仿真:
一階系統(tǒng)閉環(huán)傳遞函數(shù)為G(S)=
s1?333s=s?3,轉(zhuǎn)換為離散系統(tǒng)脈沖傳遞函數(shù)并仿真。
2)建立如圖所示二階系統(tǒng)控制模型并進(jìn)行系統(tǒng)仿真:
52s(s?2?0.4?5)25251?s(s?2?0.4?5)=s2?2?0.4?5s?52,二階系統(tǒng)閉環(huán)傳遞函數(shù)為G(S)=轉(zhuǎn)換為離散系統(tǒng)脈沖傳遞函數(shù)并仿真,改變參數(shù),觀察不同的系統(tǒng)的仿真結(jié)果。
2、控制系統(tǒng)的 Simulink 仿真
按圖建立系統(tǒng)的 Simulink 模型,對不同的輸入信號進(jìn)行仿真,改變參數(shù),觀察不同的仿真結(jié)果。
將上述系統(tǒng)離散化并仿真,觀察仿真結(jié)果
四、實驗步驟
1、根據(jù)實驗原理對控制系統(tǒng)進(jìn)行軟件仿真
2、觀察記錄輸出的結(jié)果,與理論計算值相比較
3、自行選擇參數(shù),練習(xí)仿真方法,觀察不同的仿真結(jié)果
5252s(s?2?0.4?5)s(s?2?0.4?5)?52521?1?s(s?2?0.4?5)s(s?2?0.4?5)進(jìn)行軟二階系統(tǒng)閉環(huán)傳遞函數(shù)為G(S)=件仿真如下圖:
分別進(jìn)行離散仿真:
五、實驗心得
針對這次實驗設(shè)計,我通過各種渠道,上課認(rèn)真學(xué)習(xí),請教老師、上網(wǎng)搜索,圖書館查閱,詢問同學(xué)等學(xué)習(xí)到了很多知識,一步步了解最少拍控制系統(tǒng)設(shè)計,鍛煉了自我學(xué)習(xí)能力。
盡管學(xué)習(xí)上遇到了很多困難,結(jié)果也差強(qiáng)人意。但我們在不斷處理困難的過程中磨練了處理事物的能力和耐心,也讓同學(xué)間學(xué)會了互相學(xué)習(xí),共享資源
第二篇:MATLAB與控制系統(tǒng)仿真實驗報告
《MATLAB與控制系統(tǒng)仿真》
實驗報告
2013-2014學(xué)年 第 1 學(xué)期
專業(yè): 班級: 學(xué)號: 姓名:
實驗三 MATLAB圖形系統(tǒng)一、實驗?zāi)康模?/p>
1.掌握繪制二維圖形的常用函數(shù)。2.掌握繪制三維圖形的常用函數(shù)。3.熟悉利用圖形對象進(jìn)行繪圖操作的方法。4.掌握繪制圖形的輔助操作。
二、實驗原理:
1,二維數(shù)據(jù)曲線圖
(1)繪制單根二維曲線 plot(x,y);(2)繪制多根二維曲線 plot(x,y)當(dāng)x是向量,y是有一維與x同維的矩陣時,則繪制多根不同顏色的曲線。當(dāng)x,y是同維矩陣時,則以x,y對應(yīng)列元素為橫、縱坐標(biāo)分別繪制曲線,曲線條數(shù)等于矩陣的列數(shù)。(3)含有多個輸入?yún)?shù)的plot函數(shù) plot(x1,y1,x2,y2,…,xn,yn)(4)具有兩個縱坐標(biāo)標(biāo)度的圖形 plotyy(x1,y1,x2,y2)2,圖形標(biāo)注與坐標(biāo)控制 1)title(圖形名稱); 2)xlabel(x軸說明)3)ylabel(y軸說明)4)text(x,y圖形說明)5)legend(圖例1,圖例2,…)
6)axis([xmin xmax ymin ymax zmin zmax])3, 圖形窗口的分割 subplot(m,n,p)4,三維曲線
plot3(x1,y1,z1,選項1,x2,y2,選項2,…,xn,yn,zn,選項n)5,三維曲面
mesh(x,y,z,c)與surf(x,y,z,c)。一般情況下,x,y,z是維數(shù)相同的矩陣。X,y是網(wǎng)格坐標(biāo)矩陣,z是網(wǎng)格點上的高度矩陣,c用于指定在不同高度下的顏色范圍。6,圖像處理
1)imread和imwrite函數(shù) 這兩個函數(shù)分別用于將圖象文件讀入matlab工作空間,以及將圖象數(shù)據(jù)和色圖數(shù)據(jù)一起寫入一定格式的圖象文件。
2)image和imagesc函數(shù) 這兩個函數(shù)用于圖象顯示。為了保證圖象的顯示效果,一般還應(yīng)使用colormap函數(shù)設(shè)置圖象色圖。
三、實驗儀器和設(shè)備:
計算機(jī)一臺(帶有MATLAB6.5以上的軟件環(huán)境)。
四、預(yù)習(xí)要求:
1.復(fù)習(xí)二維與三維圖形的繪圖函數(shù)。2.復(fù)習(xí)圖形輔助操作。
五、實驗內(nèi)容及步驟:
1,設(shè)y?[0.5?3sinx]cosx,在x=0~2π區(qū)間取101點,繪制函數(shù)曲線。21?x
2,已知y1=x2,y2=cos(2x),y3=y1*y2,完成下列操作:
(1)在同一坐標(biāo)系下用不同的顏色和線型繪制三條曲線;
(2)分別用條形圖、階梯圖、桿圖和填充圖繪制三條曲線。
3,已知
?x??,?x?02?e y???1In(x?1?x2),x?0??2在-5<=x<=5區(qū)間繪制函數(shù)曲線。
4,繪制函數(shù)的曲面圖和等高線
z?cosxcosye?x2?y24
其中x的21個值均勻分布在[-5,5]范圍,y的31個值均勻分布在[0,10],要求使用subplot(2,1,1)和subplot(2,1,2)將產(chǎn)生的曲面圖和登高圖畫在同一個窗口上。
5.畫出函數(shù)
z?x2?y2?sin(xy)的曲面及等高線圖。
x2y2?1繪制平面曲線,并分析參數(shù)a對其形狀的影響。6.根據(jù)2?a25?a2
四、心得體會:
通過這次實驗我能熟練掌握二維和三維圖以及其他特殊圖形的制作,弄清楚了基本的圖形操作規(guī)則,大大加深了我對matlab的興趣。
實驗二 MATLAB程序設(shè)計
一、實驗?zāi)康?/p>
1.掌握利用if語句實現(xiàn)選擇結(jié)構(gòu)的方法。
2.掌握利用switch語句實現(xiàn)多分支選擇結(jié)構(gòu)的方法。3.掌握利用for語句實現(xiàn)循環(huán)結(jié)構(gòu)的方法。4.掌握利用while語句實現(xiàn)循環(huán)結(jié)構(gòu)的方法。
二、實驗設(shè)備及條件
計算機(jī)一臺(帶有MATLAB6.5以上的軟件環(huán)境)。
三、實驗內(nèi)容
1.編寫求解方程ax2?bx?c?0的根的函數(shù)(這個方程不一定為一元二次方程,因a、b、c的不同取值而定),這里應(yīng)根據(jù)a、b、c的不同取值分別處理,有輸入?yún)?shù)提示,當(dāng)a?0,b?0,c~?0時應(yīng)提示“為恒不等式!”。并輸入幾組典型值加以檢驗。
clear,clc a=input('請輸入一個數(shù)a=');b=input('請輸入一個數(shù)b=');c=input('請輸入一個數(shù)c=');m=b^2-4*a*c;if a==0
if b==0
'為恒不等式'
end end
m=b^2-4*a*c;if m>0
x1=(-b+sqrt(m))/(2*a)
x2=(-b-sqrt(m))/(2*a)elseif m==0
x=(-b)/(2*a)else
'不存在正實根' end
2.輸入一個百分制成績,要求輸出成績等級A+、A、B、C、D、E。其中100分為A+,90分~99分為A,80分~89分為B,70分~79分為C,60分~69分為D,60分以下為E。
要求:(1)用switch語句實現(xiàn)。
(2)輸入百分制成績后要判斷該成績的合理性,對不合理的成績應(yīng)輸出出錯信息。
clear,clc for k=1:10
a(k)={89+k};b(k)={79+k};
c(k)={69+k};d(k)={59+k};end A=cell(3,6);A(1,:)={'a','b','c','d','e','f'};A(2,:)={85,76,95,100,40,65};for k=1:6
switch A{2,k}
case 100
r='A+';
case a
r='A';
case b
r='B';
case c
r='C';
case d
r='D';
otherwise
r='E';
end
A(3,k)={r};end A A =
'a'
'b'
'c'
[85]
[76]
[95]
'B'
'C'
'A'
'd'
'e'
[100]
[40]
'A+'
'E'
'f' [65] 'D'
3.利用for循環(huán)語句編寫計算n!的函數(shù)程序,取n分別為-89、0、3、5、10驗證其正確性(輸入n為負(fù)數(shù)時輸出出錯信息)。
clear,clc n=input('請輸入一個正數(shù)n=');if n<0
'輸入錯誤' elseif n==0
'n!=0' elseif n==1
'n!=1' else
y=1;
for i=1:1:n
y=y*i;
i=i+1;
end
y end 請輸入一個正數(shù)n=-89
ans =輸入錯誤 請輸入一個正數(shù)n=0
ans =n!=0 請輸入一個正數(shù)n=1
ans =n!=1 請輸入一個正數(shù)n=3
y =6 請輸入一個正數(shù)n=10
y =3628800
四、實驗心得體會:
通過本次實驗課,我能熟練運(yùn)用for循環(huán)語句,switch條件語句以及if條件語句的新用法,和在C中的區(qū)別。盡管如此,但是在實驗中依然容易把for循環(huán)跟C語言中的for語句弄混,最后經(jīng)過不懈努力下,終于弄明白了兩者之間的差別,使我能更好的運(yùn)用這些指令語句。
第三篇:控制系統(tǒng)的Matlab仿真與設(shè)計課后答案
MATLAB課后習(xí)題答案 2.1 x=[15 22 33 94 85 77 60] x(6)x([1 3 5])x(4:end)x(find(x>70))2.3 A=zeros(2,5);
A(:)=-4:5
L=abs(A)>3 islogical(L)
X=A(L)2.4 A=[4,15,-45,10,6;56,0,17,-45,0] find(A>=10&A<=20)2.5 p1=conv([1,0,2],conv([1,4],[1,1]));p2=[1 0 1 1];[q,r]=deconv(p1,p2);cq='商多項式為
';cr='余多項式為
';disp([cq,poly2str(q,'s')]),disp([cr,poly2str(r,'s')])2.6 A=[11 12 13;14 15 16;17 18 19];PA=poly(A)
PPA=poly2str(PA,'s')3.1 n=(-10:10)';y=abs(n);plot(n,y,'r.','MarkerSize',20)axis equal grid on xlabel('n')3.2 x=0:pi/100:2*pi;y=2*exp(-0.5*x).*sin(2*pi*x);plot(x,y),grid on;3.3 t=0:pi/50:2*pi;x=8*cos(t);y=4*sqrt(2)*sin(t);z=-4*sqrt(2)*sin(t);plot3(x,y,z,'p');
title('Line in 3-D Space');text(0,0,0,'origin');
xlabel('X'),ylable('Y'),zlable('Z');grid;3.4
theta=0:0.01:2*pi;
rho=sin(2*theta).*cos(2*theta);polar(theta,rho,'k');3.5
[x,y,z]=sphere(20);z1=z;
z1(:,1:4)=NaN;c1=ones(size(z1));surf(3*x,3*y,3*z1,c1);hold on z2=z;
c2=2*ones(size(z2));
c2(:,1:4)=3*ones(size(c2(:,1:4)));surf(1.5*x,1.5*y,1.5*z2,c2);colormap([0,1,0;0.5,0,0;1,0,0]);grid on hold off 第四章
function f=factor(n)if n<=1 f=1;else
f=factor(n-1)*n;end
function[s,p]=fcircle(r)s=pi*r*r;p=2*pi*r;
function k=jcsum1(n)k=0;i=0;while i<=n k=k+2^i;i=i+1;end
function k=jcsum(n)k=0;
for i=0:n k=k+2^i;end
4.1for m=100:999
m1=fix(m/100);m2=rem(fix(m/10),10);m3=rem(m,10);
if
m==m1*m1*m1+m2*m2*m2+m3*m3*m3 disp(m)
end end
4.2[s,p]=fcircle(10)4.3y=0;n=100;for i=1:n
y=y+1/i/i;end y
4.4s=0;for i=1:5
s=s+factor(i);end s
4.5sum=0;i=1;while sum<2000 sum=sum+i;i=i+1;end;n=i-2
4.6jcsum(63)jcsum1(63)4.1 for m=100:999
m1=fix(m/100);
m2=rem(fix(m/10),10);
m3=rem(m,10);
if m==m1*m1*m1+m2*m2*m2+m3*m3*m3
disp(m)
end end 4.3 y=0;n=100;
for i=1:n
y=y+1/i/i;end y 4.4 s=0;for i=1:5
s=s+factor(i);end s 4.5
sum=0;i=1;
while sum<2000
sum=sum+i;
i=i+1;end;n=i-2 4.6
i=0;k=0;while i<=63
k=k+2^i;
i=i+1;end k i
i=0;k=0;for i=0:63
k=k+2^i;end i k
第五章
function f=fxyz(u)x=u(1);y=u(2);z=u(3);f=x+y.^2./x/4+z.^2./y+2./z;
5.1A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];
b=[13,-9,6,0]';x=Ab
5.2[U,fmin]=fminsearch('fxyz',[0.5,0.5,0.5])
5.3X=linspace(0,2*pi,50);Y=sin(X);
P=polyfit(X,Y,3)AX=linspace(0,2*pi,50);Y=sin(X);Y1=polyval(P,X)
plot(X,Y,':O',X,Y1,'-*')
5.4x=0:2.5:10;h=[0:30:60]';
T=[95,14,0,0,0;88,48,32,12,6;67,64,54,48,41];xi=[0:0.5:10];hi=[0:10:60]';
temps=interp2(x,h,T,xi,hi,'cubic');
mesh(xi,hi,temps);
5.1 A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];b=[13,-9,6,0]';x=Ab 5.3 X=linspace(0,2*pi,50);Y=sin(X);P=polyfit(X,Y,3)AX=linspace(0,2*pi,50);Y=sin(X);Y1=polyval(P,X)plot(X,Y,':O',X,Y1,'-*')6.1syms x
y=finverse(1/tan(x))6.2syms x y
f=1/(1+x^2);g=sin(y);fg=compose(f,g)6.3syms x
g=(exp(x)+x*sin(x))^(1/2);dg=diff(g)
6.4F=int(int('x*exp(-x*y)','x'),'y')
6.5syms x
F=ztrans(x*exp(-x*10))6.6a=[0 1;-2-3];syms s
inv(s*eye(2)-a);
6.7 f=solve('a*x^2+b*x+c')6.8 f=solve('x+y+z=1','x-y+z=2','2*x-y-z=1')
6.9y=dsolve('D2y+2*Dy+2*y=0','y(0)=1','Dy(0)=0')ezplot(y),grid on
6.10a=maple('simplify(sin(x)^2+cos(x)^2);')
6.11f=maple('laplace(exp(-3*t)*sin(t),t,s);')6.12 syms t x
F=sin(x*t+2*t);L=laplace(F)第七章
function
[sys,x0,str,ts]=ww(t,x,u,flag)%?¨ò?á?D??μí3μ?Soˉêy A=[0,1;-0.4,-0.2];B=[0;0.2];C=[1,0];D=0;
switch flag, case 0,[sys,x0,str,ts]=mdlInitializeSizes(A,B,C,D);case 1,sys=mdlDerivatives(t,x,u,A,B,C,D);case 2,sys=mdlUpdate(t,x,u);case 3,sys=mdlOutputs(t,x,u,A,B,C,D);case 4,sys=mdlGetTimeOfNextVarHit(t,x,u);case 9,sys=mdlTerminate(t,x,u);otherwise
error(['Unhandled flag = ',num2str(flag)]);end
%=============================== function
[sys,x0,str,ts]=mdlInitializeSizes(A,B,C,D)sizes = simsizes;sizes.NumContStates = 2;sizes.NumDiscStates = 0;sizes.NumOutputs = 1;sizes.NumInputs = 1;sizes.DirFeedthrough = 1;sizes.NumSampleTimes = 1;sys = simsizes(sizes);x0 = [0;0];str = [];ts = [0 0];
%=============================== function
sys=mdlDerivatives(t,x,u,A,B,C,D)
sys = A*x+B*u;
%===============================
function sys=mdlUpdate(t,x,u)sys = [];
%=============================== function
sys=mdlOutputs(t,x,u,A,B,C,D)sys = C*x+D*u;
%=============================== function
sys=mdlGetTimeOfNextVarHit(t,x,u)
sampleTime = 1;sys = t + sampleTime;
%===============================
function sys=mdlTerminate(t,x,u)sys = [];
7.1
7.2
7.3
7.4
7.5
7.6
7.7
第八章
8.1num=[5];den=[1,2,2];sys=tf(num,den)8.1.2s = tf('s');
H = [5/(s^2+2*s+2)];H.inputdelay =2
8.1.3h=tf([0.5,0],[1,-0.5,0.5],0.1)
8.2num=2*[1,0.5];den=[1,0.2,1.01];
sys=tf(num,den)
[z,p,k]=tf2zp(num,den);zpk(z,p,k)
[A,B,C,D]=tf2ss(num,den);ss(A,B,C,D)
8.3 num=[1,5];den=[1,6,5,1];ts=0.1;
sysc=tf(num,den);sysd=c2d(sysc,ts,'tustin')8.4.0
8.4.1 %???μí3·??òí??????·?ú??DD±ào?,óD8??í¨μà,áDD′????í¨μà′?μYoˉêy r1=1;r2=2;c1=3;c2=4;G1=r1;G2=tf(1,[c1,0]);
G3=1;%ê?·?à?μ?oí??o?μ?μ?á???,2??üo?2¢,′?oˉ?a1 G4=-1;G5=1/r2;
G6=tf(1,[c2,0]);G7=-1;G8=-1;%?¨á¢?Tá??óμ?×′ì??????£Dí G=append(G1,G2,G3,G4,G5,G6,G7,G8)
%D′3??μí3μ?á??ó???ó
Q=[1 4 0 %í¨μà1μ?ê?è?ê?í¨μà4 2 1 7 %í¨μà2μ?ê?è?ê?í¨μà1,7 3 2 0
2 0 5 3 8 6 5 0 7 5 0
6 0];%?oo??ú′?oˉ?Dì??? %áD3??μí3μ?×üμ?ê?è?oíê?3???μ?±ào?
inputs=1;outputs=6;
%éú3é×éo?oó?μí3μ?×′ì??????£Dí sys=connect(G,Q,inputs,outputs)
8.4.2r1=1;r2=2;c1=3;c2=4;[A,B,C,D]=linmod('x84');[num,den]=ss2tf(A,B,C,D);sys=tf(num,den)
8.5A=[1,1,0;0,1,0;0,0,2];B=[0,0;1,0;0,-2];n=size(A)
Tc=ctrb(A,B);if n==rank(Tc)
disp('?μí3íêè??ü??');else
disp('?μí32?íêè??ü??');end 第九章
function [rtab,info]=routh(den)info=[];
vec1=den(1:2:length(den));
nrT=length(vec1);
vec2=den(2:2:length(den)-1);rtab=[vec1;vec2,zeros(1,nrT-length(vec2))];for k=1:length(den)-2, alpha(k)=vec1(1)/vec2(1);
for i=1:length(vec2),a3(i)=rtab(k,i+1)-alpha(k)*rtab(k+1,i+1);
end
if sum(abs(a3))==0 a3=polyder(vec2);
info=[info,'All elements in row ',...int2str(k+2)' are zeros;'];
elseif abs(a3(1)) info=[info,'Replaced first element;']; end rtab=[rtab;a3, zeros(1,nrT-length(a3))];vec1=vec2;vec2=a3;end 9.1num=[2,5,1];den=[1,2,3];bode(num,den);grid on;figure; nyquist(num,den); 9.2num=5*[1,5,6];den=[1,6,10,8];step(num,den);grid on;figure; impulse(num,den);grid on;9.3kosi=0.7;wn=6; num=wn^2;den=[1,2*kosi*wn,wn^2];step(num,den);grid on;figure; impulse(num,den);grid on;9.4den=[1,2,8,12,20,16,16];[rtab,info]=routh(den)a=rtab(:,1) if all(a>0) disp('?μí3ê??è?¨μ?'); else disp('?μí3ê?2??è?¨μ?'); end 9.5num=7*[1,5];den=conv([1,0,0],conv([1,10],[1,1])); [gm,pm,wg,wc]=margin(num,den) 9.1 >> sys=tf([2,5,1],[1,2,3])Transfer function: 2 s^2 + 5 s + 1---------------s^2 + 2 s + 3 >> rlocus(sys)>> nyquist(sys)>> bode(sys)9.2 >> G=tf(conv([5],[1,5,6]),[1,6,10,8]);>> step(G)>> impulse(G) sys=tf([5,25,30],[1,6,10,8]);>> step(sys)>> impulse(sys)9.4>> GH=tf(conv([7],[1,5]),conv([1,0,0],conv([1,10],[1,1]))); >> [Gm,Pm,Wcg,Wcp]=margin(GH)Gm = 0 Pm = -47.2870 Wcg = 0 Wcp = 1.4354 >> GH=tf(conv([7],[1,5]),conv([1,0,0],conv([1,10],[1,1]))); >> [Gm,Pm,Wcg,Wcp]=margin(GH)GH_close=feedback(GH,1)step(GH_close),grid on Gm = 0 Pm = -47.2870 Wcg = 0 Wcp = 1.4354 第十章 function s=bpts2s(bp,ts,delta)kosi=sqrt(1-1./(1+((1./pi).*log(1./bp)).^2)); wn=-log(delta.*sqrt(1-kosi.^2))/(kosi.*ts); s=-kosi.*wn+j.*wn.*sqrt(1-kosi.^2); function [ngc,dgc]=fa_lead(ng0,dg0,Pm,wc,w) ngv=polyval(ng0,j*wc);dgv=polyval(dg0,j*wc);g=ngv/dgv; thetag=angle(g);mg=abs(g);thetar=Pm*pi/180; tz=(1+mg*cos(thetar-thetag))/(-wc*mg*sin(thetar-thetag));tp=(cos(thetar-thetag)+mg)/(wc*sin(thetar-thetag));ngc=[tz,1];dgc=[tp,1]; function [ngc,dgc]=fg_lag_pm(ng0,dg0,w,Pm) [mu,pu]=bode(ng0,dg0,w);wgc=spline(pu,w,Pm+5-180);%2??μ?óè??ú×??à???£?èμ????μ?ê×÷?a?úí?μ????D?μ?ê ngv=polyval(ng0,j*wgc);dgv=polyval(dg0,j*wgc);g=ngv/dgv; alph=abs(1/g);T=10/alph*wgc, ngc=[alph*T,1];dgc=[T,1]; function [ngc,dgc]=fg_lag_wc(ng0,dg0,w,wc) ngv=polyval(ng0,j*wc);dgv=polyval(dg0,j*wc);g=ngv/dgv; alph=abs(1/g);T=10/(alph*wc);ngc=[alph*T,1];dgc=[T,1]; function [ngc,dgc]=fg_lead_pd(ng0,dg0,wc)ngv=polyval(ng0,j*wc);dgv=polyval(dg0,j*wc);g=ngv/dgv;mg0=abs(g); t=sqrt(((1/mg0)^2-1)/(wc^2));%·ù?μ?à?ó?aá? ngc=[t,1];dgc=[1];%Gc(s)=1+Ts function [ngc,dgc]=fg_lead_pm(ng0,dg0,Pm,w) [mu,pu]=bode(ng0,dg0,w);%?????-?μí3μ???êy?μ?ê?ìó|êy?Y [gm,pm,wcg,wcp]=margin(mu,pu,w);%?óè??-?μí3μ??à???£?èoí???D?μ?ê alf=ceil(Pm-pm+5);%?????????÷ìá1?μ?×?′ó3??°???è£? phi=(alf)*pi/180;%??×?′ó3??°??×a???a???èμ¥?? a=(1+sin(phi))/(1-sin(phi));%????a?μ dbmu=20*log10(mu);%?μí3μ???êy·ù?μ mm=-10*log10(a);%wm′|μ??????÷??êy·ù?μ wgc=spline(dbmu,w,mm);%2??μ?óè?wm£?è??awm£?wc T=1/(wgc*sqrt(a));%????T ngc=[a*T,1];dgc=[T,1]; function [ngc,dgc]=fg_lead_pm_wc(ng0,dg0,Pm,wc,w) [mu,pu]=bode(ng0,dg0,w);ngv=polyval(ng0,j*wc);dgv=polyval(dg0,j*wc); g=ngv/dgv;%?ó?-?μí3?ú?úí?μ????D?μ?ê′|μ??μ?ê?ìó|êy?YG0(jwc) theta=180*angle(g)/pi;%?-?μí3?ú?úí?μ????D?μ?ê′|μ??à???£?è£?μ¥???a?è alf=ceil(Pm-(theta+180)+5);% ×?′ó3??°?? phi=(alf)*pi/180; a=(1+sin(phi))/(1-sin(phi));dbmu=20*log10(mu);mm=-10*log10(a);wgc=spline(dbmu,w,mm);T=1/(wgc*sqrt(a)); KK=128;s1=-2+i*2*sqrt(3);a=2 ng0=[10];dg0=conv([1,0],conv([1,2],[1,8]));g0=tf(ng0,dg0); [ngc,dgc,k]=rg_lag(ng0,dg0,KK,s1,a); gc=tf(ngc,dgc)function s=kw2s(kosi,wn)s=-kosi.*wn+j*wn.*sqrt(1-kosi.^2); 10.1ng0=[1];dg0=10000*[1 0-1.1772]; g0=tf(ng0,dg0);%?ú×??a?·??ò?μ??aD£?y?μí3μ?′?μYoˉêy s=kw2s(0.7,0.5)%?úí?μ?±??·?÷μ???μ? ngc=rg_lead(ng0,dg0,s);gc=tf(ngc,1)g0c=tf(g0*gc);rlocus(g0,g0c); b1=feedback(g0,1);%?′D£?y?μí3μ?±??·′?μYoˉêy b2=feedback(g0c,1);%D£?yoó?μí3μ?±??·′?μYoˉêy figure,step(b1,'r--',b2,'b');grid on %????D£?y?°oó?μí3μ?μ¥???×?? KK=20;s1=-2+i*sqrt(6);a=1 ng0=[10];dg0=conv([1,0],[1,4]);g0=tf(ng0,dg0); [ngc,dgc,k]=rg_lag(ng0,dg0,KK,s1,a);gc=tf(ngc,dgc)g0c=tf(KK*g0*gc); b1=feedback(k*g0,1);b2=feedback(g0c,1);step(b1,'r--',b2,'b');grid on g0c=tf(KK*g0*gc);rlocus(g0,g0c); b1=feedback(k*g0,1); b2=feedback(g0c,1);figure,step(b1,'r--',b2,'b');grid on ng0=[1];dg0=conv([1,0,0],[1,5]);g0=tf(ng0,dg0);w=logspace(-3,3);KK=1;Pm=50; [ngc,dgc]=lead4(ng0,dg0,KK,Pm,w); gc=tf(ngc,dgc);g0c=tf(KK*g0*gc);bode(KK*g0,w);hold on,bode(g0c,w);grid on,hold off [gm,pm,wcg,wcp]=margin(g0c)Kg=20*log10(gm)g1=feedback(g0c,1);bode(g1),grid on, [mag,phase,w]=bode(g1);a=find(mag<=0.707*mag(1));wb=w(a(1))max(mag) b=find(mag==max(mag))wr=w(b) KK=40;Pm=50;ng0= KK *[1]; dg0=conv([1,0],conv([1,1],[1,4])); g0=tf(ng0,dg0); w=logspace(-2,4); [ngc,dgc]=fg_lead_pm(ng0,dg0,Pm, w)gc=tf(ngc,dgc),g0c=tf(g0*gc); b1=feedback(g0,1);b2=feedback(g0c,1); step(b1,'r--', b2,'b');grid on figure, bode(g0,'r--',g0c,'b',w), grid on,[gm,pm,wcg,wcp]=margin(g0c), Km=20*log10(gm) KK=200;bp=0.3;ts=0.7;delta=0.05; ng0=[1];dg0=conv([1,0],conv([0.1,1],conv([0.02 1],conv([0.01,1],[0.005 1]))));g0=tf(ng0,dg0); w=logspace(-4,3);t=[0:0.1:3];[mag,phase]=bode(KK*g0,w);[gm0,pm0,wg0,wc0]=margin(mag,phase,w),gm0=20*log10(gm0)%gm0 =-15.6769 %2?¢è·?¨?úí?μ??a?·′?μYoˉêy mr=0.6+2.5*bp; wc=ceil((2+1.5*(mr-1)+2.5*(mr-1)^2)*pi/ts), h=(mr+1)/(mr-1)w1=2*wc/(h+1), w2=h*w1 w1=wc/10;w2=25;ng1=[1/w1,1];dg1=conv([1/w2,1],conv([1,0],[1,0])); g1=tf(ng1,dg1); g=polyval(ng1,j*wc)/polyval(dg1,j*wc);K=abs(1/g);%???D?μ?ê′|·ù?μ?a1£??óK?μ g1=tf(K*g1) %3?¢è·?¨·′à??·?ú′?μYoˉêy h=tf(dg1,ng1);Kh=1/K;h=tf(Kh*h)%?úí??μ?êì?D?μ?μ1ì?D? %4?¢?é??D??ü??±ê g2=feedback(KK*g0,h);%D£?yoó£??μí3μ??a?·′?μYoˉêy b1=feedback(KK*g0,1);b2=feedback(g2,1); bode(KK*g0,'r--',g2,'b',h,'g',w);grid on figure,step(b1, 'r--',b2, 'b',t);grid on,[pos,tr,ts,tp]=stepchar(b2,delta) function [ngc,dgc]=lag2(ng0,dg0,w,KK,Pm)[mu,pu]=bode(KK*ng0,dg0,w);wgc=spline(pu,w,Pm+5-180), ngv=polyval(KK*ng0,j*wgc);dgv=polyval(dg0,j*wgc);g=ngv/dgv; alph=abs(1/g), T=10/alph*wgc, ngc=[alph*T,1];dgc=[T,1]; function [ngc,dgc]=lead4(ng0,dg0,KK,Pm,w)[mu,pu]=bode(KK*ng0,dg0,w);[gm,pm,wcg,wcp]=margin(mu,pu,w);alf=ceil(Pm-pm+5);phi=(alf)*pi/180; a=(1+sin(phi))/(1-sin(phi)), dbmu=20*log10(mu);mm=-10*log10(a); wgc=spline(dbmu,w,mm), T=1/(wgc*sqrt(a)),ngc=[a*T,1];dgc=[T,1]; function [ngc,dgc]=ra_lead(ng0,dg0,s1)ngv=polyval(ng0,s1);dgv=polyval(dg0,s1);g=ngv/dgv;thetag=angle(g);mg=abs(g);thetas=angle(s1);ms=abs(s1); tz=(sin(thetas)-mg*sin(thetag-thetas))/(mg*ms*sin(thetag));tp=-(mg*sin(thetas)+sin(thetag+thetas))/(ms*sin(thetag));ngc=[tz,1];dgc=[tp,1]; function [ngc,dgc,k]=rg_lag(ng0,dg0,KK,s1,a) ngv=polyval(ng0,s1);dgv=polyval(dg0,s1);g=dgv/ngv; k=abs(g);%?úí??÷μ???μ?′|μ??ù1ì?£??ò? beta=k/KK; [kosi1,wn1]=s2kw(s1); zc=-wn1*sin(a*pi/180)/sin(pi-atan(sqrt(1-kosi1^2)/kosi1)-(a*pi/180));%à?ó??y?ò?¨àí pc=beta*zc; ngc=beta*[1,-zc];dgc=[1,-pc]; function varargout=rg_lead(ng0,dg0,s1)if nargout==1 ngv=polyval(ng0,s1);dgv=polyval(dg0,s1);g=ngv/dgv; thetal=pi-angle(g); zc=real(s1)-imag(s1)/tan(thetal); t=-1/zc; varargout{1}=[t,1];elseif nargout==2 ngv=polyval(ng0,s1);dgv=polyval(dg0,s1); g=ngv/dgv;theta=angle(g);phi=angle(s1); if theta>0 phi_c=pi-theta; end if theta<0;phi_c=-theta end theta_z=(phi+phi_c)/2;theta_p=(phi-phi_c)/2; z_c=real(s1)-imag(s1)/tan(theta_z); p_c=real(s1)-imag(s1)/tan(theta_p); nk=[1-z_c];varargout{2}=[1-p_c];kc=abs(p_c/z_c); if theta<0 kc=-kc end varargout{1}=kc*nk;else error('ê?3?±?á?êy??2??yè·£?');end function [bp,ts]=s2bpts(s,delta)[kosi,wn]=s2kw(s); bp=exp(-kosi.*pi./sqrt(1-kosi.^2)); ts=-1./(kosi.*wn)*log(delta.*sqrt(1-kosi.^2)); function [kosi,wn]=s2kw(s)kosi=1./sqrt(1+(imag(s)/real(s)).^2); wn=-real(s)./kosi; %è?1?wn?a?o?μ£??òwnè??y£?2¢?òkosiè?·′ iwn=(wn<0);wn(iwn)=-wn(iwn);kosi(iwn)=-kosi(iwn); function [pos,tr,ts,tp]=stepchar(g0,delta) [y,t]=step(g0);[mp,ind]=max(y);dimt=length(t);yss=y(dimt); pos=100*(mp-yss)/yss;tp=t(ind);for i=1:dimt if y(i)>=1 tr=t(i); break; end end; for i=1:length(y) if y(i)<=(1-delta)*yss|y(i)>=(1+delta)*yss ts=t(i); end end 第十一章 11.1a=[0 1 0;0 0 1;-1-5-6];b=[0 0 1]'; p=[-2+4j;-2-4j;-10];K=acker(a,b,p)eig(a-b*K) 11.2a=[0 1 0;0 0 1;-6-11-6];b=[1,0,0]'; p=[-2+2*sqrt(3)*j;-2-2*sqrt(3)*j;-10]; K=acker(a,b,p)eig(a-b*K) 11.6A=[-1 0 0;0-2-3;0 0-3];B=[1 0;2 3;-3-3];C=[1 0 0;1 1 1 ]; [G,K,L]=decoupling(A,B,C) 11.8A=[0 20.6;1 0];b=[0 1]';c=[0 1];d=0; G=ss(A,b,c,d); Q=diag([1,0,0,0,0]);R=1; p=[-1.8+2.4j;-1.8-2.4j];[k,P]=lqr(A,b,Q,R);l=(acker(A',c',p))' Gc=-reg(G,k,l);zpk(Gc), eig(Gc.a), t=0:0.05:2; G_1=feedback(G*Gc,1);a1=eig(G_1.a), y_1=step(G_1,t); 第十二章 function [t,xx]=diffstate(G,H,x0,u0,N,T)xk=x0;u=u0;t=0 for k=1:N xk=G*xk+H*u;x(:,k)=xk;t=[t,k*T];end;xx=[x0,x]; 12.1 function sys=M601(t,x)u=1; sys=[x(2);x(3);-800*x(1)-80*x(2)-24*x(3)+u]; function [t,y]=ode4(A,B,C,D,x0,h,r,v,t0,tf) Ab=A-B*v*C;B=B;C=C;x=x0';y=0;t=t0; N=round((tf-t0)/h);for i=1:N k1=Ab*x+B*r; k2=Ab*(x+h*k1/2)+B*r;k3=Ab*(x+h*k2/2)+B*r;k4=Ab*(x+h*k3)+B*r;x=x+h*(k1+2*k2+2*k3+k4)/6;y=[y,C*x];t=[t,t(i)+h];end 12.1 tspan=[0,10];x0=[0,0,0]'; [t,y]=ode45('M601',tspan,x0);y1=800*y(:,1);plot(t,y1); 12.2 num=10;den=conv([1,0],conv([1,2],[1,3])); [A,B,C,D]=tf2ss(num,den);x0=[0,0,0];v=1;t0=0;tf=10;h=0.01;r=1; [t,y]=ode4(A,B,C,D,x0,h,r,v,t0,tf); plot(t,y),grid 12.3 12.4 g=[-2.8-1.4 0 0;1.4 0 0 0;-1.8-0.3-1.4-0.6;0 0 0.6 0];h=[1 0 1 0]';c=[0 0 0 1];d=0; x0=[0 0 0 0]';u=1;N=30;T=0.1; [t,xx]=diffstate(g,h,x0,u,N,T);plot(t,xx);y=c*xx;figure stairs(t,y)grid on 12.6 第十四章 14.1 clear all;load optcar.mat; t=signals(1,:);p=signals(2,:);v=signals(3,:);a=signals(4,:);theta=signals(5,:); subplot(4,1,1);plot(t,p);grid on;ylabel('????(m)');subplot(4,1,2);plot(t,v);grid on;ylabel('?ù?è(m/s)');subplot(4,1,3);plot(t,a);grid on;ylabel('?ó?ù?è(m/s2)');subplot(4,1,4);plot(t,theta);grid on;ylabel('???è(?è)'); 14.1 clear all load car.mat %??μ?è?μ?car.mat?Dμ?·???êμ?éêy?Y?á3? t=signals(1,:);x=signals(2,:);theta=signals(3,:);x1=signals(4,:);theta1=signals(5,:); plot(t,x,t,x1);ylabel('D?3μ????(m)'),grid on;% ????????á|×÷ó????ü???£Díoí??è·?£Díxμ?μ¥???×???ìó|?ú?? figure % ????????á|×÷ó????ü???£Díoí??è·?£ Díthetaμ?μ¥???×???ìó|?ú?? plot(t,theta,t,theta1);ylabel('°ú???μ(rad)'),grid on; 第二章 1>>x=[15 22 33 94 85 77 60] >>x(6)>>x([1 3 5])>>x(4:end)>>x(find(x>70))2>>T=[1-2 3-4 2-3];>>n=length(T);>>TT=T';>>for k=n-1:-1:0 >>B(:,n-k)=TT.^k;>>end >>B >>test=vander(T)3>>A=zeros(2,5);>>A(:)=-4:5 >>L=abs(A)>3 >>islogical(L)>>X=A(L)4>>A=[4,15,-45,10,6;56,0,17,-45,0] >>find(A>=10&A<=20)5>>p1=conv([1,0,2],conv([1,4],[1,1]));>>p2=[1 0 1 1];>>[q,r]=deconv(p1,p2);>>cq='商多項式為 ';cr='余多項式為 ';>>disp([cq,poly2str(q,'s')]),disp([cr,poly2str(r,'s')])6>>A=[11 12 13;14 15 16;17 18 19];>>PA=poly(A)>>PPA=poly2str(PA,'s')第三章 1>>n=(-10:10)';>>y=abs(n);>>plot(n,y,'r.','MarkerSize',20)>>axis equal >>grid on >>xlabel('n')2>>x=0:pi/100:2*pi;>>y=2*exp(-0.5*x).*sin(2*pi*x);>>plot(x,y),grid on;3>>t=0:pi/50:2*pi;>>x=8*cos(t);>>y=4*sqrt(2)*sin(t); >>z=-4*sqrt(2)*sin(t);>>plot3(x,y,z,'p');>>title('Line in 3-D Space');>>text(0,0,0,'origin');>>xlabel('X'),ylable('Y'),zlable('Z');grid;4>>theta=0:0.01:2*pi;>>rho=sin(2*theta).*cos(2*theta);>>polar(theta,rho,'k');5>>[x,y,z]=sphere(20);>>z1=z;>>z1(:,1:4)=NaN;>>c1=ones(size(z1));>>surf(3*x,3*y,3*z1,c1);>>hold on >>z2=z;>>c2=2*ones(size(z2));>>c2(:,1:4)=3*ones(size(c2(:,1:4)));>>surf(1.5*x,1.5*y,1.5*z2,c2);>>colormap([0,1,0;0.5,0,0;1,0,0]);>>grid on >>hold off 第四章 1>>for m=100:999 m1=fix(m/100);m2=rem(fix(m/10),10);m3=rem(m,10);if m==m1*m1*m1+m2*m2*m2+m3*m3*m3 disp(m)end end M文件: function[s,p]=fcircle(r)s=pi*r*r;p=2*pi*r;主程序: [s,p]=fcircle(10)3>>y=0;n=100;for i=1:n y=y+1/i/i;end >>y M文件: function f=factor(n)if n<=1 f=1;else f=factor(n-1)*n;end 主程序: >>s=0;for i=1:5 s=s+factor(i);end >>s 5>>sum=0;i=1;while sum<2000 sum=sum+i;i=i+1;end;>>n=i-2 6 for循環(huán)M文件: function k=jcsum(n)k=0;for i=0:n k=k+2^i;end 主程序: >>jcsum(63) While循環(huán)M文件: function k=jcsum1(n)k=0;i=0;while i<=n k=k+2^i;i=i+1;end 主程序: >>jcsum1(63)第五章 1>>A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];>>b=[13,-9,6,0]';>>x=Ab M文件: function f=fxyz(u)x=u(1);y=u(2);z=u(3);f=x+y.^2./x/4+z.^2./y+2./z;主程序: [U,fmin]=fminsearch('fxyz',[0.5,0.5,0.5])3>>X=linspace(0,2*pi,50);>>Y=sin(X);>>P=polyfit(X,Y,3)>>AX=linspace(0,2*pi,50);>>Y=sin(X);>>Y1=polyval(P,X)>>plot(X,Y,':O',X,Y1,'-*')4>>x=0:2.5:10;>>h=[0:30:60]';>>T=[95,14,0,0,0;88,48,32,12,6;67,64,54,48,41];>>xi=[0:0.5:10];>>hi=[0:10:60]';>>temps=interp2(x,h,T,xi,hi,'cubic');>>mesh(xi,hi,temps);第六章 1>>syms x >>y=finverse(1/tan(x))2>>syms x y >>f=1/(1+x^2);g=sin(y);>>fg=compose(f,g)3>>syms x >>g=(exp(x)+x*sin(x))^(1/2);>>dg=diff(g)4>>F=int(int('x*exp(-x*y)','x'),'y')5>>syms x >>F=ztrans(x*exp(-x*10))6>>a=[0 1;-2-3];>>syms s >>inv(s*eye(2)-a);7>>f=solve('a*x^2+b*x+c')8>>f=solve('x+y+z=1','x-y+z=2','2*x-y-z=1')9>>y=dsolve('D2y+2*Dy+2*y=0','y(0)=1','Dy(0)=0')>>ezplot(y),grid on 10>>a=maple('simplify(sin(x)^2+cos(x)^2);')11>>f=maple('laplace(exp(-3*t)*sin(t),t,s);') 12>>syms t x >>F=sin(x*t+2*t);>>L=laplace(F)第七章 第八章 1-1>>h=tf([5,0],[1,2,2])1-2>>s = tf('s');>>H = [5/(s^2+2*s+2)];>>H.inputdelay =2 1-3>>h=tf([0.5,0],[1,-0.5,0.5],0.1)2>>num=2*[1,0.5];den=[1,0.2,1.01];>>sys=tf(num,den)>>[z,p,k]=tf2zp(num,den);>>zpk(z,p,k)>>[A,B,C,D]=tf2ss(num,den);>>ss(A,B,C,D)3 >>num=[1,5];den=[1,6,5,1];ts=0.1;>>sysc=tf(num,den);>>sysd=c2d(sysc,ts,'tustin') >>r1=1;r2=2;c1=3;c2=4;>>[A,B,C,D]=linmod('x84');>>[num,den]=ss2tf(A,B,C,D);>>sys=tf(num,den)5>>A=[1,1,0;0,1,0;0,0,2];B=[0,0;1,0;0,-2];>>n=size(A)>>Tc=ctrb(A,B);if n==rank(Tc)disp('系統(tǒng)完全能控');else disp('系統(tǒng)不完全能控');end 第九章 1>>num=[2,5,1];den=[1,2,3];>>bode(num,den);grid on;>>figure;>>nyquist(num,den);2>>num=5*[1,5,6];den=[1,6,10,8];>>step(num,den);grid on;>>figure;>>impulse(num,den);grid on;3>>kosi=0.7;wn=6;>>num=wn^2;den=[1,2*kosi*wn,wn^2];>>step(num,den);grid on;>>figure;>>impulse(num,den);grid on;4 M文件: function [rtab,info]=routh(den)info=[];vec1=den(1:2:length(den));nrT=length(vec1);vec2=den(2:2:length(den)-1);rtab=[vec1;vec2, zeros(1,nrT-length(vec2))];for k=1:length(den)-2, alpha(k)=vec1(1)/vec2(1);for i=1:length(vec2), a3(i)=rtab(k,i+1)-alpha(k)*rtab(k+1,i+1); end if sum(abs(a3))==0 a3=polyder(vec2);info=[info,'All elements in row ',...int2str(k+2)' are zeros;'];elseif abs(a3(1)) rtab=[rtab;a3, zeros(1,nrT-length(a3))];vec1=vec2;vec2=a3;end 主程序: >>den=[1,2,8,12,20,16,16];>>[rtab,info]=routh(den)>>a=rtab(:,1)if all(a>0)disp('系統(tǒng)是穩(wěn)定的');else disp('系統(tǒng)是不穩(wěn)定的');end 5>>num=7*[1,5];den=conv([1,0,0],conv([1,10],[1,1]));>>[gm,pm,wg,wc]=margin(num,den)第十章 M文件: function varargout=rg_lead(ng0,dg0,s1)if nargout==1 ngv=polyval(ng0,s1);dgv=polyval(dg0,s1);g=ngv/dgv;thetal=pi-angle(g);zc=real(s1)-imag(s1)/tan(thetal);t=-1/zc;varargout{1}=[t,1];elseif nargout==2 ngv=polyval(ng0,s1);dgv=polyval(dg0,s1);g=ngv/dgv;theta=angle(g);phi=angle(s1);if theta>0 phi_c=pi-theta;end if theta<0;phi_c=-theta end theta_z=(phi+phi_c)/2;theta_p=(phi-phi_c)/2;z_c=real(s1)-imag(s1)/tan(theta_z);p_c=real(s1)-imag(s1)/tan(theta_p);nk=[1-z_c];varargout{2}=[1-p_c];kc=abs(p_c/z_c);if theta<0 kc=-kc end varargout{1}=kc*nk;else error('輸出變量數(shù)目不正確!');end 主程序: >> ng0=[1];dg0=10000*[1 0-1.1772];>>g0=tf(ng0,dg0);%滿足開環(huán)增益的為校正系統(tǒng)的傳遞函數(shù) >>s=kw2s(0.7,0.5)%期望的閉環(huán)主導(dǎo)極點 >>ngc=rg_lead(ng0,dg0,s);>>gc=tf(ngc,1)>>g0c=tf(g0*gc);>>rlocus(g0,g0c);>>b1=feedback(g0,1);%未校正系統(tǒng)的閉環(huán)傳遞函數(shù) >>b2=feedback(g0c,1);%校正后系統(tǒng)的閉環(huán)傳遞函數(shù) >>figure,step(b1,'r--',b2,'b');grid on %繪 2 M文件: function [ngc,dgc,k]=rg_lag(ng0,dg0,KK,s1,a)ngv=polyval(ng0,s1);dgv=polyval(dg0,s1);g=dgv/ngv;k=abs(g);%期望主導(dǎo)極點處的根軌跡增益 beta=k/KK;[kosi1,wn1]=s2kw(s1);zc=-wn1*sin(a*pi/180)/sin(pi-atan(sqrt(1-kosi1^2)/kosi1)-(a*pi/180));%利用正弦定理 pc=beta*zc;ngc=beta*[1,-zc];dgc=[1,-pc];主程序: >>KK=20;s1=-2+i*sqrt(6);a=1 >>ng0=[10];dg0=conv([1,0],[1,4]);>>g0=tf(ng0,dg0);>>[ngc,dgc,k]=rg_lag(ng0,dg0,KK,s1,a);>>gc=tf(ngc,dgc)>>g0c=tf(KK*g0*gc);>>b1=feedback(k*g0,1);>>b2=feedback(g0c,1);>>step(b1,'r--',b2,'b');grid on M文件: function [ngc,dgc,k]=rg_lag(ng0,dg0,KK,s1,a)ngv=polyval(ng0,s1);dgv=polyval(dg0,s1);g=dgv/ngv;k=abs(g);%期望主導(dǎo)極點處的根軌跡增益 beta=k/KK;[kosi1,wn1]=s2kw(s1);zc=-wn1*sin(a*pi/180)/sin(pi-atan(sqrt(1-kosi1^2)/kosi1)-(a*pi/180));%利用正弦定理 pc=beta*zc;ngc=beta*[1,-zc];dgc=[1,-pc];主程序: >>KK=128;s1=-2+i*2*sqrt(3);a=2 >>ng0=[10];dg0=conv([1,0],conv([1,2],[1,8]));>>g0=tf(ng0,dg0);>>[ngc,dgc,k]=rg_lag(ng0,dg0,KK,s1,a);>>gc=tf(ngc,dgc)>>g0c=tf(KK*g0*gc);>>rlocus(g0,g0c);>>b1=feedback(k*g0,1);>>b2=feedback(g0c,1);>>figure,step(b1,'r--',b2,'b');grid on 4 M文件: function [ngc,dgc]=lead4(ng0,dg0,KK,Pm,w)[mu,pu]=bode(KK*ng0,dg0,w);[gm,pm,wcg,wcp]=margin(mu,pu,w);alf=ceil(Pm-pm+5);phi=(alf)*pi/180;a=(1+sin(phi))/(1-sin(phi)), dbmu=20*log10(mu);mm=-10*log10(a);wgc=spline(dbmu,w,mm), T=1/(wgc*sqrt(a)), ngc=[a*T,1];dgc=[T,1];主程序: >>ng0=[1];dg0=conv([1,0,0],[1,5]);>>g0=tf(ng0,dg0);>>w=logspace(-3,3);>>KK=1;Pm=50;>>[ngc,dgc]=lead4(ng0,dg0,KK,Pm,w);>>gc=tf(ngc,dgc);g0c=tf(KK*g0*gc);>>bode(KK*g0,w);hold on,bode(g0c,w);grid on,hold off >>[gm,pm,wcg,wcp]=margin(g0c)>>Kg=20*log10(gm)>>g1=feedback(g0c,1);>>bode(g1),grid on, >>[mag,phase,w]=bode(g1);>>a=find(mag<=0.707*mag(1));>>wb=w(a(1)) >>max(mag)>>b=find(mag==max(mag))>>wr=w(b)5 M文件: function [ngc,dgc]=fg_lead_pm(ng0,dg0,Pm,w)[mu,pu]=bode(ng0,dg0,w);%計算原系統(tǒng)的對數(shù)頻率響應(yīng)數(shù)據(jù) [gm,pm,wcg,wcp]=margin(mu,pu,w);%求取原系統(tǒng)的相角裕度和剪切頻率 alf=ceil(Pm-pm+5);%計算控制器提供的最大超前角度,phi=(alf)*pi/180;%將最大超前角轉(zhuǎn)換為弧度單位 a=(1+sin(phi))/(1-sin(phi));%計算a值 dbmu=20*log10(mu);%系統(tǒng)的對數(shù)幅值 mm=-10*log10(a);%wm處的控制器對數(shù)幅值 wgc=spline(dbmu,w,mm);%差值求取wm,認(rèn)為wm=wc T=1/(wgc*sqrt(a));%計算T ngc=[a*T,1];dgc=[T,1];主程序: >>KK=40;Pm=50;>>ng0= KK *[1];dg0=conv([1,0],conv([1,1],[1,4]));>>g0=tf(ng0,dg0);>>w=logspace(-2,4);>>[ngc,dgc]=fg_lead_pm(ng0,dg0,Pm,w)>>gc=tf(ngc,dgc),g0c=tf(g0*gc);>>b1=feedback(g0,1);b2=feedback(g0c,1);>>step(b1,'r--', b2,'b');grid on >>figure, bode(g0,'r--',g0c,'b',w), grid on, >>[gm,pm,wcg,wcp]=margin(g0c), Km=20*log10(gm) MATLAB及控制系統(tǒng) 仿真實驗 班 級: 智能0702 姓 名: 劉保衛(wèi) 學(xué) 號: 06074053(18)實驗四 控制系統(tǒng)數(shù)學(xué)模型轉(zhuǎn)換及MATLAB實現(xiàn) 一、實驗?zāi)康?/p> 熟悉MATLAB 的實驗環(huán)境。 掌握MATLAB 建立系統(tǒng)數(shù)學(xué)模型的方法。 二、實驗內(nèi)容 (注:實驗報告只提交第2 題) 1、復(fù)習(xí)并驗證相關(guān)示例。(1)系統(tǒng)數(shù)學(xué)模型的建立 包括多項式模型(Transfer Function,TF),零極點增益模型(Zero-Pole,ZP),狀態(tài)空間模型 (State-space,SS);(2)模型間的相互轉(zhuǎn)換 系統(tǒng)多項式模型到零極點模型(tf2zp),零極點增益模型到多項式模型(zp2tf),狀態(tài)空間模 型與多項式模型和零極點模型之間的轉(zhuǎn)換(tf2ss,ss2tf,zp2ss…);(3)模型的連接 模型串聯(lián)(series),模型并聯(lián)(parallel),反饋連接(feedback) 2、用MATLAB 做如下練習(xí)。(1)用2 種方法建立系統(tǒng) 程序如下: %建立系統(tǒng)的多項式模型(傳遞函數(shù))%方法一,直接寫表達(dá)式 s=tf('s')Gs1=(s+2)/(s^2+5*s+10)%方法二,由分子分母構(gòu)造 num=[1 2];den=[1 5 10];Gs2=tf(num,den)figure pzmap(Gs1)figure pzmap(Gs1)grid on 運(yùn)行結(jié)果: 易知兩種方法結(jié)果一樣 的多項式模型。 Transfer function: s Transfer function: s + 2--------------s^2 + 5 s + 10 Transfer function: s + 2--------------s^2 + 5 s + 10(2)用2 種方法建立系統(tǒng)程序如下: %方法一 s=tf('s')Gs1=10*(s+1)/((s+1)*(s+5)*(s+10))% zpk模型 ZPK=zpk(Gs1) %方法二 % tf模型 num=[10 10];den=conv([1 1],conv([1 5],[1 10]));Gs2=tf(num,den)% zpk模型 ZPK=zpk(Gs2)figure pzmap(Gs1)figure pzmap(Gs1)grid on 運(yùn)行結(jié)果: 易知兩種方法結(jié)果一樣 的零極點模型和多項式模型。 Transfer function: s Transfer function: 10 s + 10------------------------s^3 + 16 s^2 + 65 s + 50 Zero/pole/gain: 10(s+1)------------------(s+10)(s+5)(s+1) Transfer function: 10 s + 10------------------------s^3 + 16 s^2 + 65 s + 50 Zero/pole/gain: 10(s+1)------------------(s+10)(s+5)(s+1) (3)如圖,已知G(s)和H(s)兩方框?qū)?yīng)的微分方程是: 且初始條件為零。試求傳遞函數(shù)C(s)/R(s)及E(s)/R(s)。 程序如下: %求微分方程的傳遞函數(shù)C(s)/R(s),E(s)/R(s)%求Gs=Cs/Rs n1=[20];d1=[6 10];Gs=tf(n1,d1) %求Hs=Bs/Cs n2=[10];d2=[20 5];Hs=tf(n2,d2) % C(s)/R(s)sys=feedback(Gs,Hs) % E(s)/R(s)=(Es/Cs)*(Cs/Rs)ER=sys/Gs 運(yùn)行結(jié)果: Transfer function: 20-------- % Gs=Cs/Rs 6 s + 10 Transfer function: 10-------- % Hs=Bs/Cs 20 s + 5 Transfer function: 400 s + 100--------------------- s^2 + 230 s + 250 Transfer function: 2400 s^2 + 4600 s + 1000------------------------ 2400 s^2 + 4600 s + 5000 % C(s)/R(s)% E(s)/R(s)=(Es/Cs)*(Cs/Rs)第四篇:控制系統(tǒng)的MATLAB仿真與設(shè)計課后答案
第五篇:Matlab 控制系統(tǒng) 傳遞函數(shù)模型