上海交通大学计算方法作业答案

P50 – 1

%%牛顿插值多项式

function [ c, d] = newpoly( x,y )

%这里x为n个节点的横坐标所组成的向量，y为纵坐标所组成的向量。

%c为所求的牛顿插值多项式的系数构成的向量。

n=length(x);

d=zeros(n, n);

d(: , 1)=y';

for j=2 : n

for k=j : n

d(k, j)=(d(k, j-1) - d(k-1, j-1)) / (x(k)-x(k-j+1));

end

end

c =d(n, n);

for k=(n-1) : - 1 : 1

c =conv(c, poly(x(k)));

m=length(c);

c(m)=c(m)+d(k, k);

end

>> x = 0.2:0.2:1 ;

>> y =[ 0.98,0.92,0.81,0.64,0.38] ;

>> c= newpoly(x, y )

c = -0.5208 0.8333 -1.1042 0.1917 0.9800

%%三次样条插值

x=[0.2,0.4,0.6,0.8,1.0];

y=[0.98,0.92,0.81,0.64,0.38];

x0 = [0.2,0.28,1.0,1.08];

pp=csape(x,y,'variational');

%% 三次样条函数表达式

disp(pp.coefs);

-1.3393 -0.0000 -0.2464 0.9800

0.4464 -0.8036 -0.4071 0.9200

-1.6964 -0.5357 -0.6750 0.8100

2.5893 -1.5536 -1.0929 0.6400

绘制曲线图

x2 = 0:0.01:1.2 ;

y11 = polyval(c,x2) ;

y22 = ppval(pp,x2);

x0 = [0.2,0.28,1.0,1.08];

y110 = polyval(c,x0);

y220 = ppval(pp,x0);

plot(x2,y11,'r',x0,y110,'^',x2,y22,'g',x0,y220,'h')

legend('牛顿插值','牛顿插值样点','三次样条插值','三次样条插值样点')

P50 -3

(1)x = [0,1,4,9,16,25,36,49,64] ;

y = 0:8 ;

x1 = 0:0.1:64 ;

x2 = 0:0.01:1 ;

f = lagrange(x,y)

%% 得到多项式函数表达式

L(x)= - 3.28063e-10*x^8 + 6.71268e-8*x^7 - 0.00000542921*x^6 + 0.000222972*x^5 - 0.00498071*x^4 + 0.0604294*x^3 - 0.38141*x^2 + 1.32574*x

y1 = lagrange(x,y,x1) ;

y2 = lagrange(x,y,x2) ;

(2)(2) x = [0,1,4,9,16,25,36,49,64] ;

y = 0:8 ;

x1 = 0:0.1:64 ;

x2 = 0:0.01:1 ;

%% 得到三次样条差值函数表达式

pp=csape(x,y,'not-a-knot');

disp(pp.coefs);

0.0266 -0.2998 1.2732 0

0.0266 -0.2199 0.7534 1.0000

-0.0021 0.0197 0.1529 2.0000

0.0005 -0.0112 0.1955 3.0000

-0.0000 -0.0001 0.1160 4.0000

0.0000 -0.0014 0.1026 5.0000

0.0000 -0.0005 0.0825 6.0000

0.0000 -0.0004 0.0717 7.0000

y11 = ppval(pp,x1) ;

y22 = ppval(pp,x2) ;

绘制图形

（1）在[0,64]

显然随着次数越高，多项式插值出现误差很大（2）[0,1]

在[0,1]区间上三次样条插值和多项式插值基本一致

P137-1

insucomplex_4_1.m 文件

clear ;

clc ;

%h为步长,可分别令h=1,0.1,0.01,0.001

h = [1,0.1,0.01,0.001]

for i = 1:4

h(i) ,

x=0:h(i):1;

y=sqrt(x).*log(x+eps);

%复化梯形公式

T=trapz(x,y);

T=vpa(T,7),

f=inline('sqrt(x).*log(x)',x);

%复化辛普生公式

S=quadl(f,0,1);

S=vpa(S,7),

end

>> t = -log(h) ;

>> plot(t,T,'rs',t,S,'r*')

>> lengend('复合梯形公式','复合梯形公式')