当前位置：文档之家 › 三次样条函数

三次样条函数

  • 格式：doc
  • 大小：34.50 KB
  • 文档页数：4

计算方法实验报告
1、实验题目
三次样条插函数。

2、实验内容
三次样条插值是建立在Hermite 插值的基础上的。

Hermite 插值是在一个区间上的插值，而三次样条插则是建立多个区间上插值，构造一个具有二阶光滑度的曲线，在求出给定点上对应的函数。

本实验就是建立一个能根据三次样条插值函数求根的程序。

3、算法思想
给定一个区间，并把它分成n 等份，并且给出了每个结点对就的横坐标和纵坐标。

利用程序输出给定插值点对应的值。

横坐标设为：X 0, X 1, X 2, X 3, …X n 纵坐标为Y 0, Y 1, Y 2, …Y n ,设插点为u 。

则令h k =X k+1-X k ,λk =1-+k k k h h h , μk =11--+k k k h h h , g k =3(1
11--+-+-k k k k k k k k h y y h y y λμ), 其中k=1,2,…,n-1 再根据第一类边界条件则可以确定公式6.16，再根据6.17解出方程中的m 向量，最后代入公式6.8求解。

4、源程序清单
#include <stdio.h>
#define N 21/*最大结点个数减一*/
void sanCi()
{
/*定义过程数据变量*/
float x[N],y[N],h[N]; /*横纵坐标及区间长度*/
float rr[N],uu[N],gg[N]; /*计算m 用的中间数组rr 、uu 、gg 分别对应：λ、μ、g 数组*/
float aa[N],bb[N],tt[N]; /*矩阵分解时用到的中间变量aa、bb、tt分别对应：α、β数组以及A=LU时中间矩阵*/
float mm[N]; /*最后要用到的系数m*/
int n,k,kv,chose; /* n为实际结点个数,k为下标,kv为最后确定k的值*/
float s,u; /*最后计算u对应的值*/
printf("请输入区间段数：");
scanf("%d",&n); /*输入结点个数*/
/*输入所有横坐标：*/
printf("输入所有横坐标：");
for(k=0; k<=n; k++)
scanf("%f",&x[k]);
/*输入对应纵坐标：*/
printf("输入对应纵坐标：");
for(k=0; k<=n; k++)
scanf("%f",&y[k]);
for(k=0; k<n; k++) /*计算区间长度*/
h[k] = x[k+1] - x[k];
/*构建Am=f中的行列式A*/
for(k=1; k<n; k++)
{
rr[k] = h[k]/(h[k]+h[k-1]); /*计算λ[k]*/
uu[k] = h[k-1]/(h[k]+h[k-1]); /*计算μ[k]*/
gg[k] = 3 * (uu[k]*(y[k+1]-y[k])/h[k] + rr[k]*(y[k]-y[k-1])/h[k-1]); /*计算g[k]*/
}
/* g[1] = λ[1]*f＇[0] */
/* g[n-1] = μ[n-1]*f＇[n] */
gg[1] = gg[1] - rr[1]/x[0];
gg[n-1] = gg[n-1] -uu[n-1]/x[n];
/*作A=LU分解*/
aa[1] = 2;
bb[1] = uu[1]/aa[1];
for(k=2; k<n; k++)
{
aa[k] = (float)2.0 - rr[k]*bb[k-1];
bb[k] = uu[k]/aa[k];
}
/*通过LT＝G计算中间矩阵T*/
tt[1] = gg[1]/aa[1];
for(k=2; k<n; k++)
tt[k] = (gg[k]-rr[k]*tt[k-1])/aa[k];
/*根据第一类边界条件所得*/
mm[n] = 1/x[n];
mm[0] = 1/x[0];
/*通过UM=T计算M*/
mm[n-1] = tt[n-1];
for(k=n-2; k>=1; k--)
mm[k] = tt[k]-bb[k]*mm[k+1];
do{
printf("输入插值点：");
scanf("%f",&u);
/*判断u属于哪个区间x[k],x[k+1]确定k*/
for(k=0; k<n; k++)
{
if( u>=x[k] && u<=x[k+1] )
{
kv = k;
break;
}
}
/*代入u计算值*/
s = (1-2*(u-x[kv])/(-h[kv])) * ((u-x[kv+1])/(-h[kv])) * ((u-x[kv+1])/(-h[kv])) * y[kv];
s = s + (1-2*(u-x[kv+1])/h[kv]) * ((u-x[kv])/h[kv]) * ((u-x[kv])/h[kv]) * y[kv+1];
s = s + (u-x[kv])*((u-x[kv+1])/(-h[kv]))*((u-x[kv+1])/(-h[kv])) * mm[kv];
s = s + (u-x[kv+1])*((u-x[kv])/(-h[kv]))*((u-x[kv])/(-h[kv])) * mm[kv+1];
printf("与%f对应的结点值为%f\n\n",u,s);
printf("还要继续插值吗?(1/0): ");
scanf("%d",&chose);
}while( chose != 0 );
}
void main()
{
sanCi();
}
5、程序测试结果
运行程序：
请输入区间段数：3↙
输入所有横坐标：10 11 12 13↙
输入对应纵坐标：2.3026 2.3979 2.4849 2.5649↙输入插值点：11.75↙
与11.750000对应的结点值为2.463854
还要继续插值吗?(1/0):1↙
输入插值点：10↙
与10.000000对应的结点值为2.302600
还要继续插值吗?(1/0):0↙
正常退出程序。

