数值计算作业Hermite插值
为了便于编写程序
Hermite 插值多项式H (x )的表达式如下所示:
H (x )=∑=n
i hi 1
[(x i
-x)(2a i
y i
-y i
)+y i
]
其中y i
=y (x i
),y '
i
=y '
(x i
)
h i =21)(∏≠=--n
i
j j j
i j
x x x x ,a i
=∑≠=-n
i
j j j i x x 11
用matlab 实现的算法:
1.输入X i ,Y i ,f 'i (x=1,2,3.......n),U
2.计算插值函数f 1). f=0 h=1 a=0 2). 对i=1:n
h
i
=21)(
∏≠=--n
i
j j j
i j x x x x
a i
=∑≠=-n
i
j j j i x x 11
3.f=f+h i
[(x i
-x)(2a i
y i
-y i
)+y i
] 4.输出f ,f0
编程如下:
function [f,f0]=Hermite(x,y,y_1,x0)
syms u;
f=0.0;
if(length(x)==length(y))
if (length(y)==length(y_1))
n=length(x);
else
disp; %y和y的导数的维数不相等
return
end
else
disp; % x和y的导数的维数不相等
return;
end
for i=1:n;
h=1.0;
p=0.0;
for j=1:n
if(j~=i)
h=h*(u-x(j))^2/(x(i)-x(j))^2;
p=p+1/(x(i)-x(j));
end
end
f=f+h*((x(i)-u)*(2*p*y(i)-y_1(i))+y(i));
end
f0=subs(f,'u',x0);
实例:
根据下面的数据点求出其埃尔米特插值多项式,并计算当x=1.44时y的值。
x 1 1.2 1.4 1.6 1.8 y 1 1.0954 1.1832 1.2649 1.3416 y'0.5000 0.4564 0.4226 0.3953 0.3727 解:在MATLAB命令窗口中输入以下命令:
>>x=1:0.2:1.8;
>>y=[1 1.0954 1.1832 1.2649 1.3416];
>>y_1=[0.5 0.4564 0.4226 0.3953 0.3727 ];
>>[f,f0]=Hermite(x,y,y_1,1.44)
程序运行结果如下:
f =
390625/576*(u-6/5)^2*(u-7/5)^2*(u-8/5)^2*(u-9/5)^2*(-61/3+6 4/3*u)+390625/36*(u-1)^2*(u-7/5)^2*(u-8/5)^2*(u-9/5)^2*(-11 4418258618928321/10995116277760000+337232823972231/3 5184372088832*u)+390625/16*(u-1)^2*(u-6/5)^2*(u-8/5)^2*( u-9/5)^2*(6660373488918492043/11258999068426240000+76 12884810106783/18014398509481984*u)+390625/36*(u-1)^2* (u-6/5)^2*(u-7/5)^2*(u-9/5)^2*(384779664999124623/2199023 2555520000-2855713758717179/281474976710656*u)+39062 5/576*(u-1)^2*(u-6/5)^2*(u-7/5)^2*(u-8/5)^2*(8968626627620 006931/175921860444160000-7762319875242775/2814749767 10656*u)
f0 =
1.2000
表格中的数据点是按y=x给出的,44
.1=1.2,插值结果是准确的。