三次样条插值法《数值分析》上机实验作业

《三次样条插值法《数值分析》上机实验作业》由会员上传分享，免费在线阅读，更多相关内容在教育资源-天天文库。

1、三次样条插值法《数值分析》上机实验作业昆明理工大学研究生《数值分析》上机实验作业姓名：学号：专业：一、课题名称课题七三次样条插值法1、问题提出设已知数据如下：求f(x)的三次样条插值函数S(x)。2、要求（1）满足自然边界条件s??(0.2)?s??(1.0)?0；（2）满足第一类边界条件s?(0.2)?0.20271，s?(1.0)?1.55741。（3）打印输出用追赶法解出的弯曲向量(M0,M1,M2,M3,M4)和S(0.2?0.1i)（i=0,1,2,3,4,5,6,7,8）的值。并画出y?S(x)的

2、图形。二、班级、姓名、学号三、目的和意义由于航空、造船等工程设计的需要而发展起来所谓样条插值方法，既保留了分段低次插值多项式的各种优点，又提高了插值函数的光滑性，而且具有较好的稳定性。今天，样条插值方法已成为数值逼近的一个极其重要的分支，在许多领域里得到越来越广泛地应用。其中，尤以三次样条插值函数应用最为广泛，如在高速飞机的机翼形体和船体放样等方面的应用，同时在计算机作图方面更是大有作为。它能够解决一些既有二阶光滑度，又有二阶连续导数的方程，具有良好的收敛性和稳定性。1.通过本次实验进一步了解三次样条插值函数

3、，并通过求解三弯矩方程组得出曲线函数组；2.通过MATLAB编程实现求三次样条插值函数的算法，分别考虑不同的边界条件，同时用追赶法解出弯曲向量和S(0.2?0.1i)（i=0,1,2,3,4,5,6,7,8）的值。四、计算公式(xj)?Mj(j?0,1,2,?,n)表首先我们利用S(x)的二阶导数值S??达S(x),因为在区间[xj,xj?1]上S(x)?Sj(x)是不高于三次的多项式，其二阶导数S??(x)必是线性函数，所以可表示为：S??(x)?Mjxj?1?xhj?Mj?1x?xjhj,x?[xj,x

4、j?1]对S??(x)积分两次并利用S(xj)?yj,S(xj?1)?yj?1，可定出积分常数，于是得三次样条表达式。S(x)?Mj?(yj?1?(xj?1?x)36hj6h)2j?Mj?1(x?xj)36hj?(yj?Mj6hj)2xj?1?xhjMj?1x?xjhj，j?0,1,2,?,n?1这里Mj(j?0,为了确定Mj(j?0,对S(x)1,?,n)是未知的。1,?,n)，求导得S?(x)??Mj(xj?1?x)22hj?Mj?1(x?xj)22hj?yj?1?yjhj?hj6(Mj?1?Mj)由此可

5、得S?(xj?0)?yj?1?yjhj?hj6Mj?1?hj3Mj类似的可求出S(x)在区间[xj-1,xj]上的表达式，进而得S?(xj?0)?yj?yj?1hj?1?hj?16Mj?1?hj?13Mj利用S?(xj?0)?S?(xj?0)可得μjMj?1?2Mj?λjMj?1?dj，j?1,?,n?1,（三弯矩方程）其中μj?dj?6hj?1hj?1?hj，λj?hjhj?1?hj，,xj,xj?1.j?1,?,n?1,fxj,xj?1?fxj?1,xjhj?1?hj?????6f?xj?1?其中有（n?

6、1）个未知数M0,M1,...Mn，而方程只有（n-1）个，当满足第一种边界条件时，可得另两个方程2M0?M1?6?f?x0,x1??f0??，h06?fn??f?xn?1,xn??hn?16?f?x0,x1??f0??,?n?1,gn?6?fn??f?xn?1,xn??，将h0hn?1Mn?1?2Mn?如果令?0?1,d0?上述方程综合后的一下矩阵形式：?2λ0?λ1?μ12?????μn?12λn?1??μn2???M0????M1???????Mn?1????Mn?d0?????d?1?????????

7、??dn?1???d????n?可以证明此方程组满足追赶法的条件，我们用追赶法可得M的(n?1)值，将其带入公式即得s(x)。对第二种边界条件，直接的端点方程M0?f0??,Mn?fn??并且令?0??n?0,d0?2f0??,dn?2fn??，则又得三弯矩方程同理即可求得解。五、结构程序设计1.满足自然边界条件时自定义函数：followup.m%追赶法求m%A为线性方程组的系数矩阵%b为常数向量functionm=followup(A,b)n=rank(A);fori=1:nifA(i,i)==0disp(

8、'error:对角元素中有数据为0');return;endendd=ones(n,1);a=ones(n-1,1);c=ones(n-1);fori=1:n-1a(i,1)=A(i+1,i);c(i,1)=A(i,i+1);d(i,1)=A(i,i);endd(n,1)=A(n,n);fori=2:nd(i,1)=d(i,1)-(a(i-1,1)/d(i-1,1))*c(i-1,1