合集下载

样条函数(三次样条)

样条函数(三次样条)

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

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。

三次样条函数的边界条件

三次样条函数的边界条件

三次样条函数的边界条件
一、第一种边界条件：固定边界条件（也叫夹紧边界条件）
1. 一阶导数边界条件
- 想象一下你在摆弄一个有弹性的东西，在两端你把它的“斜率”给定好了。

比如说，在区间的左端点a和右端点b，你规定了这个三次样条函数在这两点的一阶导数的值，就好像你在这两个端点抓住这个函数，让它按照你规定的倾斜程度开始和结束。

- 就好比你有一个弯曲的小棍，你在小棍的两端规定了它开始和结束弯曲的方向的陡峭程度。

2. 二阶导数边界条件
- 这个呢，就像是你在规定这个函数在边界处的“弯曲程度”。

在端点a和b，你设定了二阶导数的值。

这就好比你在控制小棍在两端的弯曲的曲率，是弯得更厉害还是比较平缓，你说了算。

二、第二种边界条件：自然边界条件
1. 二阶导数为零边界条件
- 这个边界条件很有趣哦。

它就像是你在说这个函数在边界处是比较“平滑”的，没有额外的弯曲。

就好像你在玩一个有弹性的绳子，在两端它没有被拧或者额外弯曲，是很自然地伸展着，所以二阶导数为零，就表示在边界处这个函数的弯曲没有突然的变化。

三、第三种边界条件：周期边界条件
1. 函数值和一阶导数相等边界条件
- 这就像是一个循环的东西。

想象你有一个环形的轨道，三次样条函数就沿着这个轨道走。

在区间的起点a和终点b，函数值是一样的，而且它在这两点的“倾斜程度”（一阶导数）也是一样的。

就好像你在这个环形轨道上，从一个点出发转一圈回来，函数的状态要能无缝对接，就像小火车在环形轨道上跑，回到起点的时候速度（一阶导数）和高度（函数值）都要和出发的时候一样。

第三讲三次样条函数

第三讲三次样条函数


