MATLAB04数值计算.docx
《MATLAB04数值计算.docx》由会员分享，可在线阅读，更多相关《MATLAB04数值计算.docx（14页珍藏版）》请在冰豆网上搜索。
MATLAB04数值计算
数值计算
一、概念
1、两种计算方式
（1）数值计算
（2）符号计算
2、数值问题
（1）0.1+0.2==0.3？
（2）不能准确描述邻域
（3）没有数值求极限函数
（4）没有数值求导数函数
3、
二、知识点
1、HelpBrowser
（1）operator
（2）datatype
2、数据可视化
（1）plot
（2）subplot
3、匿名函数
（1）@（x）（sin（x）+cos（x））
（2）FH=@（x）（sin（x）+cos（x））
4、浮点数体系（IEEE754）
（1）一般数
（2）±0
（3）±Inf
（4）NaN
5、微积分、矩阵和代数方程、多项式运算、拟合和最小二乘、插值和样条、Fourier分析、微分方程
6、【选学】随机数发生和统计、卷积、优化计算
三、练习题
1、【4.2.1】数值极限
【例4.2-1】设
，
，由分析知
。
试用数值法求这两个函数的极限
，
。
（1）
k=logspace（0,14,15）;
x=k*eps;
f1=（1-cos（2.*x））./x./sin（x）/2;
f2=sin（x）./x;
formatshorte
disp（[blanks（7）,'x',blanks（11）,'f1',blanks（11）,'f2']）
disp（[x',f1',f2']）
clf
n1=sum（f1<10^-10）;
subplot（2,1,1）
semilogx（x（1:
n1）,f1（1:
n1）,'ro','LineWidth',3）
holdon
semilogx（x（n1+1:
end）,f1（n1+1:
end）,'b.-','LineWidth',3）
holdoff
text（10^-14,0.2,'错误区段'）,text（10^-6,0.8,'正确区段'）
axis（[10^-20,1,0,1.1]）
gridon
title（'f2'）
subplot（2,1,2）
semilogx（x,f2,'b.-','LineWidth',3）
axis（[10^-20,1,0,1.1]）
title（'f1'）
xlabel（'x'）
gridon
shg
（2）
symst
fs1=（1-cos（2*t））/（t*sin（t））/2;
fs2=sin（t）/t;
Ls1=limit（fs1,t,0）
Ls2=limit（fs2,t,0）
2、【4.2.2】数值差分
【例4.2-2】在区间[0,3]内，计算函数
的两点前向差分和内点中心差分，并绘制相应的图形。
clf
h=0.5;
x=0:
h:
3;
f=4-（x-2）.^2;
df=diff（f）;
gf=gradient（f）;
subplot（2,1,1）
holdon
plot（x,f,'r.-'）
x1=x（1:
end-1）;
plot（x1,df,'bo-'）
legend（'f','df'）
text（1,3.4,'f3'）,text（1.5,4.2,'f4'）
text（1.05,1.1,'df3=f4-f3'）
axis（[-1,4,-2,5]）
xlabel（'x'）
boxon
holdoff
subplot（2,1,2）
holdon
plot（x,f,'r.-'）
plot（x1,df,'bo-'）
plot（x,gf,'k*-'）
legend（'f','df','gf'）
text（0,-0.5,'f1'）,text（2.9,2.4,'f7'）
text（-0.45,2.2,'gf1=f2-f1'）,text（3,-0.3,'gf7=f7-f6'）
text（0.5,2.3,'f2'）,text（1.5,4.2,'f4'）
text（1.05,1.3,'gf3=（f4-f2）/2'）
axis（[-1,4,-2,5]）
xlabel（'x'）
boxon
holdoff
shg
3、【4.2.3】数值积分
【例4.2-3】以
为被积函数，在
区间，使用quad指令计算积分
。
被积函数中的
是可调参数。
（1）
clear
formatlong
a=1.5;
y1=@（t）cos（1./（a*t.^2））;%<4>
S1=quad（y1,0.2,0.3）
4、【例4.2-4】求曲线
与
所夹区域的面积。
（1）
formatlong
f1=@（x）sqrt（x）;
f2=@（x）x.^3;%<3>
s1=quad（f1,0,1）;
s2=quad（f2,0,1）;
S=s1-s2
（2）利用quadv指令计算
f=@（x）[sqrt（x）,x.^3];
s=quadv（f,0,1）;
Sv=s
（1）-s
（2）
5、【例4.2-5】求积分
，其中
。
（1）
formatlong
x1=linspace（0,pi,100）;
x1
（1）=realmin;%<3>
y1=sin（x1）./x1;
S_t1=trapz（x1,y1）,
（2）
x2=linspace（0,pi,200）;%<6>
x2
（1）=realmin;%<7>
y2=sin（x2）./x2;
S_t2=trapz（x2,y2）
E_t1=abs（S_t2-S_t1）%<9>
（3）
ys=@（x）（sin（x）./x）;%<10>
S_simposon=quad（ys,0,pi,1e-7）
（4）
symsx
S_sym=vpa（int（sin（x）/x,0,pi））
6、【4.2.4】多重数值积分
【例4.2-8】计算
和
。
1）
symsxy
ssx01=vpa（int（int（x^y,x,0,1）,y,1,2））
ssx12=vpa（int（int（x^y,x,1,2）,y,0,1））
（2）
zz=inline（'x.^y','x','y'）;
nsx01=dblquad（zz,0,1,1,2）
nsx12=dblquad（zz,1,2,0,1）
7、【4.3.1】矩阵运算和特征参数
【例4.3-1】观察矩阵的转置操作和数组转置操作的差别。
formatrat
A=magic
（2）+j*pascal
（2）
A1=A'
A2=A.'
B1=A*A'
B2=A.*A'
C1=A*A.'
C2=A.*A.'
8、【例4.3-2】矩阵标量特征参数计算示例。
A=reshape（1:
9,3,3）
r=rank（A）
d3=det（A）
d2=det（A（1:
2,1:
2））
t=trace（A）
9、【例4.4-1】简单实阵的特征值分解及特征值阵的“复/实转换”。
A=[1,-3;2,2/3];
d=eig（A）
[V,D]=eig（A）
10、【4.5.1】求解线性方程的相关指令
【例4.5-1】“求逆”法和“左除”法解恰定方程的性能对比。
A=gallery（'randsvd',100,2e13,2）;
x=ones（100,1）;
b=A*x;
xi=inv（A）*b
tic;xd=A\b
11、【4.5.4】一般代数方程的解。
【例4.5-3】求
的零点。
（1）
S=solve（'sin（t）^2*exp（-0.1*t）-0.5*abs（t）','t'）
（2）
YC=@（t）（sin（t）.^2.*exp（-0.1*t）-0.5*abs（t））;
（3）
t=-10:
0.01:
10;
Y=YC（t）;
clf,
plot（t,Y,'r'）;
holdon
plot（t,zeros（size（t））,'k'）;
xlabel（'t'）;ylabel（'y（t）'）
holdoff
图4.5-1函数零点分布观察图
（4）
zoomon
[tt,yy]=ginput（5）;zoomoff
图4.5-2局部放大和利用鼠标取值图
tt
（5）
[t4,y4]=fzero（YC,0.1）
12、【4.7.1】多项式的运算函数
【例4.7-1】求
的“商”及“余”多项式。
（1）
formatrat
p1=conv（[1,0,2],conv（[1,4],[1,1]））;%<2>
p2=[1011];
[q,r]=deconv（p1,p2）;
cq='商多项式为';cr='余多项式为';
disp（[cq,poly2str（q,'s'）]）,disp（[cr,poly2str（r,'s'）]）%<6>
（2）
qp2=conv（q,p2）;
pp1=qp2+[zeros（1,length（qp2）-length（r））,r];%<8>
pp1==p1
13、【例4.7-2】矩阵和特征多项式，特征值和多项式根。
（1）
A=[111213;141516;171819];
PA=poly（A）%<2>
PPA=poly2str（PA,'s'）
（2）
s=eig（A）
r=roots（PA）
14、【例4.7-3】构造指定特征根的多项式。
R=[-0.5,-0.3+0.4*i,-0.3-0.4*i];
P=poly（R）
PPR=poly2str（PR,'x'）
15、【例4.7-4】多项式求值指令polyval与polyvalm的本质差别。
（1）
clear
p=[1,2,3];
disp（['给定多项式为',poly2str（p,'x'）]）
X=[1,2;3,4]
（2）
va=X.^2+2*X+3%<5>
Va=polyval（p,X）
（3）
vm=X^2+2*X+3*eye
（2）%<7>
Vm=polyvalm（p,X）
（4）
cp=poly（X）;
disp（['矩阵X的特征多项式为',poly2str（cp,'x'）]）
cpXa=polyval（cp,X）
cpX=polyvalm（cp,X）
四、函数列表
函数名
功能
plot（x,y,s,p）
画图
semilogx（x,y,s,p）
对x轴数据取以10为底的对数，再plot
subplot（x,y,z）
分割图形窗口
diff（F）
两点前向差分
gradient（F）
内点中心差分
trapz（x,y）
梯形法求积分
quad（fun,a,b,tol）
递推自适应辛普森法求误差不大于tol的函数fun的定积分近似值
quadv（fun,a,b,tol）
求被积函数数组的积分
五、习题
1、设
。
（1）求
的实数根。
（2）求
，
。
（3）求
及
使得
。
升级会员 