三次样条函数的自动求法(学院+专业+学号)

第12章三次样条众所周知，使用高阶多项式的插值常常产生病态的结果。

目前，有多种消除病态的方法。

在这些方法中，三次样条是最常用的一种。

在MATLAB中，实现基本的三次样条插值的函数有spline，ppval，mkpp和unmkpp。

在这些函数中，仅spline在《MATLAB参考指南》中有说明。

下面几节，将展示在M文件函数中实现三次样条的基本特征。

12.1 基本特征在三次样条中，要寻找三次多项式，以逼近每对数据点间的曲线。

在样条术语中，这些数据点称之为断点。

因为，两点只能决定一条直线，而在两点间的曲线可用无限多的三次多项式近似。

因此，为使结果具有唯一性。

在三次样条中，增加了三次多项式的约束条件。

通过限定每个三次多项式的一阶和二阶导数，使其在断点处相等，就可以较好地确定所有内部三次多项式。

此外，近似多项式通过这些断点的斜率和曲率是连续的。

然而，第一个和最后一个三次多项式在第一个和最后一个断点以外，没有伴随多项式。

因此必须通过其它方法确定其余的约束。

最常用的方法，也是函数spline所采用的方法，就是采用非扭结(not-a-knot)条件。

这个条件强迫第一个和第二个三次多项式的三阶导数相等。

对最后一个和倒数第二个三次多项式也做同样地处理。

基于上述描述，人们可能猜想到，寻找三次样条多项式需要求解大量的线性方程。

实际上，给定N个断点，就要寻找N-1个三次多项式，每个多项式有4个未知系数。

这样，所求解的方程组包含有4*(N-1)个未知数。

把每个三次多项式列成特殊形式，并且运用各种约束，通过求解N个具有N个未知系数的方程组，就能确定三次多项式。

这样，如果有50个断点，就有50个具有50个未知系数的方程组。

幸好，用稀疏矩阵，这些方程式能够简明地列出并求解，这就是函数spline所使用的计算未知系数的方法。

12.2 分段多项式在最简单的用法中，spline获取数据x和y以及期望值xi，寻找拟合x和y的三次样条内插多项式，然后，计算这些多项式，对每个xi的值，寻找相应的yi。

数值分析三次样条插值函数【问题】对函数f x =ex, x∈[0,1]构造等距节点的三次样条插值函数，对以下两种类型的样条函数1. 三次自然样条2. 满足S′ 0 =1,S′ 1 =e的样条并计算如下误差：max{ f x1 −S x1 ,i=1,…,N} i−i−i这里xi−1为每个小区间的中点。

对N=10,20,40比较以上两组节点的结果。

讨论你的结果。

【三次样条插值】在每一个区间[t1,t2],…,[tn−1,tn]上，S都是不同的三次多项式，我们把在[ti−1,ti]上表示S的多项式记为Si，从而，S0 x x∈[t0,t1]∈[t1,t2] S x = S1 x x…Sn−1 x x∈[tn−1,tn]通过在节点处函数值、一阶导数和二阶导数的连续性可以得到：Si−1 ti = yi= Si ti 1≤i≤ n−1Si−1′ ti = Si′ tix→ti+limS′′ x =zi=limS′′(x) x→ti−再给定z0和zn 的值就构成了4n个条件，而三次样条插值函数共4n个系数，故可以通过这4n个条件求解三次样条函数的系数，从而求得该三次样条插值函数。

特别的，当z0=zn=0 时称为自然三次样条。

文本预览：一、自然三次样条插值【自然三次样条插值算法】1.由上面的分析可知，求解三次样条函数实际上就是求解一个矩阵：u 1h 1h1u2h2h2u3…v1 z1 v2 z2 z3=v3 … z…hn−2 n−2 vn−2 z vn−1 un−1 n−1ih3…hn−3un−2hn−26…其中hi=ti+1−ti，ui=2(hi+hi−1)，ui=h(yi+1−yi)，vi=bi−bi−1 所以自然三层次样条插值的算法就是在得到端点的函数值，一次导数值和二次导数值，然后根据上述求解矩阵得到v，代入自然三次样条的表达式即可。