给定区间[0, 上 例1 给定区间 3]上 3 个点的函数 值 f(0)=0, f(1)=2, f(3)=4, 试求数 a, b, c, d, 使函数 S(x)为给定点上的三次样条插值 为给定点上的三次样条插值 函数. 函数 其中 x2 + x + d , 0≤ x ≤1 S( x) = 3 . 2 ax + bx + cx + 1, 1 ≤ x ≤ 3
答案: 答案 a = −1, b = 4, c = −2, d = 0.
给定n+1个样点 i, yi )(i=0, 1, …, n), 个样点(x 给定 个样点 确定一个三次样条插值函数需要4n个独 确定一个三次样条插值函数需要 个独 立条件. 在定义中, 已指定了4n–2个条件 即 个条件, 立条件 在定义中 已指定了 个条件
yi Mi yi +1 Mi +1 + − hi ( xi +1 − x ) + − hi ( x − xi ) 6 6 hi hi
x ∈ [ x i , x i +1 ]
( xi +1 − x )2 ( x − xi )2 S ′( x ) = − Mi + Mi +1 2hi 2hi yi +1 − yi Mi +1 − Mi hi + − hi 6
M 0 d0 M d 1 1 M 2 d2 = M M λ n −1 M n −1 d n −1 2 M n dn
对于第1型插值问题 对于第 型插值问题: 型插值问题 λ 0 = 1, d 0 = 6 ( y1 − y 0 ) h0 − y 0′ h0 , µ n = 1, d n = 6 y n′ − ( y n − y n − 1 ) hn − 1 hn − 1 . 对于第2型插值问题 型插值问题: 对于第 型插值问题 λ 0 = 0, d 0 = 2 y 0′′ , µ n = 0, d n = 2 y n′′ . 对于第3型插值问题 对于第 型插值问题: 型插值问题
三次样条插值函数

三次样条插值函数

三次样条插值函数MATLAB 编程实现三次样条插值函数为()()[)()[]1011,,,,n n n S x x x x S x S x x x x-⎧∈⎪=⎨⎪∈⎩ 利用三次埃尔米特插值函数表示三次样条插值函数，即()()()()())111111,,j j j j j j j j j j j S x y x y x m x m x x x x ααββ++++++⎡=+++∈⎣（0,1,,1j n =-）基函数满足()()()()()()21112111121121111212jj j j j j j j j j j j j j j j j j j jj j j j x x x x x x x x xx xx x x x x x xx xx x x x xx x x x x x xααββ++++++++++++⎛⎫⎛⎫--=+ ⎪⎪ ⎪⎪--⎝⎭⎝⎭⎛⎫⎛⎫--=+ ⎪⎪ ⎪⎪--⎝⎭⎝⎭⎛⎫-=-⎪ ⎪-⎝⎭⎛⎫-=-⎪ ⎪-⎝⎭由上式易得()()()()()()()()()()()()()()1331111331112211112211612612246246j j j j j j j j j j j j j j j j j j j j j j jj j j j j x x x x x x xx x x x x x x x x x x x x xx xx x x x x x x x x ααββ+++++++++++++++''=---+''=-+--+''=---+''=---则有()()()()()()()()()()()111111113333111111122221111661212242466j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j S x y x y x m x m x x x x x y x y x x x x x x x x x x x x x m x m x x x x x x x x x ααββ+++++++++++++++++++''''''''''=+++⎡⎤⎡⎤++⎢⎥⎢⎥=-+-+⎢⎥⎢⎥----⎣⎦⎣⎦⎡⎤++⎢⎥+-+-⎢⎥----⎣⎦)1,j j x x x +⎡⎤⎢⎥⎡∈⎣⎢⎥⎣⎦（0,1,,1j n =-）同理有()()()()()()()()()()()()()()()11111113333111111122221111661212242466j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j S x y x y x m x m x x x x x y x y x x x x x x x x x x x x x m x m x x x x x x x x x ααββ------------------''''''''''=+++⎡⎤⎡⎤++⎢⎥⎢⎥=-+-+⎢⎥⎢⎥----⎣⎦⎣⎦⎡⎤⎡++⎢⎥+-+-⎢⎥----⎣⎦⎣)1,j j x x x -⎤⎢⎥⎡∈⎣⎢⎥⎦（1,,j n =）根据样条函数二阶导数连续性，即()()100j j j j S x S x +''''+=-（1,,1j n =-）即()()()()()()()()()()()()()()()()111111332211111111113322111166426624j j jj j j j j j j jj j j j j j j j j j j j j j jj j j j j jj j j j j j j j x x y x x y x x x x m m x x x x x x x x x x y x x y x x x x m m x x x x x x x x ++++++++++--------------+++--------=+++----（1,,1j n =-）化简得()()()()()111111111111233j j j j j j j j j j j j j j j j j j jj j xx m x x m x x m x x x x y y y y x x x x +-+--+-++-+--+-+---=-+---（1,,1j n =-）可得线性方程组()()()()()()()()()()0121201023231213121221111110212110211032213221322122233333n n n n n n n n n n n n m m x x x x x x m x x x x x x m x x x x x x m m m x x x x y y y y x x x x x x x x y y y y x x x x y ------⨯+-+⨯⎛⎫ ⎪ ⎪ ⎪---⎛⎫⎪ ⎪--- ⎪ ⎪⎪ ⎪⎪ ⎪ ⎪--- ⎪⎝⎭ ⎪ ⎪ ⎪⎝⎭---+------+---=()()()121112112113n n n n n n n n n n n n n x x x x y y y x x x x ----------⨯⎛⎫⎪ ⎪ ⎪ ⎪⎪⎪ ⎪-- ⎪-+- ⎪--⎝⎭为了使样条插值问题有惟一解，我们在原有方程基础上增加两个边界条件。