M实验指导书.docx
《M实验指导书.docx》由会员分享，可在线阅读，更多相关《M实验指导书.docx（17页珍藏版）》请在冰豆网上搜索。
M实验指导书
湖北工业大学理学院
《Matlab程序设计与应用》
实验指导书
陈华编
信息与计算科学系
二00八年一月
附录部分程序示例
实验一Matlab运算基础
1.z1=2*sin（85*pi/180）/（1+exp
（2））
x=[21+2*i;-0.455];z2=1/2*log（x+sqrt（1+x.^2））
a=-3.0:
0.1:
3.0;z3=1/2*（exp（0.3*a）-exp（-0.3*a））.*sin（a+0.3）+log（（0.3+a）/2）
t=0:
0.5:
2.5;z4=（t>=0&t<1）.*（t.^2）+（t>=1&t<2）.*（t.^2-1）+（t>=2&t<3）.*（t.^2-2*t+1）
2.A=[1234-4;34787;3657];B=[13-1;203;3-27];
J1=A+6*B,J2=A-B+eye（3）,E1=A*B,E2=A.*B,F1=A^3,F2=A.^3,G1=A/B,G2=B\A,H1=[A,B],H2=[A（[1,3],:
）;B^2]
3.A1=[12345;678910;1112131415;1617181920;2122232425];
B1=[3016;17-69;023-4;970;41311];
C=A1*B1,D=C（3:
5,[2,3]）
4.y=100:
999;k1=find（rem（y,21）==0）;k2=length（k1）
ch='ABC123Gabc';
k3=find（ch>='A'&ch<='Z'）;
ch（k3）=[]
实验二Matlab矩阵分析与处理
1.E=eye（3）;R=rand（3,2）;O=zeros（2,3）;S=diag（[1,1]）;
A=[E,R;O,S];
B1=A^2，B2=[E,R+R*S;O,S^2]
2.H=hilb（5）;P=pascal（5）;
Hh=det（H）,Hp=det（P）
Th=cond（H）,Tp=cond（P）
实验结果表明：
由于希尔伯特矩阵H的条件数大于帕斯卡矩阵P的条件数，故帕斯卡矩阵的性能更好。
3.A=rand（5）;
det（A）,trace（A）,rank（A）,norm（A）
4.A=[-29618;20512;-885];
[V,D]=eig（A）
实验结果表明：
D中对角线上的元素是A的特征值，V中的每列是对应的特征向量。
5.
（1）A=[1/21/31/4;1/31/41/5;1/41/51/6];
b=[0.950.670.52]';x=inv（A）*b
（2）b3的微小变化引起了解的较大变化。
（3）cond（A）
A的条件数为1.3533e+003，较大，故矩阵的A的性能不好，会引起解发生较大的变化。
6.A=[1/21/31/4;1/31/41/5;1/41/51/6];
sqrtm（A）％直接作用于矩阵A
sqrt（A）％作用于矩阵A的每个元素
实验三选择结构程序设计
1.x=input（'请输入x:
'）;
if（x<0&x~=-3）
y=x^2+x-6;
elseifx>=0&x<5&x~=2&x~3
y=x^2-5*x+6
else
y=x^2-x-1;
end
y
2.score=input（'请输入成绩：
'）;
ifscore<0|score>100
disp（'error'）;
else
switchfix（score/10）
case{9,10}
disp（'A'）
case{8}
disp（'B'）
case7
disp（'C'）
case6
disp（'D'）
otherwise
disp（'E'）
end
end
3.N=input（'请输入工号N：
'）;
T=input（'请输入员工的工时数T：
'）;
ifT>120
S=120*84+（T-120）*84*（1+0.15）;
elseifT<60
S=T*84-700;
else
S=T*84;
end
disp（['应发工资为:
',num2str（S）]）;
4.R1=fix（10+（99-10）*rand
（1））
R2=fix（10+（99-10）*rand
（1））
R1+R2
5.程序一：
A=[123456;789101112;131415161718;192021222324;...
252627282930];
n=input（'请输入n的值：
'）;
try
A（n,:
）
catch
A（end,:
）
end
lasterr
程序二：
A=[123456;789101112;131415161718;192021222324;...
252627282930];
n=input（'请输入n的值：
'）;
ifn>=1&n<=5
A（n,:
）
else
A（end,:
）
disp（'error'）
end
实验四循环结构程序设计
1.法一：
n=input（‘输入n的值：
’）；
s=0;
fori=1:
n
s=s+1/i.^2;
end
a=sqrt（6*s）
法二：
n=input（‘输入n的值：
’）；
i=1:
n;
f=1/i.^2;
s=sum（f）;
a=sqrt（6*s）
2.法一（最好）：
y=1;i=1;
whiley<3
i=i+1
y=y+1/（2*i-1）;
end
n=i-1
y=y-1/（2*i-1）
法二：
n=0;y=0;
fori=1:
100
y=y+1/（2*i-1）;
if（y<3）
n=n+1;
else
break;
end
end
n
y=y-1/（2*n+1）
3.formatlong
a=input（‘请输入a的值;’）;
b=input（‘请输入b的值;’）;
x
（1）=1.0;n=1;m=1;
whilem>10^-5&n<=500
x（n+1）=a./（b+x（n））;
m=abs（x（n+1）-x（n））;
n=n+1;
end
x=x（n）
r1=（-b+sqrt（b^2+4*a））/2;
r2=（-b+sqrt（b^2+4*a））/2;
4.方法一：
f
（1）=1;f
（2）=0;f（3）=1;n=4;
whilen>3&n<=100
f（n）=f（n-1）-2*f（n-2）+f（n-3）;
n=n+1;
end
y1=max（f）,y2=min（f）,y3=sum（f）
n1=length（find（f>0））,n2=length（find（f==0））,n3=length（find（f<0））
方法二：
先建立函数文件myfun1.m
functiony=myfun1（n）
ifn==1
y=1;
elseifn==2
y=0;
elseifn==3
y=1;
else
y=myfun1（n-1）-2*myfun1（n-2）+myfun1（n-3）;％当n较大时，递归调用速度较慢，耐心等待
end
再建立命令文件ex.m，并在此调用已定义的函数myfun1
max=1;
min=1;
sum=0;
p=0;
n=0;
zero=0;
fori=1:
100
y=myfun1（i）;
sum=sum+y;
ify>0
p=p+1;
elseify==0
zero=zero+1;
else
n=n+1;
end
ify>max
max=y;
end
ifymin=y;
end
end
max,min,sum,p,n,zero
方法三：
f
（1）=1;f
（2）=0;f（3）=1;
min=0;max=1;sum=2;
n1=2;n2=1;n3=0;
forn=4:
100
f（n）=f（n-1）-2*f（n-2）+f（n-3）;
sum=sum+f（n）;
iff（n）>0
n1=n1+1;
elseiff（n）==0
n2=n2+1;
else
n3=n3+1;
end
iff（n）min=f（n）;
end
iff（n）>max
max=f（n）;
end
end
sum
max
min
n1
n2
n3
5.方法一：
i=2:
49;
y=i.*（i+1）-1;
forj=2:
sqrt（y）
n=find（rem（y,j）==0&y~=j）;
y（n）=[];
end
m=length（y）
s=sum（y）
方法二：
先建立函数文件pp.m
functiony=pp（n）
fori=2:
sqrt（n）
if（rem（n,i）==0）
y=0;
else
y=1;
end
end
再建立命令文件ex.m
num=0;
sum=0;
fori=2:
49
y=i.*（i+1）-1;
if（pp（y））
num=num+1;
sum=sum+y;
end
end
num
sum
实验五循环结构程序设计
1.function[e,l,s,c]=f（x）
e=exp（x）;l=log（x）;s=sin（x）;c=cos（x）;
2.function[a1,a2,N1,N2]=f（m1,m2,a）
g=9.8;
A=[m1*cos（a）-m1-sin（a）0;m1*sin（a）0cos（a）0;...
0m2-sin（a）0;00-cos（a）1];
B=[0m1*g0m2*g];
C=inv（A）*B';
a1=C
（1）;a2=C
（2）;N1=C（3）;N2=C（4）;
3.functiony=prime（x）
i=2:
x/2;
ifrem（x,i）~=0
y=1;
else
y=0;
end
保存为prime.m文件
法一：
j=1;
fori=10:
99
m=fix（i/10）;
n=rem（i,10）;
b=n*10+m;
ifprime（i）==1&prime（b）==1
y（j）=i;
j=j+1;
end
end
y
法二：
j=1;
fori=10:
99
m=fix（i/10）;
n=rem（i,10）;
b=n*10+m;
ifprime（i）==1&prime（b）==1
y（j）=i;
else
y（j）=0;
end
j=j+1;
end
y（find（y==0））=[];
y
4.functiony=f（x）
y=1./（（x-2）.^2+0.1）+1./（（x-3）.^4+0.01）
5.functiony=f1（n）
y=n+10*log（n^2+5）;
保存为f1.m
y=f1（40）/（f1（30）+f1（20））
functiony=f2（n）
i=1:
n;
y=sum（n*（n+1））;
保存为f2.m
y=f2（40）/（f2（30）+f2（20））
实验六高层绘图操作
1.x=linspace（0,2*pi,101）;
y=（0.5+3*sin（x）./（1+cos（x）））.*cos（x）;
plot（x,y）
2.
（1）x=linspace（0,2*pi,101）;
y1=x.^2;y2=cos（2*x）;y3=y2.*y2;
plot（x,y1,'r:
',x,y2,'y-',x,y3,'g-.'）
（2）x=linspace（0,2*pi,101）;
y1=x.^2;y2=cos（2*x）;y3=y2.*y2;
subplot（1,3,1）
plot（x,y1）
subplot（1,3,2）
plot（x,y2）
subplot（1,3,3）
plot（x,y3）
（3）x=linspace（0,2*pi,101）;
y1=x.^2;y2=cos（2*x）;y3=y2.*y2;
subplot（1,3,1）
bar（x,y1）
subplot（1,3,2）
stairs（x,y2）
subplot（1,3,3）
fill（x,y3）
3.x=linspace（-5,5,10000）;
y=（x<=0）.*（x+sqrt（pi））/exp
（2）+（x>0）.*1/2.*log（x+sqrt（1+x.^2））;
plot（x,y）
4.subplot（2,2,1）
theta=0:
0.01:
2*pi;
rho1=2*sin（2+2*theta）;
polar（theta,rho1）
subplot（2,2,2）
rho2=4*sin（2+2*theta）;
polar（theta,rho2）
subplot（2,2,3）
rho3=2*sin（4+2*theta）;
polar（theta,rho3）
subplot（2,2,4）
rho4=2*sin（2+4*theta）;
polar（theta,rho4）
5.x=linspace（-5,5,21）;
y=linspace（0,10,31）;
[x,y]=meshgrid（x,y）;
z=cos（x）.*cos（y）.*exp（-sqrt（x.^2+y.^2）/4）;
subplot（2,1,1）
surf（x,y,z）
subplot（2,1,2）,surfc（x,y,z）
6.
ezsurf（'cos（s）.*cos（t）','cos（s）.*sin（t）','sin（s）',[0,pi/2,0,3*pi/2]）
Shadinginterp
实验七低层绘图操作
1.figure（'color','r','menubar','figure','windowbuttondownfcn','disp（''LeftButtonPressed''）'）
2.x=linspace（-2,2,100）;y=x.^2.*exp（2*x）;h=plot（x,y）;
set（h,'color','r','linestyle',':
','linewidth',3）
text（1,1.^2.*exp（2*1）,'y=x^2e^{2x}'）
3.x=linspace（-2,2,10）;
[x,t]=meshgrid（x）;
z=10*exp（-0.01*x）.*sin（2000*pi*t-0.2*x+pi）;
h=surf（x,t,z）;
set（h,'edgecolor','r','linestyle',':
','linewidth',3）text（-1,-0.5,-14,'v（x,t）=10e^{-0.01x}sin（2000\pit-0.2x+\pi）'）
4.x=linspace（0,2*pi,40）;
axes（'position',[0.2,0.2,0.2,0.7]）
plot（x,sin（x））;
title（'sin（x）'）;
axes（'position',[0.4,0.5,0.2,0.1]）
plot（x,cos（x））;
title（'cos（x）'）;
axes（'position',[0.5,0.2,0.3,0.2]）
plot（x,sin（x）./（cos（x）+eps））;
title（'tan（x）'）;
axes（'position',[0.45,0.65,0.4,0.25]）
plot（x,cos（x）./（sin（x）+eps））;
title（'cot（x）'）
5.cylinder（3）;
lightingphong;materialshiny;
实验八数据处理与多项式计算
1.A=rand（1,30000）;
（1）m=mean（A）,s=std（A）,
（2）high=max（A）,low=min（A）
（3）n=length（find（A>0.5））;p=n/30000
2.score=70+7*randn（100,5）;
[high,num]=max（score）,[low,num]=min（score）
average=mean（score）,stdn=std（score）
scores=score';total=sum（scores）;
[highs,nums]=max（total）,[lows,nums]=min（total）
[zcj,xsxh]=sort（-total）;zcj=-zcj,xsxh=xsxh
3.h=6:
2:
18;
t1=[18202225302824];
t2=[15192428343230];
x=6.30:
2:
17.30;
y1=interp1（h,t1,x,'spline'）
y2=interp1（h,t2,x,'spline'）
4.x=1:
10:
101;
y=log10（x）;
p=polyfit（x,y,5）
y1=polyval（p,x）
plot（x,y,':
o',x,y1,'-*'）
5.p1=[12405];p2=[1,2];p3=[123];
（1）p=p1+[0,conv（p2,p3）]%两个多项式系数向量相加，须是同维向量
（2）x=roots（p）
（3）A=[-11.2-1.4;0.7523.5;052.5];
y1=polyval（p,A）
（4）y2=polyvalm（p,A）
实验九数值微积分与方程数值求解
实验十符号计算基础与符号微积分
实验十一级数与方程符号求解
实验十二菜单与对话框设计
升级会员 