2.根据题目中所给出的误差估计，计算在区间中点处的最大误差。

【实验】通过Mathematica编写程序得到如下结果：N=101. 计算得到zi的值为：由此可以得到各个区间的自然三次样条插值函数。

5.7 三次样条函数在制造船体和汽车外形等工艺中传统的设计方法是，首先由设计人员按外形要求，给出外形曲线的一组离散点值，施工人员准备好有弹性的样条（一般用竹条或有弹性的钢条）和压铁，将压铁放在点的位置上，调整竹条的形状，使其自然光滑，这时竹条表示一条插值曲线，我们称为样条函数。

从数学上看，这一条近似于分段的三次多项式，在节点处具有一阶和二阶连续微商。

样条函数的主要优点是它的光滑程度较高，保证了插值函数二阶导数的连续性，对于三阶导数的间断，人类的眼睛已难以辨认了。

样条函数是一种隐式格式，最后需要解一个方程组，它的工作量大于多项式拉格朗日型式或牛顿型式等显式插值方法。

定义给定区间上个节点和这些点上的函数值，。

若满足；在每个小区间上至多是一个三次多项式；在上有连续的二阶导数，则称为关于剖分的三次样条插值函数，称为样条节点。

要在每个子区间上构造三次多项式，共需要个条件，由插值条件，提供了个条件；用每个内点的关系建立条件又得到个条件；再附加两个边界条件，即可惟一确定样条函数了。

用待定系数法确定了构造样条函数的存在性和惟一性。

在具体构造样条函数时一般都不使用计算量大的待定系数法。

下面给出构造三次样插值的关系式和关系式的方法。

5.7.1 三次样条插值的M关系式引入记号。

用节点处二阶导数表示样条插值函数时称为大关系式，用一阶导数表示样条插值函数时称为小关系式。

问题5.8给定插值点，怎样构造用二阶导数表示的样条插值函数，即怎样构造关系式？假设。

由于在上为线性函数，故在上做的分段线性插值函数：令，得到（5.29）对积分两次有（5.30）将代入式（5.27）可解出故在上有（5.31）在每个小区间上具有不同的表达式，但由于在整个区间上是二阶光滑的，故有列出每一个关系式，再经计算得：（5.32）其中：由式（5.32）得到个未知数的个方程组。

现补充两个边界条件，使方程组只有惟一解。

下面分三种情况讨论边界条件。

（1）给定的值时，称为自然边界条件），此时阶方程组有个未知量，即（5.33）（2）给定的值，它们分别代入在中的表达式，得到另外两个方程：于是需要解阶的方程组：（5.34）（3）被插函数以为基本周期时，即，即；即。

样条插值是一种工业设计中常用的、得到平滑曲线的一种插值方法，三次样条又是其中用的较为广泛的一种。

1. 三次样条曲线原理假设有以下节点1.1 定义样条曲线是一个分段定义的公式。

给定n+1个数据点，共有n个区间，三次样条方程满足以下条件：a. 在每个分段区间（i = 0, 1, …, n-1，x递增），都是一个三次多项式。

b. 满足（i = 0, 1, …, n ）c. ，导数，二阶导数在[a, b]区间都是连续的，即曲线是光滑的。

所以n个三次多项式分段可以写作：，i = 0, 1, …, n-1其中ai, bi, ci, di代表4n个未知系数。

1.2 求解已知:a. n+1个数据点[xi, yi], i = 0, 1, …, nb. 每一分段都是三次多项式函数曲线c. 节点达到二阶连续d. 左右两端点处特性（自然边界，固定边界，非节点边界）根据定点，求出每段样条曲线方程中的系数，即可得到每段曲线的具体表达式。

