8.2.1 差商求导及其MATLAB程序
例 8.2.1 设.
（2）将（1）中计算的的近似值分别与精确值比较.
解 （1）编写计算的一阶导数计算的近似值和误差估计的MATLAB程序，并输入
>> x=0.79;h=[0.1,0.01,0.001,0.0001];
M=80;x1=x+h;x2=x-h; y=sin(5.*x.^2-21);
y1=sin(5.*x1.^2-21); y2=sin(5.*x2.^2-21); yq=(y1-y)./h, yh=(y-y2)./h,
wu=abs(h.*M/2), syms x,f=sin(5.*x.^2-21); yx=diff(f,x)
运行后屏幕显示利用前差公式和后差公式计算的近似值yq，yh和误差估计wu，取4位小数点计算，其中步长分别取，M=80，导函数yx
yq =
1.46596380397978 4.22848550173043 4.44250759584697 4.46320955293622
yh =
5.96885352366536 4.68672022108227 4.48833808130555 4.46779260847907
wu =
4.00000000000000 0.40000000000000 0.04000000000000 0.00400000000000
yx =
10*cos(5*x^2-21)*x
（2）计算的值.输入程序
>> x=0.79; yx =10*cos(5*x^2-21)*x,
wuq=abs(yq-yx), wuh=abs(yh-yx)
运行后屏幕显示利用前差公式和后差公式计算的近似值与精确值的绝对误差wuq，wuh和的精确值yx如下
yx =
4.46550187104484
wuq =
2.99953806706506 0.23701636931441 0.02299427519787 0.00229231810861
wuh =
1.50335165262053 0.22121835003744 0.02283621026072 0.00229073743424
8.2.2 中心差商公式求导及其MATLAB程序
利用精度为的三点公式计算的近似值和误差估计的MATLAB主程序
function [n,xi,yx,wuc]=sandian(h,xi,fi,M)
n=length(fi); yx=zeros(1,n); wuc=zeros(1,n); x1= xi(1); x2= xi(2); x3= xi(3);
y1=fi(1); y2=fi(2); y3=fi(3); xn= xi(n); xn1= xi(n-1);
xn2= xi(n-2); yn=fi(n); yn1=fi(n-1); yn2=fi(n-2);
for k=2:n-1
yx(1)=(-3*y1+4*y2-y3)/(2*h); yx(n)=(yn2-4*yn1+3*yn)/(2*h);
yx(2)=( fi(3)- fi(1))/(2*h); yx(k)=( fi(k+1)- fi(k-1))./(2*h);
wuc(1)=abs(h.^2.*M./3); wuc(n)=abs(h.^2.*M./3);
wuc(2:n-1)=abs(-h.^2.*M./6);
end
利用精度为的三点公式计算的近似值和误差估计的MATLAB主程序
function [x,yxj, wuc]=sandian3(h,xi,fi,M)
matlab求导yxj=zeros(1,3); wuc=zeros(1,3); x1= xi(1);
x2= xi(2); x3= xi(3); y1=fi(1); y2=fi(2); y3=fi(3);
for t=1:3
s(t)=((2*t-5)*y1-4*(t-2)*y2+(2*t-3)*y3)/(2*h); x=xi; y=s(t); yxj(t)=y;
if t==2
wuc(2)=abs(-h.^2*M/6);
else
wuc(1:2:3)=abs(h.^2*M/3);
end
end
例 8.2.3 设已给出的数据表8–5：
表8–5
x | 1.000 0 1.100 0 1.200 0 1.300 0 1.400 0 1.500 0 1.600 0 |
f(x) | 0.250 0 0.226 8 0.206 6 0.189 0 0.173 6 0.160 0 0.147 9 |
M= 0.750 2,试用三点公式计算下列问题：
(1)当h=0.1时，在x=1.000 0，1.100 0，1.200 0，1.300 0，1.400 0，1.500 0，1.600 0处的一阶导数的近似值，并估计误差；
(2)当h=0.2时，在x=1.000 0，1.200 0， 1.400 0，1.600 0处的一阶导数的近似值，并估计误差；
(3)当h=0.3时，在x=1.000 0，1.300 0 ，1.600 0处的一阶导数的近似值，并估计误差；
(4) 表8–5中的数据是函数在相应点的数值，将(1)至(3)计算的一阶导数的近似值与的一阶导数值比较，并求出它们的绝对误差.
解 （1）保存M文件sandian.m，sandian3.m；
（2）在MATLAB工作窗口输入如下程序
>> syms x,y=1/((1+x)^2); yx=diff(y,x,1),
yx3=diff(y,x,3),
运行后将屏幕显示的结果为
yx = yx3 =
-2/(1+x)^3 -24/(1+x)^5
（3）在MATLAB工作窗口输入如下程序
>>h=0.1; xi=1.0000:h:1.6000;
fi=[0.2500 0.2268 0.2066 0.1890 0.1736 0.1600 0.1479];
x=1:0.001:1.6; yx3 =-24./(1+x).^5; M= max(abs(yx3));
[n1,x1,yx1,wuc1]=sandian(h,xi,fi,M)
yxj1=-2./(1+xi).^3,wuyxj1=abs(yxj1- yx1)
h=0.2; xi=1.0000:h:1.6000; fi=[0.2500 0.2066 0.1736 0.1479];
x=1:0.001:1.6; yx3 =-24./(1+x).^5; M= max(abs(yx3));
[n2,x2,yx2,wuc2]=sandian(h,xi,fi,M)
yxj2=-2./(1+xi).^3,wuyxj2=abs((yxj2- yx2))
h=0.3; xi=1.0000:h:1.6000; fi=[0.2500 0.1890 0.1479];
x=1:0.001:1.6; yx3 =-24./(1+x).^5; M=max(abs(yx3));
[x3,yx3, wuc3]=sandian3(h,xi,fi,M)
yxj3=-2./(1+xi).^3,wuyxj3=abs(yxj3- yx3)
或 >> h1=0.1,
x=[1.0000,1.1000,1.2000,1.3000,1.4000,1.5000,1.6000];
f=[0.2500,0.2268,0.2066,0.1890,0.1736,0.1600,0.1479];
xi=x(1:3);f11=f(1:3); M= 0.7502;
[x11,yxj11,wuc11]=sandian3(h1,xi,f11,M),
xi= x(4:6);f12=f(4:6);
[x12,yxj12,wuc12]=sandian3(h1,xi,f12,M),
xi=x(5:7);f13=f(5:7); [x13,yxj13,wuc13]=sandian3(h1,xi,f13,M),
h2=0.2, xi= x(1:2:5);f21= f(1:2:5);
[x21,yxj21,wuc21]=sandian3(h2,xi,f21,M),
xi= x(2:2:6);f22=f(2:2:6);
[x22,yxj22,wuc22]=sandian3(h2,xi,f22,M),
xi= x(3:2:7);f23=f(3:2:7);
[x23,yxj23,wuc23]=sandian3(h2,xi,f23,M),
h3=0.3, xi= x(1:3:7);f31= f(1:3:7);
[x31,yxj31,wuc31]=sandian3(h3,xi,f31,M),
将运行的结果（略）.
8.2.3 理查森外推法求导及其MATLAB程序
（一）一般形式的理查森外推法及其MATLAB程序
利用理查森外推法计算的近似值和误差估计的MATLAB程序
function [Dy,dy,jdw,n]=diffext1(fun,x0,jdwc,max1)
h=1;j=1; n=1;jdW=1;xdW=1; x1=x0+h;x2=x0-h;
Dy(1,1)=(feval(fun,x1)- feval(fun,x2))/(2*h);
while((jdW>jdwc)&(j<max1))
j;x1=x0+2^(-j)*h;x2=x0-2^(-j)*h;
Dy(j+1,1)=(feval(fun,x1)-feval(fun,x2))/(2^(1-j)*h);
for k=1:j
k;Dy(j+1,k+1)= Dy(j+1,k)+( Dy(j+1,k)- Dy(j,k))/(4^k-1);
end
jdW=abs(Dy(j+1,j+1)-Dy(j+1,j)); j=j+1;
end
[n,n]=size(Dy);jdw=abs(Dy(n,n)-Dy(n,n-1));
发布评论