计算方法试验7

1、实验七一、阅读理解下列程序，并在计算机上运行.1. newt .m%x,y分别为插值节点列向量及节点相应的函数值 列向量，A为差商表矩阵， z为插值点，s为插值点的newton插值多项式计算值function A,s=newt(x,y,z) n=length(x);A=zeros(n);A(:,1)=y;c(1)=y； forj=2:nfor i=j:nA(i,j)=(A(i,j-1)-A(i-1,j-1)”(x(i)-x(i-j+1);endc(j)=A(j,j);%差商表每一列的第一个数即 newton插值多项式的系数ends=c(n);for k=n-1:-1:1%秦九韶算法计算newt

2、on插值多项式的值s=s*(z-x(k)+c(k);endA=x,A;以课本P32例2、P33例3测试，并计算课本P6016(X=225)i2. gaocilam (以 z=4.5 测试(此时 y=1.578721)%验证高次插值的Runge现象function y,err=gaocila(z) % z 为插值点 x=-5:5;n=11;y=f(x);for i=1:nk=1;forj=1:nif j=i k=k*(x-X(j)/(X(i)-X(j);endendy=y+k*Y(i); end err=f(z)-y;function y=f(x) % 子函数 y=1 ./(1+x .*x);3

3、. neville .m%x,y分别为插值节点向量及节点相应的函数值 列向量，A为逐步插值表,z为插值点function A=neville(x,y,z) n=length(x);A=zeros(n);A(:,1)=y; forj=2:nfor i=j:nA(i,j)=A(i-1,j-1)+(A(i,j-1)-A(i-1,j-1)*(z-x(i-j+1)/(x(i)-x(i-j+1);endendA=x,A;以课本P38例6测试.4. nevillel m%x,y分别为插值节点向量及节点相应的函数值列向量,A为逐步插值表,z为插值点，w为精度function A=neville1(x,y,z,

4、w) n=length(x);A=zeros(n);A(:,1)=y;forj=2:nfor i=j:nA(i,j)=A(i-1,j-1)+(A(i,j-1)-A(i-1,j-1)*(z-x(i-j+1)/(x(i)-x(i-j+1);if abs(A(i,j)-A(i,j-1)vw |abs(A(i,j)-A(i-1,j)vw w=1;breakendendif w=1breakendendA=x,A;以课本P38例6,对不同的精度进行测试二、编程并在计算机上调试修改运行1. 编写复化梯形公式及复化Simpson公式求值程序计算课本P75例2、P8817 2. 应用变步长梯形法编程计算 课本

5、P77例3.3. 选做题：应用Romberg算法编程计算课本P81例4.1. n ewt.m>> A,s=newt(100 121 144',10 11 12',115)100.000010.000000121.000011.00000.04760144.000012.00000.0435-0.0001s =10.7228>> x=0.3 0.4 0.5 0.6 0.7'>> y=0.29850 0.39646 0.49311 0.58813 0.68122'>> A,s=newt(x,y,0.462)A =0.3

6、0000.298500000.40000.39650.97960000.50000.49310.9665-0.0655000.60000.58810.9502-0.0815-0.053300.70000.68120.9309-0.0965-0.05000.0083s =0.4566>> A,s=newt(-2 -1 0 1 2 3',-5 1 1 1 7 25',2.25)A =-2-500000-1160000010-30001100100276310032518610010.14062. gaocila.m (以 z=4.5 测试(此时 y=1.578721 )

7、正确的程序：%佥证高次插值的Runge现象function y,err=gaocila(z)%z为插值点x=-5:5; n=11;Y=f(x);y=0;for i=1:nk=1;for j=1: nif j=ik=k*(z-x(j)/(x(i)-x(j);endendy=y+k*Y(i);enderr=f(z)-y;fun ctio n y=f(x)% 子函数y=1./(1+x.*x);>> y,err=gaocila(4.5) y =1.5787err =-1.53173. n eville.m>> x=0.3 0.4 0.5 0.6 0.7'>>

8、 y=0.2985 0.39646 0.49311 0.58813 0.68122'5>> neville(x,y,0.462) ans =Colu mns 1 through 30 0 060.3000000000000000.4000000000000000.5000000000000000.6000000000000000.7000000000000000.2985000000000000.3964600000000000.4931100000000000.5881300000000000.68122000000000000.4571952000000000.4563

9、830000000000.4570024000000000.459665800000000Colu mns 4 through 60000000.456537318000000000.4565750140000000.45655767384000000.4564963540000000.4565587576000000.4565581127628004. neville1.m>> x=0.3 0.4 0.5 0.6 0.7'>> y=0.2985 0.39646 0.49311 0.58813 0.68122'>> neville1(x,y,

10、0.462,1e-3)ans =Colu mns 1 through 30.3000000000000000.4000000000000000.5000000000000000.6000000000000000.7000000000000000.2985000000000000.3964600000000000.4931100000000000.5881300000000000.68122000000000000.4571952000000000.45638300000000000Colu mns 4 through 60 0 00 0 00 0 00 0 0>> neville1

11、(x,y,0.462,1e-4)150.2985000000000000.3964600000000000.4931100000000000.5881300000000000.6812200000000000 00 0000 0ans =Colu mns 1 through 30.3000000000000000.4000000000000000.5000000000000000.6000000000000000.700000000000000Colu mns 4 through 6000.4565373180000000.456575014000000000.4571952000000000

12、.4563830000000000.4570024000000000.459665800000000>> neville1(x,y,0.462,1e-5)0.2985000000000000.3964600000000000.4931100000000000.5881300000000000.681220000000000ans =Colu mns 1 through 30.3000000000000000.4000000000000000.5000000000000000.6000000000000000.700000000000000Colu mns 4 through 600

13、00.4571952000000000.4563830000000000.4570024000000000.4596658000000000.45653731800000000.456575014000000 0.4565576738400000.456496354000000 0.4565587576000001、复化梯形公式课本P75例2fun ctio n t1=T(a,b, n) h=(b-a)/n;T1=0;for i=1:n-1x=a+i*h;T1=2*f(x)+T1;endt1=(f(a)+T1+f(b)*(b-a)/(2* n);fun ctio n y=f(x)if x=0y

14、=s in( x)/x;elsey=1;end>> T(0,1,8)ans =0.9457P88第17题fun ctio nt1=t(a,b, n)h=(b-a)/n;T=0;for i=1:n-1x=a+i*h;T=2*f(x)+T;end t1=(1+T+f(b)*(b-a)/(2* n); fun ctio ny=f(x)y=1+exp(-x)*si n( 4*x);ans =1.2836复化Simpson公式课本P75例2fun ctio n S=s(a,b ,n) h=(b-a)/(2* n);s1=0;s2=0;for i=1:n-1x=a+2*i*h;s1=2*f(x

15、)+s1;%求节点和endfor j=0:n-1x=a+(2*j+1)*h;s2=4*f(x)+s2;%求中间节点的和endS=(f(a)+s1+s2+f(b)*(b-a)/(6* n); fun ctio ny=f(x)if x=0 y=si n(x)/x;else y=1;end>> s(0,1,4)ans =0.9461P88第17题 fun ctio n S=s(a,b ,n) h=(b-a)/(2* n);s1=0;s2=0;for i=1:n-1x=a+2*i*h;s1=2*f(x)+s1;%求节点和endfor j=0:n-1x=a+(2*j+1)*h;s2=4*f(

16、x)+s2;%求中间节点的和endS=(f(a)+s1+s2+f(b)*(b-a)/(6* n);fun ctio n y=f(x)y=1+exp(-x)*si n( 4*x);>> s(0,1,2)ans =1.30942、%变步长梯形法fun ctio ny=T(a,b,e)h=b-a;T1=h*(f(a)+f(b)/2;A=T1;while 1s=0;x=a+h/2;while 1if x<bs=s+f(x);x=x+h;elsebreak ;endendT2=T1/2+h*s/2;对应的积分值 Tn。A=A,T2;%己录区间等分数 n=2Ak,k=0,1,2if ab

17、s(T2-T1)>=eh=h/2;T1=T2; else y=T2;break ; endendAfunction y=f(x) if x=0y=si n( x)/x; elsey=1;end>> T(0,1,1e-10)A =Colu mns 1 through 3 0.920735492403948Colu mns 4 through 60.945690863582701Columns 7 through 90.946076943060063Colu mns 10 through 120.9397932848061770.9459850299343860.94608153

18、85431520.9445135216653900.9460585609627680.9460826874113470.9460830643834990.9460830702736900.946082974628235 0.946083046432447Colu mns 13 through 150.946083068871262 0.946083069993204Colu mn 160.946083070343813 ans =0.9460830703438133.选做题：fun ctio ny=T(a,b,e)h=b-a;T1=h*(f(a)+f(b)/2;A=T1;while 1s=0;x=a+h/2;while 1if x<bs=s+f(x);x=x+h;elsebreak ;endendT2=T1/2+h*s/2;A=A,T2;%己录区间等分数n=2Ak,k=0,1,2对应的积分值Tn。