三次样条插值作业题.docx
《三次样条插值作业题.docx》由会员分享，可在线阅读，更多相关《三次样条插值作业题.docx（12页珍藏版）》请在冰豆网上搜索。
三次样条插值作业题
例1设
为定义在[0,3]上的函数，有下列函数值表：
xi
0
1
2
3
yi
0
0.5
2
1.5
且
，
，试求区间[0,3]上满足上述条件的三次样条插值函数
本算法求解出的三次样条插值函数将写成三弯矩方程的形式：
其中，方程中的系数
，
，
，
将由Matlab代码中的变量Coefs_1、Coefs_2、Coefs_3以及Coefs_4的值求出。
以下为Matlab代码：
%=============================
%本段代码解决作业题的例1
%=============================
clearall
clc
%自变量x与因变量y，两个边界条件的取值
IndVar=[0,1,2,3];
DepVar=[0,0.5,2,1.5];
LeftBoun=0.2;
RightBoun=-1;
%区间长度向量，其各元素为自变量各段的长度
h=zeros（1,length（IndVar）-1）;
fori=1:
length（IndVar）-1
h（i）=IndVar（i+1）-IndVar（i）;
end
%为向量μ赋值
mu=zeros（1,length（h））;
fori=1:
length（mu）-1
mu（i）=h（i）/（h（i）+h（i+1））;
end
mu（i+1）=1;
%为向量λ赋值
lambda=zeros（1,length（h））;
lambda
（1）=1;
fori=2:
length（lambda）
lambda（i）=h（i）/（h（i-1）+h（i））;
end
%为向量d赋值
d=zeros（1,length（h）+1）;
d
（1）=6*（（DepVar
（2）-DepVar
（1））/（IndVar
（2）-IndVar
（1））-LeftBoun）/h
（1）;
fori=2:
length（h）
a=（DepVar（i）-DepVar（i-1））/（IndVar（i）-IndVar（i-1））;
b=（DepVar（i+1）-DepVar（i））/（IndVar（i+1）-IndVar（i））;
c=（b-a）/（IndVar（i+1）-IndVar（i-1））;
d（i）=6*c;
end
d（i+1）=6*（RightBoun-（DepVar（i+1）-DepVar（i））/（IndVar（i+1）-IndVar（i）））/h（i）;
%为矩阵A赋值
%将主对角线上的元素全部置为2
A=zeros（length（d）,length（d））;
fori=1:
length（d）
A（i,i）=2;
end
%将向量λ的各元素赋给主对角线右侧第一条对角线
fori=1:
length（d）-1
A（i,i+1）=lambda（i）;
end
%将向量d的各元素赋给主对角线左侧第一条对角线
fori=1:
length（d）-1
A（i+1,i）=mu（i）;
end
%求解向量M
M=A\d';
%求解每一段曲线的函数表达式
fori=1:
length（h）
Coefs_1=M（i）/（6*h（i））;
Part_1=conv（Coefs_1,...
conv（[-1,IndVar（i+1）],...
conv（[-1,IndVar（i+1）],[-1,IndVar（i+1）]）））;
S_1=polyval（Part_1,[IndVar（i）:
0.01:
IndVar（i+1）]）;
Coefs_2=M（i+1）/（6*h（i））;
Part_2=conv（Coefs_2,...
conv（[1,-IndVar（i）],...
conv（[1,-IndVar（i）],[1,-IndVar（i）]）））;
S_2=polyval（Part_2,[IndVar（i）:
0.01:
IndVar（i+1）]）;
Coefs_3=（DepVar（i）-M（i）*h（i）^2/6）/h（i）;
Part_3=conv（Coefs_3,[-1,IndVar（i+1）]）;
S_3=polyval（Part_3,[IndVar（i）:
0.01:
IndVar（i+1）]）;
Coefs_4=（DepVar（i+1）-M（i+1）*h（i）^2/6）/h（i）;
Part_4=conv（Coefs_4,[1,-IndVar（i）]）;
S_4=polyval（Part_4,[IndVar（i）:
0.01:
IndVar（i+1）]）;
S=S_1+S_2+S_3+S_4;
plot（[IndVar（i）:
0.01:
IndVar（i+1）],S,'LineWidth',1.25）
%在样条插值曲线的相应位置标注该段曲线的函数表达式
text（i-1,polyval（Part_1,3）,...
['\itS',num2str（i）,'（x）=',num2str（Coefs_1）,'（',num2str（IndVar（i+1））,'-x）^{3}+',...
num2str（Coefs_2）,'（x-',num2str（IndVar（i））,'）^{3}+',num2str（Coefs_3）,...
'（',num2str（IndVar（i+1））,'-x）+',num2str（Coefs_4）,'（x-',num2str（IndVar（i））,'）'],...
'FontName','TimesNewRoman','FontSize',14）
holdon
end
%过x=1和x=2两个横轴点作垂线%
line（[1,1],[2.5,-0.5],'LineStyle','--'）;
line（[2,2],[2.5,-0.5],'LineStyle','--'）;
%为x轴和y轴添加标注
xlabel（'\itx','FontName','TimesNewRoman',...
'FontSize',14,'FontWeight','bold'）;
ylabel（'\its（x）','FontName','TimesNewRoman',...
'Rotation',0,'FontSize',14,'FontWeight','bold'）;
最终，三次样条插值函数s（x）表达式为：
曲线的图像如图所示：
例2已知函数值表：
xi
1
2
4
5
yi
1
3
4
2
试求在区间[1,5]上满足上述函数表所给出的插值条件的三次自然样条插值函数
本算法求解出的三次样条插值函数将写成三弯矩方程的形式：
其中，方程中的系数
，
，
，
将由Matlab代码中的变量Coefs_1、Coefs_2、Coefs_3以及Coefs_4的值求出。
以下为Matlab代码：
%=============================
%本段代码解决作业题的例2
%=============================
clearall
clc
%自变量x与因变量y的取值
IndVar=[1,2,4,5];
DepVar=[1,3,4,2];
%区间长度向量，其各元素为自变量各段的长度
h=zeros（1,length（IndVar）-1）;
fori=1:
length（IndVar）-1
h（i）=IndVar（i+1）-IndVar（i）;
end
%为向量μ赋值
mu=zeros（1,length（h））;
fori=1:
length（mu）-1
mu（i）=h（i）/（h（i）+h（i+1））;
end
mu（i+1）=0;
%为向量λ赋值
lambda=zeros（1,length（h））;
lambda
（1）=0;
fori=2:
length（lambda）
lambda（i）=h（i）/（h（i-1）+h（i））;
end
%为向量d赋值
d=zeros（1,length（h）+1）;
d
（1）=0;
fori=2:
length（h）
a=（DepVar（i）-DepVar（i-1））/（IndVar（i）-IndVar（i-1））;
b=（DepVar（i+1）-DepVar（i））/（IndVar（i+1）-IndVar（i））;
c=（b-a）/（IndVar（i+1）-IndVar（i-1））;
d（i）=6*c;
end
d（i+1）=0;
%为矩阵A赋值
%将主对角线上的元素全部置为2
A=zeros（length（d）,length（d））;
fori=1:
length（d）
A（i,i）=2;
end
%将向量λ的各元素赋给主对角线右侧第一条对角线
fori=1:
length（d）-1
A（i,i+1）=lambda（i）;
end
%将向量d的各元素赋给主对角线左侧第一条对角线
fori=1:
length（d）-1
A（i+1,i）=mu（i）;
end
%求解向量M
M=A\d';
%求解每一段曲线的函数表达式
fori=1:
length（h）
Coefs_1=M（i）/（6*h（i））;
Part_1=conv（Coefs_1,...
conv（[-1,IndVar（i+1）],...
conv（[-1,IndVar（i+1）],[-1,IndVar（i+1）]）））;
S_1=polyval（Part_1,[IndVar（i）:
0.01:
IndVar（i+1）]）;
Coefs_2=M（i+1）/（6*h（i））;
Part_2=conv（Coefs_2,...
conv（[1,-IndVar（i）],...
conv（[1,-IndVar（i）],[1,-IndVar（i）]）））;
S_2=polyval（Part_2,[IndVar（i）:
0.01:
IndVar（i+1）]）;
Coefs_3=（DepVar（i）-M（i）*h（i）^2/6）/h（i）;
Part_3=conv（Coefs_3,[-1,IndVar（i+1）]）;
S_3=polyval（Part_3,[IndVar（i）:
0.01:
IndVar（i+1）]）;
Coefs_4=（DepVar（i+1）-M（i+1）*h（i）^2/6）/h（i）;
Part_4=conv（Coefs_4,[1,-IndVar（i）]）;
S_4=polyval（Part_4,[IndVar（i）:
0.01:
IndVar（i+1）]）;
S=S_1+S_2+S_3+S_4;
plot（[IndVar（i）:
0.01:
IndVar（i+1）],S,'LineWidth',1.25）
%在样条插值曲线的相应位置标注该段曲线的函数表达式
text（i,polyval（Part_1,5）,...
['\itS',num2str（i）,'（x）=',num2str（Coefs_1）,'（',num2str（IndVar（i+1））,'-x）^{3}+',...
num2str（Coefs_2）,'（x-',num2str（IndVar（i））,'）^{3}+',num2str（Coefs_3）,...
'（',num2str（IndVar（i+1））,'-x）+',num2str（Coefs_4）,'（x-',num2str（IndVar（i））,'）'],...
'FontName','TimesNewRoman','FontSize',14）
holdon
end
%过x=2和x=4两个横轴点作垂线%
line（[2,2],[4.5,0.5],'LineStyle','--'）;
line（[4,4],[4.5,0.5],'LineStyle','--'）;
%为x轴和y轴添加标注
xlabel（'\itx','FontName','TimesNewRoman',...
'FontSize',14,'FontWeight','bold'）;
ylabel（'\its（x）','FontName','TimesNewRoman',...
'Rotation',0,'FontSize',14,'FontWeight','bold'）;
最终，三次自然样条插值函数s（x）表达式为：
曲线的图像如图所示：
例3课后习题与思考题第7题
xi
0.25
0.30
0.39
0.45
0.53
yi
0.5000
0.5477
0.6245
0.6708
0.7280
试求在区间[0.25,0.53]上满足上述函数表所给出的插值条件的三次自然样条插值函数s（x）
求解出的三次样条插值函数将写成三弯矩方程的形式：
本题采用和例2基本相同的Matlab代码，只改变初始条件。
最终，三次自然样条插值函数s（x）表达式为：
曲线的图像如图所示：
例4课后习题与思考题第6题
求
的二次插值式
，使：
并计算
的近似值并估计误差。
本题采用拉格朗日二次插值法进行计算：
以下为Matlab代码：
%=============================
%本段代码解决课本第2章习题与思考题第6题
%=============================
clearall
clc
%自变量x与因变量y的取值
IndVar=[100,121,144];
DepVar=[10,11,12];
%构造拉格朗日插值函数
Coefs_1=DepVar
（1）/...
（（IndVar
（1）-IndVar
（2））*（IndVar
（1）-IndVar（3）））;
Part_1=conv（Coefs_1,...
conv（[1,-IndVar
（2）],[1,-IndVar（3）]））;
f_1=polyval（Part_1,[IndVar
（1）:
0.01:
IndVar（3）]）;
Coefs_2=DepVar
（2）/...
（（IndVar
（2）-IndVar
（1））*（IndVar
（2）-IndVar（3）））;
Part_2=conv（Coefs_2,...
conv（[1,-IndVar
（1）],[1,-IndVar（3）]））;
f_2=polyval（Part_2,[IndVar
（1）:
0.01:
IndVar（3）]）;
Coefs_3=DepVar（3）/...
（（IndVar（3）-IndVar
（1））*（IndVar（3）-IndVar
（2）））;
Part_3=conv（Coefs_3,...
conv（[1,-IndVar
（1）],[1,-IndVar
（2）]））;
f_3=polyval（Part_3,[IndVar
（1）:
0.01:
IndVar（3）]）;
f=f_1+f_2+f_3;
plot（[IndVar
（1）:
0.01:
IndVar（3）],f,'LineWidth',1.25）;
%在样条插值曲线的相应位置标注该段曲线的函数表达式
text（110,polyval（Part_1,110）+polyval（Part_2,110）+polyval（Part_3,110）,...
['\itP_2（x）=',num2str（Coefs_1）,'（x-',num2str（IndVar
（2））,...
'）（x-',num2str（IndVar（3））,'）+',num2str（Coefs_2）,'（x-',num2str（IndVar
（1））,...
'）（x-',num2str（IndVar（3））,'）+',num2str（Coefs_3）,'（x-',num2str（IndVar
（1））,...
'）（x-',num2str（IndVar
（2））,'）'],...
'FontName','TimesNewRoman','FontSize',14）
%为x轴和y轴添加标注
xlabel（'\itx','FontName','TimesNewRoman',...
'FontSize',14,'FontWeight','bold'）;
ylabel（'\itP_2（x）','FontName','TimesNewRoman',...
'Rotation',0,'FontSize',14,'FontWeight','bold'）;
%求115平方根的估计值
fun115=polyval（Part_1,115）+polyval（Part_2,115）+polyval（Part_3,115）
%求115平方根的准确值
sqrt（115）
%求误差，单位为%
est=100*（fun115-sqrt（115））/fun115
最终，拉格朗日插值函数
的表达式为：
从计算结果可以看出，
的准确值为10.7238，而通过拉格朗日插值法求出的近似值为10.7228，误差为0.0098%
的图像为：
升级会员 