s=0; for i=1:n
t=ones(1,length(xx)); for j=1:n if j~=i
t=t.*(xx-x(j))/(x(i)-x(j)); end end s=s+t*y(i); end yy=s;
在Matlab窗口中执行： x=[-2.00 0.00 1.00 2.00]; y=[17.00 1.00 2.00 17.00]; xx=0.6; zuoye10(x,y,xx) 结果如下： ans =
0.25600000000000
13．编制Newton插值法的通用Matlab程序，并求f(0.5)的近似值. 已知的数值如下表所示
xi00.20.40.60.8
f(xi)0.19950.39650.58810.7721.09461
16
解：
Newton插值法的通用Matlab程序：
x=[0 0.2 0.4 0.6 0.8]';y=[0.1995 0.3965 0.5881 0.7721 0.9461]';xx=0.5; m=length(x);n=length(y);
if m~=n, error('向量x与y的长度必须一致');end for j=2:n，for i=j:n
newton_chazhi(i,j)=(newton_chazhi(i-1,j-1)-newton_chazhi(i,j-1))/(x(1+i-j)-x(i)); end end
yy= newton_chazhi(1,1);
t=1;for i=2:n，t=t*(xx-x(i-1));yy=yy+ newton_chazhi(i,i)*t; end
%计算误差近似值
err=newton_chazhi(n,n);for i=1:n,err=err*(xx-x(i));end,err yy
运行结果为：
err =-2.343750000006213e-06 yy=0.681187500000000
??y??f(x,y)??y|x?x0?y020．编写解常微分方程?的四阶龙格库塔法通用程序.
解：
四阶龙格库塔法通用程序：
function [x,y]=zuoye20(dyfun,xspan,y0,h)
17
x=xspan(1):h:xspan(2); y(1)=y0;
for n=1:length(x)-1
k1=feval(dyfun,x(n),y(n));
k2=feval(dyfun,x(n)+h/2,y(n)+h/2*k1); k3=feval(dyfun,x(n)+h/2,y(n)+h/2*k2); k4=feval(dyfun,x(n+1),y(n)+h*k3); y(n+1)=y(n)+h*(k1+2*k2+2*k3+k4)/6; end x=x'; y=y';
21. 利用上题的程序求解初值问题:
?dy2?2??xy,?dx3??y(0)?1x?[0,1.2]
32y?1?x的数值解,并将结果与准确解进行比较. 取h?0.4.
解：
Matlab程序：
function [x,y]=zuoye21(dyfun,xspan,y0,h) x=xspan(1):h:xspan(2); y(1)=y0;
for n=1:length(x)-1
18
k1=feval(dyfun,x(n),y(n));
k2=feval(dyfun,x(n)+h/2,y(n)+h/2*k1); k3=feval(dyfun,x(n)+h/2,y(n)+h/2*k2); k4=feval(dyfun,x(n+1),y(n)+h*k3); y(n+1)=y(n)+h*(k1+2*k2+2*k3+k4)/6; end x=x' y=y'
在Matlab窗口中执行： dyfun=inline('2*x*y^(-2)/3'); xspan=[0,1.2]; y0=1; h=0.4;
nark18(dyfun,xspan,y0,h) 结果如下： x =
0 0.40000000000000 0.80000000000000 1.20000000000000 y =1.00000000000000 1.05075062243354 1.17933175912031 1.34631542500313
19
计算值与准确值比较如下：
22. 编写四阶龙格库塔法程序，并求解洛伦兹型系统（不许调用现成的龙格库塔法程序软件包）：
????x??y??yz?x???rx?qy?sxz?y?z????bz??xy
解：为了便于求解，这里首先为各个参数赋予初值，实际中根据需要修改程序中的参数即可。
取??0.25，??0.06，b?0.4，??0.5，r?120，q?1.3，s?1.5，???20.取初始点为(0.005,0.4596,?0.1146). 程序如下：
%预定义homework22.m文件 function f=homework22(t,y)
20