插值和连续性：, 其中i = 0, 1, …, n-1微分连续性：, 其中i = 0, 1, …, n-2样条曲线的微分式：将步长带入样条曲线的条件：a. 由(i = 0, 1, …, n-1)推出b. 由(i = 0, 1, …, n-1)推出c. 由(i = 0, 1, …, n-2)推出由此可得：d. 由(i = 0, 1, …, n-2)推出设，则a. 可写为：，推出b. 将ci, di带入可得：c. 将bi, ci, di带入(i = 0, 1, …, n-2)可得：端点条件由i的取值范围可知，共有n-1个公式，但却有n+1个未知量m 。

要想求解该方程组，还需另外两个式子。

所以需要对两端点x0和xn的微分加些限制。

选择不是唯一的，3种比较常用的限制如下。

a. 自由边界(Natural)首尾两端没有受到任何让它们弯曲的力，即。

具体表示为和则要求解的方程组可写为：b. 固定边界(Clamped)首尾两端点的微分值是被指定的，这里分别定为A和B。

三次样条曲线表达式三次样条曲线是一种常用的插值方法，用于在给定的数据点之间进行平滑插值。

三次样条曲线由一组三次多项式组成，在每个相邻数据点之间使用不同的三次多项式来实现平滑曲线的连接。

这些三次多项式在相邻数据点处具有连续的一阶和二阶导数，以确保曲线的平滑性和连续性。

要表示三次样条曲线的表达式，首先需要确定数据点的坐标和插值方法。

假设我们有n个数据点 (x_0, y_0), (x_1, y_1), ..., (x_n, y_n)，我们可以使用三次样条插值方法来计算出每个相邻数据点之间的三次多项式表达式。

具体而言，三次样条曲线由n-1个三次多项式组成，每个三次多项式在相邻数据点之间定义。

假设第i个三次多项式在区间[x_{i-1}, x_i] 上的表达式为 S_i(x)，其中 i = 1, 2, ..., n-1。

在每个区间 [x_{i-1}, x_i] 上，三次样条曲线的表达式为：S_i(x) = a_i(x x_i)^3 + b_i(x x_i)^2 + c_i(x x_i) + d_i.其中 a_i, b_i, c_i, d_i 是待定系数。

为了确定这些系数，我们需要满足一些插值条件，例如在每个数据点处函数值相等，以及相邻区间的三次多项式在相邻数据点处的一阶和二阶导数相等。

通过解这些条件，可以得到每个区间上的三次多项式的系数，从而得到整个三次样条曲线的表达式。

最终的三次样条曲线表达式将是所有区间上三次多项式的组合。

总之，三次样条曲线的表达式包括了每个相邻数据点之间的三次多项式表达式，通过满足插值条件来确定每个区间上的系数，从而得到整个曲线的表达式。

这样的表达式能够实现对给定数据点之间的平滑插值，从而得到连续且平滑的曲线。

matlab 三次样条函数
《Matlab三次样条函数》
Matlab三次样条函数是由三次多项式函数组成的连续函数，根
据指定函数值，可以计算得到每一段函数的三次多项式系数。

经常用于数值计算中。

1、定义
Matlab中的三次样条函数（Cubic Spline）是一种连续函数，
由于其在每一段函数上都是三次多项式函数的拼接，因此又被称为三次拉格朗日插值法。

2、原理
Matlab中的三次样条函数可以用以下公式表示：
y=ax^3+bx^2+cx+d
其中，x为变量，y为函数值，a、b、c、d四个参数可以由指定的函数值计算出来。

首先，需要指定函数数据，一般是步长一致的函数值，即：
{x_1,y_1,x_2,y_2,x_3,y_3,...}
然后，根据以上函数值可以计算出每一段三次多项式函数的系数，即：
a_1,b_1,c_1,d_1
a_2,b_2,c_2,d_2
...
最后，函数总值可以由以上每一段的三次多项式函数拼接而成，
即：
y(x)=
begin{align}
t&a_1x^3+b_1x^2+c_1x+d_1, xin[x_1,x_2]
t&a_2x^3+b_2x^2+c_2x+d_2, xin[x_2,x_3]
t& ....
end{align}
3、应用
Matlab三次样条函数经常用于数值计算中，比如：（1）求解积分；
（2）求解微分方程；
（3）近似计算目标函数；
（4）曲线拟合等。

