实验二.docx
《实验二.docx》由会员分享，可在线阅读，更多相关《实验二.docx（28页珍藏版）》请在冰豆网上搜索。
实验二
实验二
1：
建立m文件，求z=
所表示的三维曲面（图1.3-4）。
x,y的取值范围是[8,8]。
输入代码为：
clear;x=-8:
0.5:
8;
y=x';
X=ones（size（y））*x;
Y=y*ones（size（x））;
R=sqrt（X.^2+Y.^2）+eps;%<5>
Z=sin（R）./R;%<6>
surf（X,Y,Z）;%
colormap（cool）%
xlabel（'x'）,ylabel（'y'）,zlabel（'z'）
运行结果为：
2：
多项式四则运算。
建立m文件，输入以下代码并执行。
format
disp（''）
disp（'多项式四则运算'）,disp（''）,pause
disp（'相乘:
a=[2,4,6,8],b=[3,6,9],c=conv（a,b）'）,pause;disp（''）
a=[2,4,6,8],b=[3,6,9],c=conv（a,b）,pause
disp（'相加:
d=a+[0,b]'）,pause,d=a+[0,b],pause;disp（''）
disp（'相除:
[q,r]=deconv（c,a）'）,pause,[q,r]=deconv（c,a）,pause;disp（''）
disp（'a1=a+1'）,pause,a1=a+1,pause;disp（''）
disp（'[q1,r1]=deconv（c,a1）'）,pause,[q1,r1]=deconv（c,a1）,pause;disp（''）
disp（'conv（q1,a1）+r1'）,pause,conv（q1,a1）+r1,pause;disp（''）
disp（'多项式求导数:
e=polyder（c）'）,pause,e=polyder（c）,pause;disp（''）
disp（'多项式方程求根'）,disp（''）,pause;disp（''）
disp（'ra=roots（a）,rb=roots（b）;rc=roots（c）;[[ra;rb],rc]'）,pause;disp（''）
formatcompact,ra=roots（a）,rb=roots（b）,rc=roots（c）,pause,format
disp（''）
disp（'多项式求值:
（polyval）'）,pause;
disp（''）
disp（'线性间隔:
w=linspace（0,10）;'）,pause,w=linspace（0,10）;disp（''）
disp（'A=polyval（a,j*w）;plot（w,abs（A））'）,pause;disp（''）
A=polyval（a,j*w）;h=plot（w,abs（A））;set（h,'linewidth',2）,pause
disp（'B=polyval（b,j*w）;plot（w,abs（B））'）,pause;
disp（''）
B=polyval（b,j*w）;plot（w,abs（B））,pause
disp（'subplot（2,2,1）;plot（w,abs（B./A））,subplot（2,2,3）;plot（w,angle（B./A））'）,pause;disp（''）
subplot（2,2,1）;h=plot（w,abs（B./A））;set（h,'linewidth',2）
subplot（2,2,3）;h=plot（w,angle（B./A））;,set（h,'linewidth',2）,pause
disp（'对数等间隔:
w1=logspace（a,b,n）%在10^a到10^b之间按等比分为n份'）,pause;disp（''）
disp（''）
disp（'w1=logspace（-1,1）'）;pause,w1=logspace（-1,1）;pause
disp（''）
disp（'F=polyval（b,j*w1）./polyval（a,j*w1）;'）,pause;
disp（''）
disp（'subplot（2,2,2）,loglog（w1,abs（F））'）,pause;
disp（''）
F=polyval（b,j*w1）./polyval（a,j*w1）;subplot（2,2,2）
loglog（w1,abs（F）,'linewidth',2）,pause
disp（''）
disp（'subplot（2,2,4）;semilogx（w1,angle（F））'）,pause;
disp（''）
subplot（2,2,4）;semilogx（w1,angle（F）,'linewidth',2）,pause
disp（''）
disp（'曲线拟合:
'）,pause;
disp（''）
disp（'原始数据:
'）;
disp（''）
x=0:
0.1:
1;y=[-.447,1.978,3.28,6.16,7.08,7.34,7.66,9.56,9.48,9.30,11.2];
subplot（2,3,1）,plot（x,y,'o'）,pause
disp（'a1=polyfit（x,y,1）;xi=linspace（0,1）;y1=polyval（a1,xi）;'）,pause;
disp（''）
a1=polyfit（x,y,1）;xi=linspace（0,1）;y1=polyval（a1,xi）;pause
disp（'plot（x,y,''o'',xi,y1,''b''）,pause'）
disp（''）
subplot（2,3,2）,plot（x,y,'o',xi,y1,'b'）,pause
disp（'a2=polyfit（x,y,2）;y2=polyval（a2,xi）;plot（x,y,''o'',xi,y2,''m''）'）,pause
disp（''）
a2=polyfit（x,y,2）;y2=polyval（a2,xi）;
subplot（2,3,3）,plot（x,y,'o',xi,y2,'r'）,pause
disp（'a3=polyfit（x,y,3）;y3=polyval（a3,xi）;plot（x,y,''o'',xi,y3,''r''）'）,pause
disp（''）
a3=polyfit（x,y,3）;y3=polyval（a3,xi）;
subplot（2,3,4）,plot（x,y,'o',xi,y3,'r'）,pause
disp（'a9=polyfit（x,y,9）;y9=polyval（a9,xi）;plot（x,y,''o'',xi,y9,''b''）'）,pause
disp（''）
a9=polyfit（x,y,9）;y9=polyval（a9,xi）;
subplot（2,3,5）,plot（x,y,'o',xi,y9,'b'）,pause
disp（'a10=polyfit（x,y,10）;y10=polyval（a10,xi）;plot（x,y,''o'',xi,y10,''b''）'）,pause
disp（''）
a10=polyfit（x,y,10）;y10=polyval（a10,xi）;
subplot（2,3,6）,plot（x,y,'o',xi,y10,'b'）,pause
disp（'一维插值:
interp1（x,y,xi,''method''）'）,pause;disp（''）
disp（'z1=interp1（x,y,xi）;plot（x,y,''o'',xi,z1）'）,pause
disp（''）
z1=interp1（x,y,xi）;subplot（1,2,1）,plot（x,y,'o',xi,z1）,pause
disp（'z2=interp1（x,y,xi,''spline''）;plot（x,y,''o'',xi,z2,''b''）'）,pause
disp（''）
z2=interp1（x,y,xi,'spline'）;subplot（1,2,2）,plot（x,y,'o',xi,z2,'b'）,pause
disp（'二维插值:
interp2（x,y,z,xi,yi,''method''）'）,pause;
disp（''）
disp（'width=1:
5;depth=1:
3;temps=[8281808284;7963616581;8484828586];'）,pause
disp（''）
width=1:
5;depth=1:
3;temps=[8281808284;7963616581;8484828586];pause
%disp（'沿深度的中心线的温度分布:
'）,pause,echoon
%wi=1:
0.2:
5;d=2;t1=interp2（width,depth,temps,wi,d）;,pause
%t2=interp2（width,depth,temps,wi,d,'cubic'）;plot（wi,t1,'--',wi,t2,'g'）,pause
disp（'沿宽度和深度二维分布:
'）,pause,disp（''）
disp（'di=1:
0.2:
3;wi=1:
0.2:
5;tc=interp2（width,depth,temps,wi,di'',''cubic''）;'）,pause
disp（''）
di=1:
0.2:
3;wi=1:
0.2:
5;tc=interp2（width,depth,temps,wi,di','cubic'）;pause
subplot（1,1,1）
disp（''）
disp（'h=mesh（wi,di,tc）,set（h,''linewidth'',1.5'）
disp（''）
h=mesh（wi,di,tc）;set（h,'linewidth',1.5）,pause
disp（'线性微分方程的解:
'）,disp（''）,pause
disp（'脉冲响应:
'）,pause
disp（''）
disp（'a=[1,5,4,7];b=[3,0.5,4];[r,p,k]=residue（b,a）'）,pause
disp（''）
a=[1,5,4,7];b=[3,0.5,4];[r,p,k]=residue（b,a）,pause
disp（'t=0:
0.2:
10;yi=r
（1）*exp（p
（1）*t）+r
（2）*exp（p
（2）*t）+r（3）*exp（p（3）*t）;plot（t,yi）'）,pause
disp（''）
t=0:
0.2:
10;yi=r
（1）*exp（p
（1）*t）+r
（2）*exp（p
（2）*t）+r（3）*exp（p（3）*t）;
subplot（1,2,1）,h=plot（t,yi）;set（h,'linewidth',2）,pause
disp（''）
disp（'阶跃响应:
'）
disp（''）
disp（'a=[1,5,4,7,0];b=[3,0.5,4];[r,p,k]=residue（b,a）'）,pause
disp（''）
a=[1,5,4,7,0];b=[3,0.5,4];[r,p,k]=residue（b,a）,pause
disp（''）
disp（'ys=r
（1）*exp（p
（1）*t）+r
（2）*exp（p
（2）*t）+r（3）*exp（p（3）*t）+r（4）;plot（t,ys）'）,pause
disp（''）
ys=r
（1）*exp（p
（1）*t）+r
（2）*exp（p
（2）*t）+r（3）*exp（p（3）*t）+r（4）;
subplot（1,2,2）,h=plot（t,ys）;set（h,'linewidth',2）,pause
disp（''）
disp（'在控制工具箱中:
'）
disp（''）
disp（'yi=impulse（b,a）,及ys=step（b,a）'）,pause
echooff,clf
多项式四则运算
相乘:
a=[2,4,6,8],b=[3,6,9],c=conv（a,b）
a=
2468
b=
369
c=
624609610272
相加:
d=a+[0,b]
d=
271217
相除:
[q,r]=deconv（c,a）
q=
369
r=
000000
a1=a+1
a1=
3579
[q1,r1]=deconv（c,a1）
q1=
2.00004.66677.5556
r1=
0007.55567.11114.0000
conv（q1,a1）+r1
ans=
624609610272
多项式求导数:
e=polyder（c）
e=
3096180192102
多项式方程求根
ra=roots（a）,rb=roots（b）;rc=roots（c）;[[ra;rb],rc]
ra=
-1.6506
-0.1747+1.5469i
-0.1747-1.5469i
rb=
-1.0000+1.4142i
-1.0000-1.4142i
rc=
-1.6506
-1.0000+1.4142i
-1.0000-1.4142i
-0.1747+1.5469i
-0.1747-1.5469i
多项式求值:
（polyval）
线性间隔:
w=linspace（0,10）;
A=polyval（a,j*w）;plot（w,abs（A））
B=polyval（b,j*w）;plot（w,abs（B））
subplot（2,2,1）;plot（w,abs（B./A））,subplot（2,2,3）;plot（w,angle（B./A））
对数等间隔:
w1=logspace（a,b,n）%在10^a到10^b之间按等比分为n份
w1=logspace（-1,1）
F=polyval（b,j*w1）./polyval（a,j*w1）;
subplot（2,2,2）,loglog（w1,abs（F））
subplot（2,2,4）;semilogx（w1,angle（F））
曲线拟合:
原始数据:
a1=polyfit（x,y,1）;xi=linspace（0,1）;y1=polyval（a1,xi）;
plot（x,y,'o',xi,y1,'b'）,pause
a2=polyfit（x,y,2）;y2=polyval（a2,xi）;plot（x,y,'o',xi,y2,'m'）
a3=polyfit（x,y,3）;y3=polyval（a3,xi）;plot（x,y,'o',xi,y3,'r'）
a9=polyfit（x,y,9）;y9=polyval（a9,xi）;plot（x,y,'o',xi,y9,'b'）
a10=polyfit（x,y,10）;y10=polyval（a10,xi）;plot（x,y,'o',xi,y10,'b'）
一维插值:
interp1（x,y,xi,'method'）
z1=interp1（x,y,xi）;plot（x,y,'o',xi,z1）
z2=interp1（x,y,xi,'spline'）;plot（x,y,'o',xi,z2,'b'）
二维插值:
interp2（x,y,z,xi,yi,'method'）
width=1:
5;depth=1:
3;temps=[8281808284;7963616581;8484828586];
沿宽度和深度二维分布:
di=1:
0.2:
3;wi=1:
0.2:
5;tc=interp2（width,depth,temps,wi,di','cubic'）;
h=mesh（wi,di,tc）,set（h,'linewidth',1.5
线性微分方程的解:
脉冲响应:
a=[1,5,4,7];b=[3,0.5,4];[r,p,k]=residue（b,a）
r=
3.2288
-0.1144+0.0730i
-0.1144-0.0730i
p=
-4.4548
-0.2726+1.2235i
-0.2726-1.2235i
k=
[]
t=0:
0.2:
10;yi=r
（1）*exp（p
（1）*t）+r
（2）*exp（p
（2）*t）+r（3）*exp（p（3）*t）;plot（t,yi）
阶跃响应:
a=[1,5,4,7,0];b=[3,0.5,4];[r,p,k]=residue（b,a）
r=
-0.7248
0.0767+0.0764i
0.0767-0.0764i
0.5714
p=
-4.4548
-0.2726+1.2235i
-0.2726-1.2235i
0
k=
[]
ys=r
（1）*exp（p
（1）*t）+r
（2）*exp（p
（2）*t）+r（3）*exp（p（3）*t）+r（4）;plot（t,ys）
在控制工具箱中:
yi=impulse（b,a）,及ys=step（b,a）
3：
三维绘图和屏幕控制
clf
disp（'三维曲线'）,echoon
pause,z=0:
0.1:
4*pi;x=cos（z）;y=sin（z）;plot3（x,y,z）;shg
pause,x=-8:
0.5:
8;y=x';X=ones（size（y））*x;Y=y*ones（size（x））;
pause,subplot（1,2,1）,plot（X,Y）,
pause,subplot（1,2,2）,plot（Y,X）,
pause,subplot（1,1,1）,R=sqrt（X.^2+Y.*Y）;z=sin（R）./R;mesh（z）
pause,R=sqrt（X.^2+Y.*Y）+eps;z=sin（R）./R;figure
（1）,mesh（z）
pause,R=abs（X）+abs（Y）+eps;z1=sin（R）./R;figure
（2）,mesh（z1）
pause,text（10,30,1.1,'R=abs（X）+abs（Y）,z=sin（R）./R'）,view（20,0）,
pause,disp（'屏幕分割'）
pause,disp（'坐标设定'）
pause,subplot（2,2,1）,R=sqrt（X.^2+Y.*Y）;z=sin（R）./R;meshc（z）
pause,title（'meshc（z）,shadingflat'）,shadingflat
pause,subplot（2,2,2）,R=sqrt（X.^2+Y.*Y）+eps;z=sin（R）./R;meshz（z）
pause,title（'meshz（z）,shadinginterp'）,shadinginterp
pause,subplot（2,2,3）,R=abs（X）+abs（Y）+eps;z1=sin（R）./R;surfc（z1）
pause,title（'surfc（z1）,shadingflat'）,shadingflat,%colormap（gray）
pause,subplot（2,2,4）;surfl（z1）,view（20,0）
pause,title（'surfl（z1）,view（20,0）'）
4：
直流电路的稳态计算
clear,formatcompact
R1=2;R2=4;R3=12;R4=4;R5=12;R6=4;R7=2;%为给定元件赋值
%解问题
（1）
display（'解问题
（1）'）
a11=R1+R2+R3;a12=-R3;a13=0;%将系数矩阵各元素赋值
a21=-R3;a22=R3+R4+R5;a23=-R5;
a31=0;a32=-R5;a33=R5+R6+R7;
b1=1;b2=0;b3=0;
%输入解
（1）的已知条件
us=input（'us='）;
%列出系数矩阵A
A=[a11,a12,a13;a21,a22,a23;a31,a32,a33]
B=[b1;0;0];I=A\B*us;%I=[ia;ib;ic]
ia=I
（1）;ib=I
（2）;ic=I（3）;
%解出所需变量
i3=ia-ib,u4=R4*ib,u7=R7*ic
%利用电路的线性性质及问题
（1）的解
display（'解问题
（2）'）
u42=input（'给定u42='）;
%由问题
（1）得出待求量与us的比例系数
k1=i3/us;k2=u4/us;k3=u7/us;
us2=u42/k2,i32=k1/k2*u42,u72=k3/k2*u42%按比例方法求出所需变量
解问题
（1）
us=10
A=
18-120
-1228-12
0-1218
i3=
0.3704
u4=
2.2222
u7=
0.7407
解问题
（2）
给定u42=6
us2=
27.0000
i32=
1.0000
u72=
2
5:
一阶电路开关的暂态计算
clear,formatcompact
%figure
（1）,fg521,pause
r1=3;us=18;is=3;r2=12;r3=6;C=1;%给出原始数据
%算出初值ir20及uc0
uc0=-12;ir20=uc0/r2;ir30=uc0/r3;
ic0=is-ir20-ir30;
%算出终值ir2f及ucf
ir2f=is*r3/（r2+r3）;
ir3f=is*r2/（r2+r3）;
ucf=ir2f*r2;icf=0;
%注意时间数组的设置，在t=0及10附近设两个点
t=[-2,-1,0-eps,0+eps,1:
9,10-eps,10+eps,11:
20];%
%找出时间与数组下标的关系，t=10+eps对应下标15
figure
（2）,plot（t）,grid,pause
%t<0时的值
uc（1:
3）=-12;ir2（1:
3）=3;
%求充电时常数
T=r2*r3/（r2+r3）*C;
uc（4:
14）=ucf+（uc0-ucf）*exp（-t（4:
14）/T）;%
%用三要素法求输出
ir2（4:
14）=ir2f+（ir20-ir2f）*exp（-t（4:
14）/T）;
%求t=10+eps时的各初值
uc（15）=uc（14）;ir2（15）=is;
%求uc和ir2在新区间终值ucf2和ir2f
uc
升级会员 