三次样条曲线的生成算法的研究
三次样条曲线的生成算法
本文由天空乐园河南自考网整理分享
摘要
三次样条函数曲线具有的最高多项式插值精度是三次多项式函数，对其进行推广构造的三次参数样条曲线应至少具有同样的插值精度。
本文讨论了构造三次参数样条曲线中节点选取问题，相邻两节点之间的跨度规范化为1，提出了构造2GC 三次参数样条曲线的新方法。文中首先讨论了2GC 三次参数样条曲线需满足的连续性方程，然后讨论了平面有序五点确定一组三次多项式函数曲线和平面有序六点唯一确定一条三次多项式函数曲线。在此基础上，提出了为给定数据点选取节点值的新方法。新方法构造的2GC 三次参数样条曲线具有三次多项式函数的插值精度。最后以具体数据点对新方法和已有的四种节点选取方法构造的插值曲线的精度做了比较。
关键词：三次样条曲线；曲线拟合；计算机图形学
自1946年美国数学家I. J. Schoenberg 提出样条函数[1]以来，样条函数以其构造简单、易于计算又有很好的力学背景等特点而被广泛用于科学计算、工程设计和计算机辅助设计等领域，成为最重要的曲线和曲面构造方法之一。在样条函数的应用中，三次样条函数由于具有极小模性质、最佳逼近性质和很强的收敛性[2,3,4]等而成为最主要的方法应用于构造插值曲线和曲面。
用样条函数方法构造三次插值曲线，曲线的连续性基本可满足实际应用的要求。当曲线的端点条件确定之后，曲线的精度和形状是由曲线需满足的连续性方程唯一决定的。在小挠度的情况下，插值曲线的精度和形状都是非常理想的。对大挠度曲线和任意平面数据点，则需推广三次样条函数方法构造三次参数样条曲线，此时需知道每个数据点处的参数值（节点值）。在实际应用中，这些参数值一般是无法预先给定的，所以构造三次参数样条曲线的第一步是对给定数据点参数化，即为每个数据点指定节点值。如果指定的节点值是精确的，给定适当的端点条件，可使构造的插值曲线的代数精度达到三次参数多项式。构造三次参数样条曲线，当曲线的端点条件确定之后，能够决定曲线插值精度的量只有节点。因此构造三次参数样条曲线的关键是如何选择节点。目前常用的节点选取方法有4种，均匀参数化法、累加弦长参数化法、向心参数化法[5]和修正弦长参数化法[6]。这些方法虽然在实际中得到了较为广泛的应用，但从逼近的角度看，它们的插值精度较低，其插值多项式的最高精度是线性的。最近一个确定节点的方法[7]具有二次多项式插值精度，如果用来构造三次参数样条曲线，这个精度也是较低的。
三次样条函数曲线具有的最高多项式插值精度是三次多项式函数，对其进行推广构造的三次参数样条曲线应至少具有同样的插值精度。从这一目标出发，本文讨论了构造三次参数样条曲线中节点选取问题，相邻两节点之间的跨度规范化为1，提出了构造2GC 三次参数样条曲线的新方法。文中首先讨论了2GC 三次参数样条曲线需满足的连续性方程，然后讨论了平面有序五点确定一组三次多项式函数曲线和平面有序六点唯一确定一条三次多项式函数曲线。在此基础上，提出了为给定数据点选取节点值的新方法。新方法构造的2GC 三次参数样条曲线具有三次多项式函数的插值精度。最后以具体数据点对新方法和已有的四种节点选取方法构造的插值曲线的精度做了比较。
平面自由曲线——不能用一个标准代数方程精确表示。实际中应用很多，如轮船船身放样。
将放样过程抽象为:平面上给定若干点（型值点），找一个代数方程，逼近或插值上述型值点。
理论上，n个点，可以找到一个n-1次多项式来逼近，但n太大时，多项式次数太高，计算复杂，难以控制。
工程上，降低次数，且分段定义。
样条函数自提出以来，以其构造简单，易于计算，及很好的力学背景等特点被广泛用于科学计算，工程设计和计算机辅助设计等领域，从而成为最重要的曲线和曲面构造方法之一。
三次样条曲线在使用中存在局限性，且表示方法缺乏几何不变性。即当平面直角坐标系中得型值点发生旋转等几何变形时，其曲线的形状也发生变形，严重时甚至不能保证满足X1 在周期性三次样条曲线拟合控制多边形时，取得了很小的偏差；基于累加弦长的三次参数样条曲线插值在数控系统中取得了较好的效果，但是以累加弦长为参数的三次参数样条曲线插值和基样条的函数插值在各分段曲线两端曲率的符号相同的情况下都有可能产生这段曲线上的拐点，造成曲线不光顾。因此一些准测提出检查多余的拐点，YE J等人修正了Kjellander的方法，并从累加弦长参数化和光顾函数两方面消除了三次参数样条的震荡和回折。 在曲线拟合中，插值过程可具体使用线性（liner）插值，三系样条（spline）插值，立方(cubic)插值等方法，在曲线插值法中最常用的是线性插值法，它是估计2个主干点之间数值的最简单，最易实现的方法，但采用线性插值法会有以下缺点：1曲线不能显示连接主干点间的凸状弧线；2从曲线导出远期曲线时会形成人为的“尖头”。 因此，通常采用样条法来构造曲线，它通过构造多项式（1个或1组不同阶多项式）来形成1条把所有主干点连接起来的平滑曲线，一般常选择3次曲线（根据3次插值样条函数所得的曲线）进行拟合。3次样条曲线具有良好的数学特征，而且用3次曲线去拟合时，其结果要比线性插值估计更接近于工程实际情况，但是在工程应用中，我们利用三次样条插值方法，对相同的控制点，只可以得到1条光滑曲线，如果我们想基于相同的控制点，得到多条不同曲线，依靠传统算法，是无法实现的，这就限制了三次样条在工程中，特别是印刷领域中得应用。 在印刷领域，特别市在工艺前端，传统的方法是利用曲线进行分色操作，由于色彩的特殊性，即人眼对色彩的感觉不尽相同。如果只有1条色度曲线，那么工艺人员就无法对色彩效果进行有效的对比。因此印刷工艺的特殊性要求能够根据相同控制点，得出多条曲线，实现不同的印刷色彩效果，从中选出最佳的色度曲线。在这一点上，传统的方法是通过修改基本的控制点，生成新的控制曲线实行，本文提出改进的3次样条算法，实现了在相同的控制点上，生成了曲线不同的新曲线。增加了生成曲线的条数，从而使得印刷前端的工艺操作人员，对控制图像的色度曲线有更多的选择。 三次参数样条曲线的构造 设平面上给定了n个数据点目标是构造一条对n个数据点插值的三次参数样条曲线P (t)。设ti(待定)是与i P相对应的节点，令则区间上的三次参数曲线Pi(t)可定义如下： 其中 为[0,1]区间上的三次埃尔米特基函数，为节点ti处的切矢。P t ()应满足的连续性方程是[6,8]： 其中; 令； 则(2.2)可写成如下形式： 其中； 如果Si给定，则可得到(2.4)中的Ni,i=1,2,…,n,由 知，对i =1,2, …,n-1,(2.1)可写成如下形式： 其中0≤s≤1。 显然，由(2.6)定义的样条曲线是连续的。方程组(2.4)中有n ? 2方程，n个未知量,解方程组需增加两个端点条件，方法如下： (1)对封闭或周期曲线 所以， (2)两端点处的二阶导数为零 此时，对(2.4)中的i =2和n-1 ，相应的方程为 所求方程中的未知量为 (3)给定端点切向条件M1 和Mn ， 其中F1和F2的值确定将在第四节中讨论。 设，是给定数点中连续五点，与Pj 相 对应的节点是t j(待定)。为讨论方便，我们做如下变换。设不在一条直线上，把坐标变换 分别施加到Pi和ti 。则在ovw坐标系中的 坐标为 相应的节点值分别为和其中is 由 (2.3)定义。对插值的三次参数曲线 是 其中； 如果给定数据点是一条三次多项式函数曲线上的点，则vi(s)和wi(s)中二次项系数之比等于三次项系数之比[8]，即 化简得： 其中满足由于Q (s)是唯一的，从而有 直接验证知，两式是等价的，并且可写成如下形式 由( 3.4)知，这时一个关于Si 的一元五次方程，用公式法不能求出精确解。可用如下方法求精确解。考虑Pk ，k=i-1,i+3，把变换(3.1)施加到Pk ，由对称性得 其中由(3.4)定义。联立(3.6)和(3.7)可求出精确解Si 。因此平面六点可唯一定义一条三次多项式函数曲线。 对边界数据点Pi，当 3 = i ，相应于(3.7)的方程是 联立(3.6) 和(3.8) 可求出精确解S3对P2 , 由(3.4) 知，相应于P2 的参 数是所以，同理可求出 。 计算题 可用两种方法计算Si ：1) 对(3.6)直接用数值方法求；2) 联立(3.6)和(3.7)用公式法求。所求的Si应满足下列条件： 这样的Si可能多于一个。所以，平面有序五点不能唯一决定一条三次多项式 函数曲线。确定Si的方法是所有满足(3.9)的加权平均。 下面讨论权函数的取值。在所有中，只有一个值是所求的。考虑，如果给定数据点是三次多项式函数上的点，由(3.3)和(3.4)知，中应有一个Si,l满足 如果给定数据点和P k不是三次多项式上的点，则, 。在这种情况下，希望值小的Si,l 对形成Si的影响大。为此定义 Si由下式定义： 显然，如果Si,j ,满足Ei,j=0，则Si = Si,j。 记(3.6)的右端为F(Si)。计算表明，在大多数情况下，F(Si)有符合(3.9)的两个实根或两个重根，如图1所示。在图1 中，F(Si)虽然有三个实根，但左边的根不符合条件(3.9)。如果(3.6)中没有满足(3.9)的Si，则选取F(Si) 在[0,1]区间上的极小值点(F(Si)在[0,1]区间上的极值点的值大于零)或极大值点(F(Si)在[0,1]区间上的极值点的值小于零)作为Si。例如，对图2 所示的情况选取极大值点为Si。如果F(Si) 在[0,1]区间上既没有符合(3.9)的实根，也没有符合条件的极小或极大值点，如图 3 所示，则Si的选取应极小化下式的值 即，对 应使(3.3)定义的两个三次参数曲线的三次项系数的差最小。 例子： 本节我们以实例对新方法、累加弦长法、向心法，修正弦长法和二次精度法做比较。用于比较的数据点取自一条给定曲线。用五种方法分别对给定数据点构造三次参数样条曲线，以样条曲线的插值精度对五种方法进行比较。定义数据点的曲线是 其中， 三次参数曲线F(t) 具有如下性质：K = 1,2,3,4时是凸的，K = 5,6,7,8 时有两个拐点，K = 9时有一个尖点，K = 10 ,11,12 时有一个圈，图4 是K = 3,6,9,12 时F (s)在区间[0,1] 上的图形。 用于比较的区间是[0,1]。区间被分成20 个子区间定义数据点Pi=F(Ti) , i=0,1,2,…,20 ,Ti由下式定义 其中0≤λ≤0.25 。λ的取值使数据点相邻两弦长满足 五种方法用绝对误差曲线E(s) 比较,E(s) 定义如下 其中P(s) 表示五种方法构造的三次参数样条曲线之一，F(t) 为(5.1)或(5.3)定义的曲线，P i(s)表示P(s)在区间[Si,Si+1]上的部分，|P(s)-F(t)|表示点P(s)到F(t)的距离。 (5.2)中取λ= 0.25 时，五种方法产生的最大绝对误差见表1。 表 1 五种方法对(5.1)插值产生的最大误差 E(s) 新方法二次精度累加弦长修正弦长向心 K=1 9.07e-6 1.68e-5 3.91e-4 4.97e-4 8.07e-4 K=2 5.08e-6 9.24e-6 2.23e-5 2.77e-4 1.29e-3 K=3 1.73e-7 4.89e-6 2.42e-5 5.50e-4 1.87e-3 K=4 4.73e-5 1.47e-5 4.50e-5 9.40e-4 2.37e-3 K=5 1.01e-4 3.88e-4 5.43e-5 1.36e-3 2.70e-3 K=6 3.60e-4 6.21e-4 2.60e-4 1.68e-3 2.73e-3 K=7 4.33e-4 1.64e-3 1.44e-3 1.61e-3 2.06e-3 K=8 2.19e-4 3.08e-4 4.25e-3 1.07e-3 1.37e-3 K=9 2.68e-4 9.34e-4 1.11e-3 4.57e-4 1.84e-3 K=10 2.75e-4 4.25e-4 6.32e-3 1.93e-3 2.31e-3 K=11 1.72e-4 3.83e-4 3.15e-3 3.82e-3 3.70e-3 K=12 1.46e-4 3.79e-4 1.23e-3 5.17e-3 5.55e-3 下面用一个椭圆定义数据点对五种方法构造的闭合曲线进行比较。椭圆的方程是 X = 3cos(2πt) (5.3) Y = 2sin(2πt) 区间[0,1]也被分成20 个小区间定义数据点。本例中也对用精确节点值构造的三次参数插值曲线进行比较。对应于(5.2)中不同的λ，六种方法产生的最大绝对误差见表2。计算知，本例中λ= 0时新方法和二次精度法指定的节点都是精确的。表2说明，对本例，新方法构造的插值曲线的精度比用精确节点构造的 插值曲线的精度高。 表 2 五种方法对椭圆插值产生的最大误差 E(s) 精确参数新方法二次精度累加弦长修正弦长向心 λ=.0 7.68e-5 7.68e-5 7.68e-5 6.65e-4 1.87e-4 3.56e-4 λ=.05 1.36e-4 1.16e-4 2.08e-4 7.55e-4 3.46e-3 4.68e-3 λ=.10 2.20e-4 1.68e-4 3.85e-4 8.41e-4 7.48e-3 1.02e-2 λ=.15 3.25e-4 2.26e-4 6.03e-4 9.23e-4 1.24e-2 1.65e-2 λ=.20 4.51e-4 2.88e-4 8.62e-4 1.00e-3 1.83e-2 2.37e-2 λ=.25 6.00e-4 3.53e-4 1.17e-3 1.07e-3 2.50e-2 3.18e-2 我们还把区间[0,1]分成30、40 等个小区间对五种方法进行比较，其结果类似于表1和表2。 总结： 在平面上唯一确定一条三次多项式函数曲线需要6 个有序数据点。如果给定的 6 个有序数据点是三次多项式函数F 上的点，则6 个数据点可精确定义F 。 对平面上给定的n个数据点，本文提出了一个对其进行参数化的新方法。用指定节点构造的三次参数样条曲线是2GC 连续的，其多项式准确集包括所有三次和小于三次的函数多项式。实例计算也表明，新方法构造的三次参数样条曲线对给定数据点，尤其是对凸的数据点，具有较高的插值精度。 继续的工作是：研究是否存在一个对数据点参数化的方法，使构造的三次参数样条曲线的多项式准确集包括所有三次和小于三次的参数多项式。推广本文思想构造对空间数据点插值的2GC 三次参数样条曲线。 在对叠化的编辑模型分析的基础上，针对有些非叠化区域由于物体及摄像机运动的影响，也符合叠化编辑模型的特征，利用叠化所特有的连续一致性对候选收索区域进行确认，减少了因物体及摄像机运动所产生的误检，提高了算法的查准率。同时，本算法是在检测出切变镜头边界的基础上检测叠化区域，避免了所选取视频序列的边界恰为叠化区域造成的误检，提高了算法的查全率，实验中和其他算法的检测结果进行了比较，结果表明本文给出的算法查全率和查准率都较高。 参考文献 【1】Differential Synthetic Aperture Radar Interferometry DINSAR) for 3D Coastal Geomorphology ReconstructionMa ged Marghany?and Mazlan hashim 【2】Incorporating Prior Knowledge in Cubic Spline Approximation Application to the Identification of Reaction Kinetic Models János Madár 1, János Abonyi 1, Hans Roubos 2, Ferenc Szeifert 1 【3】自适应三次样条插值逼近算法研究——戴任权 【4】曲线曲面的B-3样条拟合及实现——丛伟 【5】计算机图形学抛物及三次样条曲线 【6】图形学_05_曲线之生成算法 【7】计算机图形学第五章曲线与曲面的生成 【8】基于DXF文件格式的三次参数样条曲线的生成——中国矿业大学——张剑英 【9】基于3次样条插值控制色度曲线算法分析——商云辉 【10】Hermite三次样条插值的车型机器人路径规划研究——彭辉 【11】高精度三次参数样条曲线的构造———张彩明 【12】自然边界条件下的三次样条计算——张汉茹 由开封市同力基础工程有限公司整理。 本文由天空乐园河南自考网分享。 例1 设)(x f 为定义在[0,3]上的函数，有下列函数值表： 且2.0)('0=x f ，1)('3-=x f ，试求区间[0,3]上满足上述条件的三次样条插值函数)(x s 本算法求解出的三次样条插值函数将写成三弯矩方程的形式： ) ()6()() 6()(6)(6)(211123 13 1j j j j j j j j j j j j j j j j x x h h M y x x h h M y x x h M x x h M x s -- + -- + -+ -= +++++其中，方程中的系数 j j h M 6， j j h M 61+，j j j j h h M y )6(2- ， j j j j h h M y ) 6(211++- 将由Matlab 代码中的变量Coefs_1、Coefs_2、Coefs_3以及Coefs_4的值求出。 以下为Matlab 代码： %============================= % 本段代码解决作业题的例1 %============================= clear all 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); for i = 1 : length(IndVar) - 1 h(i) = IndVar(i + 1) - IndVar(i); end % 为向量μ赋值 CENTRAL SOUTH UNIVERSITY 数值分析实验报告 三次样条插值方法的应用 一、问题背景 分段低次插值函数往往具有很好的收敛性，计算过程简单，稳定性好，并且易于在在电子计算机上实现，但其光滑性较差，对于像高速飞机的机翼形线船体放样等型值线往往要求具有二阶光滑度，即有二阶连续导数，早期工程师制图时，把富有弹性的细长木条（即所谓的样条）用压铁固定在样点上，在其他地方让他自由弯曲，然后沿木条画下曲线，称为样条曲线。样条曲线实际上是由分段三次曲线并接而成，在连接点即样点上要求二阶导数连续，从数学上加以概括就得到数学样条这一概念。下面我们讨论最常用的三次样条函数及其应用。 二、数学模型 样条函数可以给出光滑的插值曲线（面），因此在数值逼近、常微分方程和偏微分方程的数值解及科学和工程的计算中起着重要的作用。 设区间[]b ,a 上给定有关划分b x x n =<<<= 10x a ，S 为[]b ,a 上满足下面条件的函数。 ● )(b a C S ,2∈； ● S 在每个子区间[]1,+i i x x 上是三次多项式。 则称S 为关于划分的三次样条函数。常用的三次样条函数的边界条件有三种类型： ● Ⅰ型 ()()n n n f x S f x S ''0'',==。 ● Ⅱ型 ()()n n n f x S f x S ''''0'''',==，其特殊情况为()()0''''==n n x S x S 。 ● Ⅲ型 ()() 3,2,1,0,0==j x S x S n j j ，此条件称为周期样条函数。 鉴于Ⅱ型三次样条插值函数在实际应用中的重要地位，在此主要对它进行详细介绍。 三、算法及流程 按照传统的编程方法，可将公式直接转换为MATLAB 可是别的语言即可；另一种是运用矩阵运算，发挥MATLAB 在矩阵运算上的优势。两种方法都可以方便地得到结果。方法二更直观，但计算系数时要特别注意。这里计算的是方法一的程序，采用的是Ⅱ型边界条件，取名为spline2.m 。 Matlab 代码如下： function s=spline2(x0,y0,y21,y2n,x) %s=spline2(x0,y0,y21,y2n,x) %x0,y0 are existed points,x are insert points,y21,y2n are the second 第三章 样条插值和曲线拟合 1．x y = 有如下的函数表 8。 解 先作差商表 4 167 1210 13 9 3 42015 11008 16012 4 60 13 1611 1 10 0-?- -- 故：8.2)48(5 1 2)8(1=-+=p 819047619.2) 98)(48(210 1 )48(512)8(2=----+=p 844444.2)98)(48)(18(3 4201) 48)(18(601 )18(311)8(3=---?+----+=p 6222.2)1(4781008 1478601) 18(86 1 )08(10)8(4=-???-??+---?+=p 已知 828427.28=，因此选定 )8(,16,9,42321p x x x ===最接近8。 利用Neville 方法得: xi 8-xi f(xi) 2.8284271 8 0 8 1 7 1 -1.33333333 3.3333333 2.4 4 4 2 2.866666667 2.6222222 2.8 2.8444444 9 -1 3 2.819047619 2.8571429 16 -8 4 f(8)= 2.828427125 xi 8-xi f(xi) 8 0 8 1 7 1 -1 1/3 3 1/3 2 2/5 4 4 2 2 13/15 2 28/45 2 4/5 2 38/45 9 -1 3 2 86/105 2 6/7 16 -8 4 已知 828427.28=，故选定)8(,16,9,42321 p x x x ====2.819047619最接近8. 三次样条插值在工程拟合中的应用 摘要: 介绍了工程实验、勘测、设计中常见的列表函数之数值插值方法、程序实现及工程应用, 应用此法可方便地将任何列表函数计算到工程设计、施工所需要的精确程度, 给 出了各参数随主要参数变化而变化的光滑曲线, 并将其应用推广到一般情况. 关键词: 列表函数; 数值拟合; 三次样条插值; MA TLAB 程序设计与应用 在实际工程中, 广泛存在这样的问题: 根据设计要求和具体的工程条件, 在初始设计阶段会勘测得到若干组该工程的控制参数, 但这些参数之间彼此离散、不够密集, 利用它们来施工则不能满足施工的精度要求. 为了解决这一问题, 需要对已知的参数数据进行分析处理, 进行必要的插值、拟合, 以达到施工所需要的数据精度.本文以工程实例为基础, 对实际工程中插值方法的选取、插值的实现和插值曲线的拟合加以讨论, 提出能得到较合乎实际的插值方法, 给出一般工程人员就能实现的计算方法以及能得到光滑曲线的拟合方法. 1 工程应用实例 表1 所示的为某双曲拱坝体形原始参数[ 1 对于这一类工程列表参数有一个显著的特点:尽管不同工程的参数多寡不同, 但都是由n 行k 列的离散的列表数据给出, 虽然同一行代表某工程特定位置的几个参数(或高程参数, 或上游 半径参数?) , 但相邻两行由于位置距离太大, 两行各参数之间究竟存在什么数值关系, 对工 程设计、施工有何影响, 这是工程技术人员需要弄清楚的[ 2 ].以双曲拱坝为例, 它沿整个高程的变化是一个连续光滑的空间曲面. 从施工需要来看, 这些数据太稀疏, 难以满足设计、施工放样与钢筋配置等要求, 如果照此施工, 则有可能达不到工程精度、降低工程效率; 从计算机图形模拟来看, 要生成这个曲面仅由这一列表函数是得不到光滑曲面的, 是不可取的. 所以, 为使计算精确, 满足工程施工过程中任何断面位置、任意水平位置、任意高程位置所必需的施工数据与设计图纸, 保证工程施工的高品质,就要求作精确的数据处理.进一步分析可知, 在这 些参数表中, 各行的参数都随某一主要参数的变化而变化, 如上游半径参数随高程的变化而变化?, 它们的这种函数关系,在数值分析中有许多的方法可以求得. 但是哪种方法能更好、更合 乎实际地给出平滑曲线呢? 下面所选的插值方法能够较好地满足这一要求. 1设计目的、要求 对龙格函数2 2511 )(x x f += 在区间[-1,1]上取10=n 的等距节点，分别作多项式插值、三次样条插值和三次曲线拟合，画出)(x f 及各逼近函数的图形，比较各结果。 2设计原理 （1） 多项式插值：利用拉格朗日多项式插值的方法，其主要原理是拉格朗日多项 式，即： 01,,...,n x x x 表示待插值函数的1n +个节点， 0()()n n j k k j j k L x y l x y ===∑，其中0,1,...,j n =； 011011()...()()...() ()()...()...()...() k k n k k k k k k k n x x x x x x x x l x x x x x x x x x -+-+----= ---- （2） 三次样条插值：三次样条插值有三种方法，在本例中，我们选择第一边界条 件下的样条插值，即两端一阶导数已知的插值方法： 00'()'S x f = '()'n n S x f = （3）三次曲线拟合：本题中采用最小二乘法的三次多项式拟合。最小二乘拟合是 利用已知的数据得出一条直线或者曲线，使之在坐标系上与已知数据之间的距离的平方和最小。在本题中，n= 10，故有11个点，以这11个点的x 和 y 值为已知数据，进行三次多项式拟合，设该多项式为 23432xi i i i p a a x a x ax =+++，该拟合曲线只需2[]xi i p y -∑的值最小即可。 3采用软件、设备 计算机、matlab 软件 4设计内容 1、多项式插值： 在区间[] -上取10 1,1 n的等距节点，带入拉格朗日插值多项式中，求出各个节点的插值， = 并利用matlab软件建立m函数，画出其图形。 在matlab中建立一个lagrange.m文件，里面代码如下： %lagrange 函数 function y=lagrange(x0,y0,x) n=length(x0);m=length(x); for i=1:m z=x(i); s=0.0; for k=1:n p=1.0; for j=1:n if j~=k p=p*(z-x0(j))/(x0(k)-x0(j)); end end s=p*y0(k)+s; end y(i)=s; end 建立一个polynomial.m文件，用于多项式插值的实现，代码如下： %lagrange插值 x=[-1:0.2:1]; y=1./(1+25*x.^2); x0=[-1:0.02:1]; y0=lagrange(x,y,x0); y1=1./(1+25*x0.^2); plot(x0,y0,'--r') %插值曲线 hold on %原曲线 plot(x0,y1,'-b') 运行duoxiangshi.m文件，得到如下图形： MATLAB 程序设计期中考查 在许多问题中，通常根据实验、观测或经验得到的函数表或离散点上的信息，去研究分析函数的有关特性。其中插值法是一种最基本的方法，以下给出最基本的插值问题——三次样条插值的基本提法： 对插值区间[]b a ,进行划分：b x x x a n ≤ 学习报告—— 三次样条函数插值问题的讨论 班级：数学二班 学号：152111033 姓名：刘楠楠 样条函数： 由一些按照某种光滑条件分段拼接起来的多项式组成的函数；最常用的样条函数为三次样条函数，即由三次多项式组成，满足处处有二阶连续导数。 一、三次样条函数的定义： 对插值区间[,]a b 进行划分，设节点011n n a x x x x b -=<< <<=，若 函数2()[,]s x c a b ∈在每个小区间1[,]i i x x +上是三次多项式，则称其为三次样条函数。如果同时满足()()i i s x f x = (0,1,2)i n =，则称()s x 为()f x 在 [,]a b 上的三次样条函数。 二、三次样条函数的确定： 由定义可设：101212 1(),[,] (),[,]()(),[,] n n n s x x x x s x x x x s x s x x x x -∈??∈?=???∈?其中()k s x 为1[,]k k x x -上的三次 多项式，且满足11(),()k k k k k k s x y s x y --== (1,2,,k n = 由2()[,]s x C a b ∈可得：''''''()(),()(),k k k k s x s x s x s x -+-+== 有''1()(),k k k k s x s x -++= ''''1()(),(1 ,2,,1)k k k k s x s x k n -+ +==-， 已知每个()k s x 均为三次多项式，有四个待定系数，所以共有4n 个待定系数，需要4n 个方程才能求解。前面已经得到22(1)42n n n +-=-个方程，因此要唯一确定三次插值函数，还要附加2个条件，一般上，实际问题通常对样条函数在端点处的状态有要求，即所谓的边界条件。 1、第一类边界条件：给定函数在端点处的一阶导数，即 ''''00(),()n n s x f s x f == 2、第二类边界条件：给定函数在端点处的二阶导数，即 给定数据表如下： 试求三次样条插值S（X），并满足条件： i)S’(0.25)=1.0000, S’(0.53)-0.6868; ii) S”(0.25)= S”（0.53）=0； 解： 由给定数据知： h0 =0.3-0.25 - 0.05 , h 1=0.39-0.30-0.09 h 2=0.45-0.39-0.06, h 3=0.53-0.45-0.08 由μ i=h i/(h i1+h i), λ i= h i/(h i1+h i) 得: μ1= 5/14 ; λ 1= 9/14 μ2= 3/5 ; λ 2= 2/5 μ3= 3/7 ; λ 3=4/7 0.25 0.5000 ﹨ ﹨ 1.0000 ∕﹨ 0.25 0.5000 ∕ -0.9200-f[x 0,x 0, x 1 ] ﹨∕ 0.9540 ∕﹨ 0.30 0.5477 -0.7193-f[x 0,x 1,x 2 ] ﹨∕ 0.8533 ∕﹨ 0.39 0.6245 -0.5440-f[x1,x2,x 3 ] ﹨∕ 0.7717 ∕﹨ 0.45 0.6708 -0.4050-f[x 2,x 3,x 4 ] ﹨∕ 0.7150 ∕﹨ 0.53 0.7280 -0.3525-f[x 3,x 4,x 5 ] ﹨∕ 0.6868 ∕ 0.53 0.7280 i)已知一节导数边界条件，弯矩方程组 ┌┐┌┐ │ 2 1 │┌M 0 ┐│-0.9200 ︳ ︳5/14 2 9/14 ︳︳M ︳︳-0.7193 ︳ 1 ︳3/5 2 2/5 ︳︳M 2 ︳_6 ︳-0.5440︳ ︳ 3/7 2 4/7 ︳︳M ︳︳-0.4050 ︳ 3 计算方法实验报告 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 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 数学实验（三次样条） 数学实验（三次样条插值） 实验1： 某汽车制造商用三次样条插值设计车门的曲线，其中一段的数据如下： i x 0 1 2 3 4 5 6 7 8 9 10 i y 0.0 0.79 1.53 2.19 2.71 3.03 3.27 2.89 3.06 3.19 3.29 i y 0.8 0.2 用三次样条插值求)(10x S ，用软件绘制)(10x S 的图像，即车门的曲线。 1、计算)(10x S ： 程序： // splineaaaa.cpp : 定义控制台应用程序的入口点。 //**三次样条差值** //***第一步，利用差商，代替导数，求差商； //***第二步，利用追赶法求解三对角方程组，得到M[i]; //***第三步，将求得值带入三次样条函数，求得S(x); #include 主程序M文件（King.m）： function varargout = King(varargin) % KING MATLAB code for King.fig % KING, by itself, creates a new KING or raises the existing % singleton*. % % H = KING returns the handle to a new KING or the handle to % the existing singleton*. % % KING('CALLBACK',hObject,eventData,handles,...) calls the local % function named CALLBACK in KING.M with the given input arguments. % % KING('Property','Value',...) creates a new KING or raises the % existing singleton*. Starting from the left, property value pairs are % applied to the GUI before King_OpeningFcn gets called. An % unrecognized property name or invalid value makes property application % stop. All inputs are passed to King_OpeningFcn via varargin. % % *See GUI Options on GUIDE's Tools menu. Choose "GUI allows only one % instance to run (singleton)". % % See also: GUIDE, GUIDATA, GUIHANDLES % Edit the above text to modify the response to help King % Last Modified by GUIDE v2.5 06-May-2013 21:39:33 % Begin initialization code - DO NOT EDIT gui_Singleton = 1; gui_State = struct('gui_Name', mfilename, ... 'gui_Singleton', gui_Singleton, ... 'gui_OpeningFcn', @King_OpeningFcn, ... 'gui_OutputFcn', @King_OutputFcn, ... 'gui_LayoutFcn', [] , ... 'gui_Callback', []); 样条插值是一种工业设计中常用的、得到平滑曲线的一种插值方法，三次样条又是其中用的较为广泛的一种。 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, …, n b. 每一分段都是三次多项式函数曲线 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。则可以推出 2015级《数值分析》课外课堂大作业 论文题目：基于多项式插值与三次样条插值曲线拟合的比较姓名：XXX 学号：XXXXXXXXXXX 学院：XXXXXXXXXXXXXXX 专业方向: XXXXXXXXXXXXXXX 联系方式：（QQ号） （手机号） 导师姓名： 完成人（亲笔）签字 二0一五年十二月 基于多项式插值与三次样条插值曲线拟合的比较 摘要：在数值计算中经常要计算函数，当函数只在有限点集上给定函数值要包含改点集的区间上用公式给出函数的简单表达式，这就涉及在已知区间上用简单函数逼近已知复杂函数问题。本文为了解决这类问题就采用多项式插值与三次样条插值两种插值法并利用MATLAB数值分析软件进行编程，实现相应数据的曲线拟合以获得最佳曲线模型与相应数据的曲线拟合，选出最优的插值法以解决所给数据的曲线拟合问题。 关键词：函数；多项式插值；三次样条插值；曲线拟合；MATLAB Abstract:In numerical analysis ,the function value is often calculated .when the function is only given a function point set ,the simple expression of the function is given by the interval .which involves the use of a simple function to approximate the known complex function .in order to solve this problem ,we use polynomial interpolation and cubic spline interpolation tow kind of interpolation method and use MATLAB numerical analysis software to program ,to achieve the curve fitting of the corresponding date to obtain the best cure fitting ,and to choose the best interpolation method to solve the problem of curve fitting to the date. Keyword: Function ; Polynomial interpolation ; Cubic spline interpolation ; Fitting of a curve ; MATLAB 三次B 样条曲线插补算法的VC 实现 孔凡国,郝尚华,钟廷志 （五邑大学 机电系，广东 江门 529020） 摘要：结合三次B 样条曲线的数学性质，理论分析了三次B 样条曲线的插补过程，推导出了三次B 样条曲线的插补运算公式。在对三次B 样条曲线进行研究的基础上，介绍三次B 样条曲线插补算法实现方法，并在VC 环境中对B 样条曲线的插补过程进行了实现，实现了数控插补的动态运算。 关键词：三次B 样条；插补；VC VC Implement of the Cubic B-spline curve Interpolation KONG Fan-guo ，HAO Shang-hua ，ZHONG Ting-zhi （Department of Mechatronics ，Wuyi University ，Jiangmen 529020，China ） Abstract ：Combined with mathematics characters of the cubic B-spline curve, the interpolation process of the cubic B-spline curve is analyzed and calculation interpolation formula is deduced. Finally the cubic B-spline curve interpolation operator is programmed under the platform of VC to realize dynamic simulation. Keyword ：cubic B-spline ；interpolation ；VC 在数控加工中用一小段直线或圆弧去拟合实际曲线，这种拟合方法就是“插补”。 它实质上是根据有限的信息完成“数据密化”的工作。插补的计算方法和计算精度影响到整个数控系统的精度和效率，因此插补算法对整个数控系统的性能指标至关重要，可以说插补是整个数控系统控制软件的核心。 1 三次B 样条曲线插补原理 B 样条曲线是对Bezier 曲线的改进，它不仅保留了Bezier 曲线的优点，而且具有局部控制的能力，B 样条曲线方程可为： 给定n+1个控制点i P （i =0，1，……，n ），也称为特征多边形的顶点，k 次（k ＋1阶）B 样条曲线的表达式是： ,0()()n i i k i P u N u ==∑P 1k n ≤≤ （1） 在上式中当k=3，i =0，1，2，3时，可得三次B 样条曲线方程是： 3,30()()i i i P u N u ==∑P 01u ≤≤ （2） 用矩阵形式可表示为： 013223133136301()1 01303061410P u u u u u --????????-??????=≤≤??????-???????? P P P P （3） 第三章样条插值与曲线拟合 学习目标：掌握分段线性插值、分段Hermite插值、样条插值 方法以及贝齐尔曲线拟合曲 线的方法。重点是分段线性 插值、分段Hermite插值、样 条插值。 1901年龙格（Runge ）给出一个例子： ，定义在区间[-1,1]上，这是一个很光滑的函数，它的任意阶导数都存在，对它在[-1,1]上作等距节点插值时，插值多项式的情况见图1 §1 多项式插值的龙格现象 22511)(x x f += 俄罗斯数学家伯恩斯坦（C.H.Bernstein ）在1916年还给出如下定理。 定理1 函数 在[-1,1]上取n 个等距节点 ，构造n-1次插值多项式 ，当n 增大时，除了-1,0,1三点外，在[-1,1]中任何点处都不收敛于 。 x x f = )(1,11=-=n x x )(1 x P n -x 上述介绍的现象和定理告诉我们用高次插值多项式是不妥当的，从数值计算上来看也是这样，前一章介绍过的差分的误差传播会随阶数的提高越来越严重，因此，实践上作插值时一般只用一次、二次，最多用三次插值多项式。而提高插值精度的方法，可采用分段插值： 譬如在[a,b]上定义的连续函数 ，在[a,b]上取节点 构造一个分段一次多项式 ，即 在 上为 由一次插值的余项知在 上， )(x f b x x x x a n n =<<<<=-121 )(x L )(x L ],[1+i i x x i i i i i i i i x x x x x f x x x x x f --+--++++1111)()() )()(,,()()()(11++--=-=i i i i x x x x x x x f x L x f x R ],[1+i i x x 228 )(h M x R ≤ 3．6三次样条插值 一、教学目标及基本要求 通过对本节的学习，使学生掌握三次样条插值方法。 二、教学内容及学时分配 本章主要介绍线性方程求根的迭代法的加速方法。要求 1．了解数值分析的三次样条函数及有关概念。 2．正确理解三次样条差值的基本思想、数学原理、算法设计。 3．了解插值是数值逼近的重要方法之一，正确理解三次样条插值的基本思想、计算公式、算法设计、程序框图设计和源程序。 4．掌握三次样条差值原理和程序设计方法。 三、教学重点难点 1．教学重点：三次样条函数、三次样条插值。 2. 教学难点：三次样条插值。 四、教学中应注意的问题 多媒体课堂教学为主。适当提问，加深学生对概念的理解，迭代加速的算法实现。 五、教案正文 一 样条函数的概念 分段线性插值在节点处没有连续的一阶导函数，其光滑性较差。对于飞机的机翼的型线及船舶型往往要求有二阶光滑度（即在节点处要求二阶导函数连续）。 样条函数的概念来源于工程设计的实践。所谓“样条”（spline）是早期工程设计中的一种绘图工具，它是富有弹性的细长条。绘图时，用压铁迫使样条通过指定的型值点，并保证样条的光滑外形。在绕度不大的情况下，样条的曲线即为三次样条函数。 二 几何意义 三 构造三次样条函数的理论分析 如上图所示，通过已知的六个点，构造5个三次多项式函数分别是：红色、蓝色、黑色、紫色和绿色5根曲线。为确定一根曲线，就需要确定4个待定系数，所以总共需要4*5=20个待定系数。 另外，分析需要的约束条件。 每一根函数都要过已知的左右两个点，则有5*2=10个约束条件。 此外，每两个相邻曲线在相邻点处要求充分光滑，即在连接点处左右两个函数在该点具有1次和2次的导函数连续，图中有4个“中间点”，故又有4*2=8个约束条件。 若在整个图形的两端在加2个约束条件，整个3次样条函数就确定了。如： ①左右两端点上的1阶导函数已知； ②左右两端点上的2阶导函数已知，如()()00n S x S x ′′′′==（称为自然边界条件）； ③若原来的函数f(x)是以xn-x0为周期的周期函数，则y0=yn，且()()000n S x S x ′′′′+=?。 数值计算方法作业 实验4.3 三次样条差值函数 实验目的： 掌握三次样条插值函数的三弯矩方法。 实验函数: dt e x f x t ?∞ -- = 2 2 21 )(π x 0.0 0.1 0.2 0.3 0.4 F(x) 0.5000 0.5398 0.5793 0.6179 0.7554 求f(0.13)和f(0.36)的近似值 实验内容： （1） 编程实现求三次样条插值函数的算法，分别考虑不同的边界条件； （2） 计算各插值节点的弯矩值； （3） 在同一坐标系中绘制函数f(x)，插值多项式，三次样条插值多项式的曲线 比较插值结果。 实验4.5 三次样条差值函数的收敛性 实验目的： 多项式插值不一定是收敛的，即插值的节点多，效果不一定好。对三次样条插值函数如何呢？理论上证明三次样条插值函数的收敛性是比较困难的，通过本实验可以证明这一理论结果。 实验内容： 按照一定的规则分别选择等距或非等距的插值节点，并不断增加插值节点的个数。 实验要求： （1） 随着节点个数的增加，比较被逼近函数和三样条插值函数的误差变化情 况，分析所得结果并与拉格朗日插值多项式比较； （2） 三次样条插值函数的思想最早产生于工业部门。作为工业应用的例子，考 实验名称 实验 4.3三次样条插值函数（P126） 4.5三次样条插值函数的收敛性（P127） 实验时间 姓名 班级 学号 成绩 虑如下例子：某汽车制造商根据三次样条插值函数设计车门曲线，其中一 段数据如下： k x 0 1 2 3 4 5 6 7 8 9 10 k y 0.0 0.79 1.53 2.19 2.71 3.03 3.27 2.89 3.06 3.19 3.29 k y ' 0.8 0.2 算法描述： 拉格朗日插值： 错误！未找到引用源。 其中错误！未找到引用源。是拉格朗日基函数，其表达式为：() ∏ ≠=--=n i j j j i j i x x x x x l 0) ()( 牛顿插值： ) )...()(](,...,,[....))(0](,,[)0](,[)()(1102101210100----++--+-+=n n n x x x x x x x x x x f x x x x x x x f x x x x f x f x N 其中????? ?? ?? ?????? --=--= --= -)/(]),...,[],...,[(]...,[..],[],[],,[)()(],[01102110x x x x x f x x x f x x x f x x x x f x x f x x x f x x x f x f x x f n n n n i k j i k j k j i j i j i j i 三样条插值： 所谓三次样条插值多项式Sn(x)是一种分段函数，它在节点Xi(a 已知一组数据点，编写一程序求解三次样条插值函数满足 并针对下面一组具体实验数据 求解，其中边界条件为. 解：Matlab计算程序为： clear clc x=[0.25 0.3 0.39 0.45 0.53] y=[0.5000 0.5477 0.6245 0.6708 0.7280] n=length(x); for i=1:n-1 h(i)=x(i+1)-x(i); end for i=1:n-2 k(i)=h(i+1)/(h(i)+h(i+1)); u(i)=h(i)/(h(i)+h(i+1)); end for i=1:n-2 gl(i)=3*(u(i)*(y(i+2)-y(i+1))/h(i+1)+k(i)*(y(i+1)-y(i))/h(i)); end g0=3*(y(2)-y(1))/h(1); g00=3*(y(n)-y(n-1))/h(n-1); g=[g0 gl g00]; g=transpose(g) k1=[k 1]; u1=[1 u]; Q=2*eye(5)+diag(u1,1)+diag(k1,-1) m=transpose(Q\g) syms X; for i=1:n-1 p1(i)=(1+2*(X-x(i))/h(i))*((X-x(i+1))/h(i))^2*y(i); p2(i)=(1-2*(X-x(i+1))/h(i))*((X-x(i))/h(i))^2*y(i+1); p3(i)=(X-x(i))*((X-x(i+1))/h(i))^2*m(i); p4(i)=(X-x(i+1))*((X-x(i))/h(i))^2*m(i+1); a = [1.0000 0.5000 0.2500 0.0500 0.0100 0.0050 0.0010]; b = [100.0000 97.8000 94.5000 79.0000 32.8000 23.0000 13.2000]; xx= [ 1.0000 0.0500 0.0020]; yy=interp1(a,b,xx,'pchip'); %plot绘图 plot(a,b,'o',a,b); hold on xxx=1:-0.01:0.001; yspline=interp1(a,b,xxx,'spline'); ypchip=interp1(a,b,xxx,'pchip'); ycubic=interp1(a,b,xxx,'cubic'); plot(xxx,yspline, '--bo'); hold on plot(xxx,ypchip,'-rs'); hold on plot(xxx,ycubic, '-kx'); grid on xlabel('土壤粒径（毫米）'); ylabel('颗粒累积百分数（%）'); title('土壤颗粒级配曲线') hold off 三种三次样条插值函数比较 批量转换程序： a = [2 1.0000 0.5000 0.2500 0.0500 0.0100 0.0050 0.0010]; b = [100.4 100 99 97.6 87.2 36.2 24.8 10 100.5 100 99.1 97.9 88 38 25 10.2 100.6 100 99.4 98.6 90 41.3 28.5 12.5 100.5 100 99.4 98.7 89.6 41 28.9 14.1三次样条插值作业题
三次样条插值方法的应用
样条插值和曲线拟合
三次样条插值在工程拟合中的应用
三次样条拟合范例
三次样条插值的MATLAB实现
关于三次样条插值函数的学习报告(研究生)资料
试求三次样条插值S(X)
三次样条函数
三次样条差值拟合车门曲线
三次样条曲线插补模拟
样条函数(三次样条)
基于多项式插值与三次样条插值曲线拟合的比较
三次B样条曲线插补算法的VC实现
第三章 样条插值和曲线拟合
三次样条插值
数值分析作业-三次样条插值..
Matlab程序三次样条插值函数
三次样条函数程序