沈阳航空航天大学数学软件课程设计（设计程序）题目三次样条插值函数班级 / 学号学生姓名指导教师沈阳航空航天大学课程设计任务书课程名称数学软件课程设计院（系）理学院专业信息与计算科学班级学号姓名课程设计题目三次样条插值函数课程设计时间: 2010 年12月20日至2010 年12月31日课程设计的内容及要求：1.三次样条插值函数给出函数在互异点处的值分别为。

(1)掌握求三次样条插值函数的基本原理；(2)编写程序求在第一边界条件下函数的三次样条插值函数；(3)在区间上取n=10，20，分别用等距节点对函数作三次样条插值函数，利用(1)的结果画出插值函数的图形，并在该图形界面中同时画出的图形。

[要求]1.学习态度要认真，要积极参与课程设计，锻炼独立思考能力；2.严格遵守上机时间安排；3.按照MATLAB编程训练的任务要求来编写程序；4.根据任务书来完成课程设计论文；5.报告书写格式要求按照沈阳航空航天大学“课程设计报告撰写规范”；6.报告上交时间：课程设计结束时上交报告；7.严谨抄袭行为。

指导教师年月日负责教师年月日学生签字年月日沈阳航空航天大学课程设计成绩评定单课程名称数学软件课程设计院（系）理学院专业信息与计算科学课程设计题目三次样条插值函数学号姓名指导教师评语：课程设计成绩指导教师签字年月日目录一正文 (1)1问题分析 (1)1.1 题目 (1)1.2 分析 (1)2 研究方法原理 (1)2.1 求三次样条插值多项式，算法组织 (1)3 算例结果 (3)二总结 (7)参考文献 (8)附录 (9)源程序: (9)程序1 (9)程序2 (10)程序3 (12)程序 4 (12)一 正文1问题分析 1.1 题目三次样条插值函数给出函数在互异点处的值分别为。

（1）掌握求三次样条插值函数的基本原理；（2）编写程序求在第一边界条件下函数的三次样条插值函数； （3）在区间上取n=10，20，分别用等距节点对函数作三次样条插值函数，利用(1)的结果画出插值函数图形，并在该图形界面中同时画出的图形。

一、实验目的和要求（1） 掌握三次样条方法的原理，加深对插值方法的理解。

（2） 理解三次样条在高次差值的逼近效果比多项式的优越性！（3） 能构造出正确的算法结构并能实现具体问题的计算。

二、实验内容和原理实验内容：三次样条插值 解决具体问题：对于函数 ，取等距节点5,0,1,,10i i i x =-+= ，设已给出节点上的函数值以及左右两个端点的一阶导数值，按上述算法进行样条插值。

实验原理：（1）、三次样条的定义：设在[a,b]上有n+1个节点，f(xi)= yi,i=0,1,2,…,n 。

若S3(x)满足1) S3(x)在每个[xi,xi+1](i=0,1,…,n-1)上是不高于三次的多项式；2) S3(x)， ， 在[a,b]上连续；3) S3(xi)=yi(i=0,1,…,n)。

则称S3(x)为f(x)关于节点x0,x1,…,xn 的三次样条插值函数。

（2）、三次样条插值函数的边界条件：待定系数个数：4n已知条件：补充条件：这两个条件通常在区间[a,b]的两个端点给出，称为边界条件（3）、三次样条插值函数的求法：21()1f x x =+)('3x s )(''3x s 3333333(),(0,1,)(0)(0),(1,2,1)(0)(0),(1,2,1)(0)(0),(1,2,1)i i i i ii i i S x y i n S x S x i n S x S x i n S x S x i n ==⎧⎪-=+=-⎪⎨''-=+=-⎪⎪''''-=+=-⎩ 3011011i i+1i i i+120212021()()()()()...(1)h =x -x , x x x ()(1)(21)()(23)()(1)()(1)i i i i i i i i i i i i i i x x x x x x x x S x y y h m h m h h h h x x x x x x x x x x x x ϕϕψψϕϕψψ++----=+++≤≤=-+=-+=-=-得到系数矩阵，利用追赶法求解利用式（1）求解（4）、追赶法：追的过程（消元过程）赶的过程（回代过程）1,2,1][32)(322 111111111112121111111-=-++-+=++++⇔-+-=++++-----+----+⨯+--+----n i h y y h h h h y y h h h m h h h m m h h h h y y h y y h m m h m m ii i i i i i i i i i i i i i i i i i i i h h h h i i i i i ii i i i i i ii i i1,2,142624611211121-=++--=+------++n i h m m h y y h m m h y y i i i i i i i i i i i i ，)3......(..........2)1(2)1(2)1()1(211121212232232212011211⎪⎪⎪⎩⎪⎪⎪⎨⎧'-=+-=++-=++---=+-----------n n n n n n n n n n n n y αm m αm αm m αm αm m αy αm αm ββββ ⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡---=---210212102122221n n n A αααααα三、算法框图与程序框图（1）、算法框图：（2）、程序框图三次样条追赶法四、实验数据#include <iostream> #include<iomanip>using namespace std; double x[100]={0}; double y[100]={0}; int N;int n;double M0,Mn; double c[100]={0}; double p[100]={0};double h[100]={0};double d[100]={0};double a[100]={0};double X;double G0,G1,W0,W1,S3;int main(){//输入数据cout<<"输入节点的数目";cin>>N;n=N-1;cout<<"输入各个节点以及节点的值";for(int i=0;i<=n;i++){cout<<endl<<setw(15)<<"x["<<i<<"]=";cin>>x[i];cout<<setw(15)<<"y["<<i<<"]=";cin>>y[i];}cout<<"输入两个端点的一阶导数值";cout<<endl<<setw(15)<<"M0=";cin>>M0;cout<<setw(15)<<"Mn=";cin>>Mn;//用（39）式计算值for(int i=0;i<n;i++)h[i]=x[i+1]-x[i];for(int j=1;j<n;j++){c[j]=h[j-1]/(h[j-1]+h[j]);p[j]=3*((1-c[j])*(y[j]-y[j-1])/h[j-1]+c[j]*(y[j+1]-y[j])/h[j]); // cout<<p[j]<<" ";}//用追赶法求解//表示d[i];d[0]=M0;d[1]=p[1]-(1-c[1])*M0;for(int i=2;i<=n-2;i++){d[i]=p[i];}d[n-1]=p[n-1]-c[n-1]*Mn;d[n]=Mn;//表示a[i];for(int j=2;j<=n-1;j++){a[j]=1-c[j];}//按照框图6-4计算,其中b[i]=2;double b=2;c[1]=c[1]/b;d[1]=d[1]/b;int k=0;double t=0;for(k=2;k<=n-2;k++){t=b-c[k-1]*a[k];c[k]=c[k]/t;d[k]=(d[k]-d[k-1]*a[k])/t;}d[n-1]=(d[n-1]-d[n-2]*a[n-1])/(b-c[n-2]*a[n-1]);for(int i=n-2;i>=1;i--){d[i]=d[i]-c[i]*d[i+1];}cout<<"各节点的导数值如下："<<endl;for(int i=0;i<N;i++)cout<<"M["<<i<<"]="<<d[i]<<endl;cout<<endl<<"输入要求值点的X值：X="; cin>>X;int E=0;while(E<=n){if(E>=n){cout<<"输入有误";return 0;}if(x[E]<X&&X<x[E+1]){break;}elseE=E+1;}double R= (X-x[E])/h[E];G0=(R-1)*(R-1)*(2*R+1);G1=R*R*(-2*R+3);W0=R*(R-1)*(R-1);W1=R*R*(R-1);// cout<<n;// cout<<endl<<E;//cout<<"hahahhah";//cout<<G0<<endl<<G1<<endl<<W0<<endl<<W1<<endl<<(X-x[E])/h[E]<<endl<<d[E]<<endl <<d[E+1];S3=G0*y[E]+G1*y[E+1]+h[E]*W0*d[E]+h[E]*W1*d[E+1];cout<<endl<<S3;}/*4-1-11133561*//*11-5 0.03846-4 0.05882-3 0.10000-2 0.20000-1 0.500000 1.000001 0.500002 0.200003 0.100004 0.058825 0.038460.01479-0.01479*/。

§8 三次样条插值问题的提出:上面讨论的分段低次插值函数都有一致收敛性，但光滑性较差，对于像高速飞机的机翼形线，船体放样等型值线往往要求有二阶光滑度，即有二阶连续导数，早期工程师制图时，把富有弹性的细长木条（所谓样条）用压铁固定在样点上，在其他地方让它自由弯曲，然后画下长条的曲线，称为样条曲线。

它实际上是由分段三次曲线并接而成，在连接点即样点上要求二阶导数连续，从数学上加以概括就得到数学样条这一概念。

下面我们讨论最常用的三次样条函数。

三次样条函数:定义:函数],[)(2b a C x S ∈，且在每个小区间],[1+j j x x 上是三次多项式，其中b x x x a n =<<<=L 10是给定节点，则称)(x S 是节点n x x x ,,,10L 上的三次样条函数。

若在节点j x 上给定函数值),,1,0)((n j x f y j j L ==，并成立 ),,1,0()(n j y x S j j L ==，则称)(x S 为三次样条插值函数。

从定义知要求出)(x S ，在每个小区间],[1+j j x x 上要确定4个待定系数，共有n 个小区间，故应确定n 4个参数。

根据)(x S 在],[b a 上二阶导数连续，在节点)1,,2,1(−=n j x j L 处应满足连续性条件)0()0(+=−j j x S x S ，''(0)(0).j j S x S x −=+,).0()0(+′′=−′′j j x S x S 共有33−n 个条件，再加上)(x S 满足插值条件),,1,0()(n j y x S j j L ==，共有24−n 个条件，因此还需要2个条件才能确定)(x S 。

通常可在区间],[b a 端点n x b x a ==,0上各加一个条件（称为边界条件），可根据实际问题的要求给定。

常见的有以下三种：1° 已知两端的一阶导数值，即 n n f x S f x S ′=′′=′)(,)(00.2° 两端的二阶导数已知，即 ''00(),()n n S x f S x f ′′′′′′==， 其特殊情况 0)()(0=′′=′′n x S x S ， 称为自然边界条件。

三次样条函数的自动求法
夏省祥;于正文
【期刊名称】《山东建筑大学学报》
【年(卷),期】2003(018)004
【摘要】在Matlab的The Spline Toolbox中,没有给出三次样条函数表达式的求法,可在教学过程中,或在实际问题中,我们需要知道样条插值函数的分段表达式.在现行数值分析教材中,一般都是通过解方程确定三次样条插值函数的表达式,但这种方法的工作量很大.在本文中,我们用Matlab语言编制了三个小程序,给出在三种边界条件下,三次样条插值函数表达式的自动求法.
【总页数】4页(P86-89)
【作者】夏省祥;于正文
【作者单位】山东建筑工程学院,数理系,山东,济南,250014;山东建筑工程学院,数理系,山东,济南,250014
【正文语种】中文
【中图分类】O24
【相关文献】
1.法显西行求法青岛登岸及佛经译传刍议——纪念法显西行求法登陆归国1600周年 [J], 翟广顺
2.利用函数思想解释数列通项公式求法——以《一类数列通项公式的求法》一课教学为例 [J], 刘铁龙
3.在建构和演进中寻求法治进路,在平衡中寻求法治之美——评《国家机关权责平
衡问题研究》 [J], 方慧
4.椭圆曲率半径的数学求法与物理求法的比较 [J], 李景琴
5.关于二元三次样条函数空间的维数 [J], 罗炯兴;王华桥
因版权原因，仅展示原文概要，查看原文内容请购买。