数值分析实验报告-插值、三次样条Word版
- 格式:doc
- 大小:359.50 KB
- 文档页数:6
实验报告:牛顿差值多项式&三次样条问题:在区间[-1,1]上分别取n=10、20用两组等距节点对龙格函数21()25f x x作多项式插值及三次样条插值,对每个n 值,分别画出插值函数及()f x 的图形。
实验目的:通过编程实现牛顿插值方法和三次样条方法,加深对多项式插值的理解。
应用所编程序解决实际算例。
实验要求:1. 认真分析问题,深刻理解相关理论知识并能熟练应用;2. 编写相关程序并进行实验;3. 调试程序,得到最终结果;4. 分析解释实验结果;5. 按照要求完成实验报告。
实验原理:详见《数值分析 第5版》第二章相关内容。
实验内容:(1)牛顿插值多项式1.1 当n=10时:在Matlab 下编写代码完成计算和画图。
结果如下:代码:clear allclcx1=-1:0.2:1;y1=1./(1+25.*x1.^2);n=length(x1);f=y1(:);for j=2:nfor i=n:-1:jf(i)=(f(i)-f(i-1))/(x1(i)-x1(i-j+1));endendsyms F x p ;F(1)=1;p(1)=y1(1);for i=2:nF(i)=F(i-1)*(x-x1(i-1));p(i)=f(i)*F(i);endsyms PP=sum(p);P10=vpa(expand(P),5);x0=-1:0.001:1;y0=subs(P,x,x0);y2=subs(1/(1+25*x^2),x,x0);plot(x0,y0,x0,y2)grid onxlabel('x')ylabel('y')P10即我们所求的牛顿插值多项式,其结果为:P10(x)=-220.94*x^10+494.91*x^8-9.5065e-14*x^7-381.43*x^6-8.504e-14*x^5+123.36*x^4+2.0202e-1 4*x^3-16.855*x^2-6.6594e-16*x+1.0并且这里也能得到该牛顿插值多项式的在[-1,1]上的图形,并和原函数进行对比(见Fig.1)。
《数值分析课程设计-三次样条插值》报告掌握三次样条插值函数的构造方法,体会三次样条插值函数对被逼近函数的近似。
三次样条插值函数边界条件由实际问题对三次样条插值在端点的状态要求给出。
以第1 边界条件为例,用节点处二阶导数表示三次样条插值函数,用追赶法求解相关方程组。
通过Matlab 编制三次样条函数的通用程序,可直接显示各区间段三次样条函数体表达式,计算出已给点插值并显示各区间分段曲线图。
引言分段低次样条插值虽然计算简单、稳定性好、收敛性有保证且易在电子计算机上实现,但只能保证各小段曲线在连接处的连续性,不能保证整件曲线的光滑性。
利用样条插值,既可保持分段低次插值多项式,又可提高插值函数光滑性。
故给出分段三次样条插值的构造过程算法步骤,利用Matlab软件编写三次样条插值函数通用程序,并通过数值算例证明程序的正确性。
三次样条函数的定义及特征定义:设[a,b] 上有插值节点,a=x1<x2<…xn=b,对应函数值为y1,y2,⋯yn。
若函数S(x) 满足S(xj) = yj ( j = 1,2, ⋯,n ), S(x) 在[xj,xj+1] ( j =1,2,⋯,n-1)上都是不高于三的多项式(为了与其对应j 从1 开始,在Matlab 中元素脚标从1 开始)。
当S(x) 在 [a,b] 具有二阶连续导数。
则称S(x) 为三次样条插值函数。
要求S(x) 只需在每个子区间[xj,xj+1] 上确定 1 个三次多项式,设为:Sj(x)=ajx3+bjx2+cjx+dj, (j=1,2,⋯,n-1) (1)其中aj,bj,cj,dj 待定,并要使它满足:S(xj)=yj, S(xj-0)=S(xj+0), (j=2,⋯,n-1) (2)S'(xj-0)=S'(xj+0), S"(xj-0)=S"(xj+0), (j=2,⋯,n-1) (3)式(2)、(3)共给出n+3(n-2)=4n-6 个条件,需要待定4(n-1) 个系数,因此要唯一确定三次插值函数,还要附加2个边界条件。
数值分析课程设计报告书院系名称:学生姓名:专业名称:班级:时间:实验一 三次样条插值的三弯矩法一、实验目的已知数据i x ,()i i y f x =,0,,i n =及边界条件()n j x y j j 1,0),(2=,求)(x f 的三次样条插值函数)(x S .要求输出用追赶法解出的弯矩向量0[,,]n M M M =及()(),0,,,0,1,2k i S t i m k ==的值.画出)(x S y =的图形,图形中描出插值点(,)i i x y 及(,())i i t S t 分别用‘o ’和‘*’标记.二、实验原理1.用追赶法求解第二类边界条件的三弯矩方程:0010012111121111[,,]21[,,]26[,,]212[,,]n n n n n n n n n n f x x x M f x x x M M f x x x M f x x x μλμλ------⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦ 其中1111,,j jj j j j j j j j j h h h x x h h h h μλ-+--===-++.2.得出样条函数表达式:332211111()()()()()6666j j j j j j j j j j j j j j j jx x x x M h x x M h x x S x M M y y h h h h +++++----=++-+-. 3.计算(k)(),0,,,0,1,2i S t i m k ==.三、实验结果所用数据:x=[-2.223,-1.987,-1.8465,-1.292,-1.2266,-1.1056,-0.8662,-0.6594,-0.2671,-0.0452,0.5385,1.2564,1.4398,1.5415,1.7646,1.9678,2.236];y=[0.83995,1.1696,1.3141,1.6992,1.7312,1.7847,1.8708,1.9262,1.9881,1.9997,1.9511,1.7169,1.618,1.5543,1.3871,1.191,0.81662];d2s1= -4.5000;d2sn= -4.8967; %第二种边界条件t=[-2.223,-1.9443,-1.6656,-1.3869,-1.1083,-0.82956,-0.55088,-0.27219,0.0065,0.28519,0.56387,0.84256,1.1212,1.3999,1.6786,2.236]; ;(指定计算点)计算结果:-2.5-2-1.5-1-0.500.51 1.52 2.50.811.21.41.61.82四、实验分析通过实验结果我们,知道三弯矩法求出满足初始条件的三次样条函数,与其他插值函数的构造相比,三次样条插值法的计算量要小得多。
数值分析实验报告
一、实验背景
本实验主要介绍了数值分析的各种方法。
在科学计算中,为了求解一
组常微分方程或一些极限问题,数值分析是一种有用的方法。
数值分析是
一种运用计算机技术对复杂模型的问题进行数学分析的重要手段,它利用
数学模型和计算机程序来解决复杂的数学和科学问题。
二、实验内容
本实验通过MATLAB软件,展示了以下几种数值分析方法:
(1)拉格朗日插值法:拉格朗日插值法是由法国数学家拉格朗日发
明的一种插值方法,它可以用来插值一组数据,我们使用拉格朗日插值法
对给定的点进行插值,得到相应的拉格朗日多项式,从而计算出任意一个
点的函数值。
(2)最小二乘法:最小二乘法是一种常用的数据拟合方法,它可以
用来拟合满足一定函数的点的数据,它的主要思想是使得数据点到拟合曲
线之间的距离的平方和最小。
(3)牛顿插值法:牛顿插值法是一种基于差商的插值方法,它可以
用来插值一组数据,可以求得一组数据的插值函数。
(4)三次样条插值:三次样条插值是一种基于三次样条的插值方法,它可以用来对一组数据进行插值,可以求得一组数据的插值函数。
三、实验步骤
1.首先启动MATLAB软件。
第1篇一、实验目的1. 理解并掌握插值法的基本原理和常用方法。
2. 学习使用拉格朗日插值法、牛顿插值法等数值插值方法进行函数逼近。
3. 分析不同插值方法的优缺点,并比较其精度和效率。
4. 通过实验加深对数值分析理论的理解和应用。
二、实验原理插值法是一种通过已知数据点来构造近似函数的方法。
它广泛应用于科学计算、工程设计和数据分析等领域。
常用的插值方法包括拉格朗日插值法、牛顿插值法、样条插值法等。
1. 拉格朗日插值法拉格朗日插值法是一种基于多项式的插值方法。
其基本思想是:给定一组数据点,构造一个次数不超过n的多项式,使得该多项式在这些数据点上的函数值与已知数据点的函数值相等。
2. 牛顿插值法牛顿插值法是一种基于插值多项式的差商的插值方法。
其基本思想是:给定一组数据点,构造一个次数不超过n的多项式,使得该多项式在这些数据点上的函数值与已知数据点的函数值相等,并且满足一定的差商条件。
三、实验内容1. 拉格朗日插值法(1)给定一组数据点,如:$$\begin{align}x_0 &= 0, & y_0 &= 1, \\x_1 &= 1, & y_1 &= 4, \\x_2 &= 2, & y_2 &= 9, \\x_3 &= 3, & y_3 &= 16.\end{align}$$(2)根据拉格朗日插值公式,构造插值多项式:$$P(x) = \frac{(x-x_1)(x-x_2)(x-x_3)}{(x_0-x_1)(x_0-x_2)(x_0-x_3)}y_0 + \frac{(x-x_0)(x-x_2)(x-x_3)}{(x_1-x_0)(x_1-x_2)(x_1-x_3)}y_1 + \frac{(x-x_0)(x-x_1)(x-x_3)}{(x_2-x_0)(x_2-x_1)(x_2-x_3)}y_2 + \frac{(x-x_0)(x-x_1)(x-x_2)}{(x_3-x_0)(x_3-x_1)(x_3-x_2)}y_3.$$(3)计算插值多项式在不同点的函数值,并与实际值进行比较。
武汉理工大学学生实验报告书实验课程名称数值分析开课学院计算机科学与技术学院指导老师姓名学生姓名学生专业班级2010—2010学年第一学期实验课程名称:数值分析第二部分:实验调试与结果分析(可加页)一、调试过程(包括调试方法描述、实验数据记录,实验现象记录,实验过程发现的问题等)(1)用拉格朗日插值法计算时,输入及运行结果如下:拉格朗日插值法牛顿插值法(2)利用二次插值计算时,输入及运行结果如下:拉格朗日插值法牛顿插值法(3)用艾尔米特插值法计算时,f(x)的插值多项式H5(x)=(1+4*x)*(x-0.5)*(x-0.5)*(x-2)*(x-2)+(3.90807-6.03838*x)*(x-2)*(x-2)*x*x+(2.34573-4.16674*x)*x*x*(x-0.5)*(x-0.5)(4)各插值算法的精度差异比较经过比较,拉格朗日插值法要比牛顿插值法算法的计算量多一些,拉格朗日插值法后一次计算时用到了前一次计算的结果,提高了运算的效率,但拉格朗日插值法在构造艾尔米特插值法时很方便,将坐标点和对应的导数结合起来的精度比线性插值的精度又要高一些。
但从实验数据来看,在坐标不是很多的情况下,已知的点越多精度也就相对较高。
对于实验要求的第二组数据用拉格朗日插值法(或者牛顿插值法)实验结果如下:一下分别是二阶、三阶、四阶、五阶插值得到的结果以上只是实验结果的一部分,改变插值的位置时,得到的实验结果精度也是有所不同的。
由以上结果分析可知,插值次数并不是越多越好,多了反而会让结果更加偏离真实结果,这充分说明了高次插值存在“病态性质”,在已知点很多的情况下应该采用分段低次插值,将拉格朗日插值法和牛顿插值法运用到分段低次插值法当中,这样得到的结果可能胡更加精确。
数值分析插值实验报告引言插值是数值分析中常用的一种技术,通过已知点的函数值来推测未知点的函数值。
在实际应用中,我们经常需要根据有限的数据点来估计连续函数的值,这时插值就起到了关键作用。
本实验旨在通过插值方法来推测未知数据点的函数值,并对比不同插值方法的精度和效果。
实验目的1.了解插值的基本概念和方法;2.掌握常见的插值方法,如拉格朗日插值、牛顿插值等;3.对比不同插值方法的精度和效果,分析其优缺点。
实验步骤1.数据采集:选取一组已知数据点,作为插值的基础。
这些数据点可以是从实际场景中测量得到的,也可以是人为设定的。
2.插值方法选择:根据实验要求和数据特点,选择适合的插值方法。
常见的插值方法包括拉格朗日插值、牛顿插值、分段线性插值等。
3.插值计算:根据选定的插值方法,利用已知数据点进行计算,并得到插值结果。
4.结果分析:比较插值结果与实际数据的差异,并评估插值方法的精度和效果。
可以使用误差分析等方法进行评估。
5.优化调整:根据实验结果和需求,对插值方法进行优化调整,以提高插值的准确性和可靠性。
实验结果与讨论通过实验,我们得到了不同插值方法的结果,并进行了对比和分析。
根据实验数据和误差分析,我们可以得出以下结论:1.拉格朗日插值方法具有较高的插值精度,在一定程度上能够准确地模拟实际数据。
2.牛顿插值方法相对于拉格朗日插值方法而言,对于大量数据点的计算速度更快,但在少量数据点的情况下,两者的精度差异较小。
3.分段线性插值方法适用于数据点较为离散的情况,能够提供较为平滑的插值结果。
4.插值方法的选择应根据具体需求和数据特点进行,没有一种插值方法适用于所有情况。
实验总结通过本次实验,我们对插值方法有了更深入的了解,并掌握了常见的插值方法的原理和应用。
实验结果表明,插值方法在数值分析中起到了重要的作用,能够准确地推测未知点的函数值。
然而,在实际应用中,我们还需要考虑数据的特点、插值方法的适用性以及计算效率等因素。
关于三次样条插值函数的学习报告三次样条插值函数是一种广泛应用于数值分析领域的插值方法,用于逼近一组已知数据点构成的函数。
在这篇学习报告中,我将介绍三次样条插值函数的定义、原理、应用及其优缺点,并通过实际例子说明其如何在实际问题中使用。
一、三次样条插值函数的定义三次样条插值函数是指用分段三次多项式对一组已知数据点进行插值的方法。
具体来说,对于已知数据点$(x_0,y_0),(x_1,y_1),...,(x_n,y_n)$,三次样条插值函数会在每相邻两个数据点之间构造一个三次多项式,使得这些多项式在相应的数据点上满足插值条件,并且在相邻两个多项式之间满足一定的连续性条件。
二、三次样条插值函数的原理三次样条插值函数的原理是利用三次多项式在每个数据点上的取值和导数值来确定三次多项式的系数,从而构造出满足插值条件和连续性条件的插值函数。
具体来说,对于每个相邻的数据点$(x_i,y_i),(x_{i+1},y_{i+1})$,我们可以构造一个三次多项式$S_i(x)$,满足以下条件:1.$S_i(x_i)=y_i$,$S_i(x_{i+1})=y_{i+1}$,即在数据点上满足插值条件;2.$S_i'(x_{i+1})=S_{i+1}'(x_{i+1})$,$S_i''(x_{i+1})=S_{i+1}''(x_{i+1})$,即在数据点上满足连续性条件。
通过求解上述条件,可以得到每个相邻数据点之间的三次多项式$S_i(x)$,从而得到整个插值函数。
三、三次样条插值函数的应用三次样条插值函数在数值分析领域有广泛的应用,尤其在曲线拟合、数据逼近等问题中起到重要作用。
例如,当我们需要根据已知的离散数据点绘制平滑的曲线图形时,可以使用三次样条插值函数来进行插值,从而得到更加连续和光滑的曲线。
另外,在信号处理、图像处理等领域也常常会用到三次样条插值函数。
例如,在数字图像处理中,我们需要对像素点进行插值以得到更高分辨率的图像,三次样条插值函数可以很好地满足这个需求,使图像更加清晰和真实。
南京理工大学课程考核论文课程名称:高等数值分析论文题目:基于matlab的函数插值方法性能比较姓名:xxx学号:xxxxxxxxxx成绩:摘要函数插值是指已知某函数的在若干离散点上的函数值或者导数信息,通过求解该函数中待定形式的插值函数以及待定系数,使得该函数在给定离散点上满足约束。
本文首先介绍了五种插值方法:线性插值、lagrange插值、newdun插值、三次样条插值和三次B样条插值,并从五种插值法的基本思想和具体实例仿真入手,探讨了五种插值法的优缺点。
通过对五种插值法的对比研究及实际应用的总结,从而使我们在以后的应用中能够更好、更快的解决问题。
关键字插值法对比matlab目录摘要 (2)0 引言 (4)1插值问题的提出、发展史及简单应用 (4)1.1插值问题的提出 (4)1.2插值法的发展史 (4)1.3插值法的简单应用 (4)2 五种插值法的定义 (5)2.1线性插值法 (5)2.2Lagrange插值法 (5)2.3Newton插值法 (6)2.4 三次样条插值法 (6)2.5B样条插值 (6)3五种插值法的matlab仿真实现 (8)4五种插值方法性能对比 (11)5结束语 (12)参考文献 (12)0 引言近半世纪由于计算机的广泛使用和造船、航空、精密机械加工等世纪问题的需要,使插值法在理论上和实践上得到进一步发展,尤其是20世纪40年代后期发展起来的样条插值等,更获得广泛应用,称为计算机图形学的基础。
插值是指已知某函数的在若干离散点上的函数值或者导数信息,通过求解该函数中待定形式的插值函数以及待定系数,使得该函数在给定离散点上满足约束。
插值函数又叫作基函数,如果该基函数定义在整个定义域上,叫作全域基,否则叫作分域基。
1插值问题的提出、发展史及简单应用1.1插值问题的提出许多实际问题都用函数来表示某种内在规律的数量关系,其中相当一部分函数是通过实验或观测得到的。
虽然()x f 在某个区间[]b a ,上是存在的,有的还是连续的,但却只能给出[]b a ,上一系列点i x 的函数值()() 2,1,0==i x f y i i ,这只是一张函数表.有的函数虽有解析表达式,但由于计算复杂,使用不方便,通常也造一个函数表,如大家熟悉的三角函数表、对数表、平方根和立方根表.为了研究函数的变化规律,往往需要求出不在表中的函数值.因此,我们希望根据给定的函数表做一个既能反映函数()x f 的特性,又便于计算简单函数()x p ,用()x p 近似()x f 。
实验报告:牛顿差值多项式&三次样条问题:在区间[-1,1]上分别取n=10、20用两组等距节点对龙格函数21()25f x x作多项式插值及三次样条插值,对每个n 值,分别画出插值函数及()f x 的图形。
实验目的:通过编程实现牛顿插值方法和三次样条方法,加深对多项式插值的理解。
应用所编程序解决实际算例。
实验要求:1. 认真分析问题,深刻理解相关理论知识并能熟练应用; 2. 编写相关程序并进行实验; 3. 调试程序,得到最终结果; 4. 分析解释实验结果; 5. 按照要求完成实验报告。
实验原理:详见《数值分析 第5版》第二章相关内容。
实验内容:(1)牛顿插值多项式1.1 当n=10时:在Matlab 下编写代码完成计算和画图。
结果如下:代码:clear allclcx1=-1:0.2:1;y1=1./(1+25.*x1.^2);n=length(x1);f=y1(:);for j=2:nfor i=n:-1:jf(i)=(f(i)-f(i-1))/(x1(i)-x1(i-j+1));endendsyms F x p;F(1)=1;p(1)=y1(1);for i=2:nF(i)=F(i-1)*(x-x1(i-1));p(i)=f(i)*F(i);endsyms PP=sum(p);P10=vpa(expand(P),5);x0=-1:0.001:1;y0=subs(P,x,x0);y2=subs(1/(1+25*x^2),x,x0);plot(x0,y0,x0,y2)grid onxlabel('x')ylabel('y')P10即我们所求的牛顿插值多项式,其结果为:P10(x)=-220.94*x^10+494.91*x^8-9.5065e-14*x^7-381.43*x^6-8.504e-14*x^5+123.36* x^4+2.0202e-14*x^3-16.855*x^2-6.6594e-16*x+1.0并且这里也能得到该牛顿插值多项式的在[-1,1]上的图形,并和原函数进行对比(见Fig.1)。
数值计算实验—实验报告2一、实验项目:第二类边界条件三次样条插值二、实验目的和要求a.通过本实验深入地理解三次样条插值多项式的基本原理b.通过数值算例更好的领会三次样条插值多项式具有较高的准确性三、实验内容1.用调试好的程序解决如下问题:点中点处的函数值,并将计算结果与sinx在相应点的数值相比较。
n=8;p1=0.4794;pn=0.9463;u=[0.6,0.8,1.0,1.2,1.4,1.6,1.8];p=7;x=[0.5,0.7,0.9,1.1,1.3,1.5,1.7,1.9];y=[0.4794,0.6442,0.7833,0.8912,0.9636,0.9975,0.9917,0.9463];for i=1:n-1h(i)=x(i+1)-x(i);enda2(1)=1;g(1)=3*(y(2)-y(1))/h(1)-p1*h(1)/2;for k=2:n-1a1(k-1)=h(k)/(h(k)+h(k-1));a2(k)=h(k-1)/(h(k)+h(k-1));g(k)=3*a2(k)*(y(k+1)-y(k))/h(k)+3*a1(k-1)*(y(k)-y(k-1))/h(k-1); enda1(n-1)=1;g(n)=3*(y(n)-y(n-1))/h(n-1)+pn * h(n-1)/2;%追赶法求三转角方程b1(1)=2;m(1)=g(1)/2;b2(1)=a2(1)/b1(1);for i=2:nb1(i)=2-a1(i-1)*b2(i-1);if(i~=n)b2(i)=a2(i)/b1(i);endm(i)=(g(i)-a1(i-1)*m(i-1))/b1(i);endfor i=n-1:-1:1m(i)=m(i)-b2(i)*m(i+1);endfor j=1:pfor i=1:nif((u(j)>=x(i))&&(u(j)<x(i+1)))k=i;break;endends(j)=0;s(j)=s(j)+(h(k)+2*(u(j)-x(k)))*(u(j)-x(k+1))^2*y(k)/(h(k))^3;s(j)=s(j)+(h(k)-2*(u(j)-x(k+1)))*(u(j)-x(k))^2*y(k+1)/(h(k))^3;s(j)=s(j)+(u(j)-x(k))*(u(j)-x(k+1))^2*m(k)/(h(k))^2;s(j)=s(j)+(u(j)-x(k+1))*(u(j)-x(k))^2*m(k+1)/(h(k))^2;end(2).运行结果3. 根据Lagrange插值多项式基本原理编制程序,并计算下面的数值算例:=-5+kh,其中h=10/n,n=10,20,40.给定函数f(x)=1/(1+x^2)(-5≤x≤5),取等距节点xk边界条件为S''(x0)=f''(x0),S''(x n)=f''(x n).用上述算法计算S10(x),S20(x), S40(x),并与函数f(x)以及10次Lagrange插值多项式L10(x)在给定点处的函数值进行比较。
数值计算方法作业实验4.3 三次样条差值函数实验目的:掌握三次样条插值函数的三弯矩方法。
实验函数:dt ex f xt ⎰∞--=2221)(πx 0.0 0.1 0.2 0.3 0.4 F(x) 0.5000 0.5398 0.57930.61790.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 ky ' 0.80.2算法描述:拉格朗日插值:错误!未找到引用源。
其中错误!未找到引用源。
是拉格朗日基函数,其表达式为:()∏≠=--=ni j j j i ji 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 ji j i j i三样条插值:所谓三次样条插值多项式Sn(x)是一种分段函数,它在节点Xi(a<X0<X1……<Xn<b)分成的每个小区间[x i-1,x i ]上是三次多项式,其在此区间上的表达式如下:],[),6()6(]6)([6)(6)()(111113131i i ii i i i i i i i i i i i i i i i i i x x x h yM h M h h y x M M h h y y h x x Mi h x x M x S -------∈-+-+---+-+-=式中Mi=)(i x S ''.因此,只要确定了Mi 的值,就确定了整个表达式,Mi 的计算方法如下:令⎪⎪⎪⎪⎩⎪⎪⎪⎪⎨⎧=---+=+=+=+--++++++],,[6)(6111111111i i i i i i i i i i i i i i i i i i i ix x x f h y y h y y h h d h h h h h h λμ则Mi 满足如下n-1个方程:1,...2,1,211-==+++-n i d M M M i i i i i i λμ 常用的边界条件有如下几类:(1) 给定区间两端点的斜率m 0,m n ,即n n n m y x S m y x S ='='='=')(,)(000 (2) 给定区间两端点的二阶导数M0,Mn,即n n n M y x S M y x S =''=''=''='')(,)(000 (3) 假设y=f(x)是以b-a 为周期的周期函数,则要求三次样条插值函数S (x )也为周期函数,对S (x )加上周期条件2,1,0),0()0()(0)(=-=+p x S x S n p p对于第一类边界条件有⎪⎪⎩⎪⎪⎨⎧--=+--=+--)(62)(6211001110n n n n n n i h y y mn h M M m h y y h M M对于第二类边界条件有⎩⎨⎧=+=+-n n n n d M M d M M 221100μλ其中n n n n nnn M u x x f m h d M m x x f h d )1(2]),[(6)1(2)],[(6100001010-+-=-+-=-μλλ那么解就可以为⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛----n n n n n n n d d d d d M M M M M 1210121011...2...............2............................1..2.1......0..2μλμλμλ 对于第三类边界条件,)0()0(,,000-=+==n n n x S x S M M y y ,由此推得0010012d M M M n =-++μλ,其中]),1[],[(6,,101010110n n nn n n x x f x x f h h d h h h h h h --+=+=+=μλ,那么解就可以为: ⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫ ⎝⎛-------1221012101221100...2.............2..............................2..,,.......,..22n n n n n n n d d d d d M M M M M n μλλμλμμλ 程序代码: 1拉格朗日插值函数Lang.mfunction f=lang(X,Y,xi) %X 为已知数据的横坐标 %Y 为已知数据的纵坐标 %xi 插值点处的横坐标%f 求得的拉格朗日插值多项式的值 n=length(X); f=0; for i=1:n l=1; for j=1:i-1l=l.*(xi-X(j))/(X(i)-X(j)); end ; for j=i+1:nl=l.*(xi-X(j))/(X(i)-X(j)); end ;%拉格朗日基函数 f=f+l*Y(i); endfprintf('%d\n',f) return2 牛顿插值函数newton.mfunction f=newton(X,Y,xi) %X 为已知数据的横坐标 %Y 为已知数据的纵坐标 %xi 插值点处的横坐标%f 求得的拉格朗日插值多项式的值 n=length(X);newt=[X',Y'];%计算差商表for j=2:nfor i=n:-1:1if i>=jY(i)=(Y(i)-Y(i-1))/(X(i)-X(i-j+1));else Y(i)=0;endendnewt=[newt,Y'];end%计算牛顿插值f=newt(1,2);for i=2:nz=1;for k=1:i-1z=(xi-X(k))*z;endf=f+newt(i-1,i)*z;endfprintf('%d\n',f)return3三次样条插值第一类边界条件Threch.mfunction S=Threch1(X,Y,dy0,dyn,xi)% X为已知数据的横坐标%Y为已知数据的纵坐标%xi插值点处的横坐标%S求得的三次样条插值函数的值%dy0左端点处的一阶导数% dyn右端点处的一阶导数n=length(X)-1;d=zeros(n+1,1);h=zeros(1,n-1);f1=zeros(1,n-1);f2=zeros(1,n-2);for i=1:n%求函数的一阶差商h(i)=X(i+1)-X(i);f1(i)=(Y(i+1)-Y(i))/h(i);endfor i=2:n%求函数的二阶差商f2(i)=(f1(i)-f1(i-1))/(X(i+1)-X(i-1));d(i)=6*f2(i);endd(1)=6*(f1(1)-dy0)/h(1);d(n+1)=6*(dyn-f1(n-1))/h(n-1);%¸赋初值A=zeros(n+1,n+1);B=zeros(1,n-1);C=zeros(1,n-1);for i=1:n-1B(i)=h(i)/(h(i)+h(i+1));C(i)=1-B(i);endA(1,2)=1;A(n+1,n)=1;for i=1:n+1A(i,i)=2;endfor i=2:nA(i,i-1)=B(i-1);A(i,i+1)=C(i-1);endM=A\d;syms x;for i=1:nSx(i)=collect(Y(i)+(f1(i)-(M(i)/3+M(i+1)/6)*h(i))*(x-X(i))...+M(i)/2*(x-X(i))^2+(M(i+1)-M(i))/(6*h(i))*(x-X(i))^3);digits(4);Sx(i)=vpa(Sx(i));%三样条插值函数表达式endfor i=1:ndisp('S(x)=');fprintf('%s (%d,%d)\n',char(Sx(i)),X(i),X(i+1));endfor i=1:nif xi>=X(i)&&xi<=X(i+1)S=Y(i)+(f1(i)-(M(i)/3+M(i+1)/6)*h(i))*(xi-X(i))+M(i)/2*(xi-X(i))^2+(M (i+1)-M(i))/(6*h(i))*(xi-X(i))^3;endenddisp('xi S');fprintf('%d,%d\n',xi,S);return4 三次样条插值第二类边界条件Threch2.mfunction [Sx]=Threch2(X,Y,d2y0,d2yn,xi)X为已知数据的横坐标%Y为已知数据的纵坐标%xi插值点处的横坐标%S求得的三次样条插值函数的值%d2y0左端点处的二阶导数% d2yn右端点处的二阶导数n=length(X)-1;d=zeros(n+1,1);h=zeros(1,n-1);f1=zeros(1,n-1);f2=zeros(1,n-2);for i=1:n%求一阶差商h(i)=X(i+1)-X(i);f1(i)=(Y(i+1)-Y(i))/h(i);endfor i=2:n%求二阶差商f2(i)=(f1(i)-f1(i-1))/(X(i+1)-X(i-1));d(i)=6*f2(i);endd(1)=2*d2y0;d(n+1)=2*d2yn;%赋初值A=zeros(n+1,n+1);B=zeros(1,n-1);C=zeros(1,n-1);for i=1:n-1B(i)=h(i)/(h(i)+h(i+1));C(i)=1-B(i);endA(1,2)=0;A(n+1,n)=0;for i=1:n+1A(i,i)=2;endfor i=2:nA(i,i-1)=B(i-1);A(i,i+1)=C(i-1);endM=A\d;syms x;for i=1:nSx(i)=collect(Y(i)+(f1(i)-(M(i)/3+M(i+1)/6)*h(i))*(x-X(i))... +M(i)/2*(x-X(i))^2+(M(i+1)-M(i))/(6*h(i))*(x-X(i))^3);digits(4);Sx(i)=vpa(Sx(i));endfor i=1:ndisp('S(x)=');fprintf('%s (%d,%d)\n',char(Sx(i)),X(i),X(i+1));endfor i=1:nif xi>=X(i)&&xi<=X(i+1)S(i)=Y(i)+(f1(i)-(M(i)/3+M(i+1)/6)*h(i))*(xi-X(i))+M(i)/2*(xi-X(i))^2 +(M(i+1)-M(i))/(6*h(i))*(xi-X(i))^3;endenddisp('xi S');fprintf('%d,%d\n',xi,S);return5插值节点处的插值结果main3.mclearclcX=[0.0,0.1,0.2,0.3,0.4];Y=[0.5000,0.5398,0.5793,0.6179,0.7554];xi=0.13;%xi=0.36;disp('xi=0.13');%disp('xi=0.36');disp('拉格朗日插值结果');lang(X,Y,xi);disp('牛顿插值结果');newton(X,Y,xi);disp('三次样条第一类边界条件插值结果');Threch1(X,Y,0.40,0.36,xi);%0.4,0.36分别为两端点处的一阶导数disp('三次样条第二类边界条件插值结果');Threch2(X,Y,0,-0.136,xi);%0,-0.136分别为两端点处的二阶导数6将多种插值函数即原函数图像画在同一张图上main2.mclearclcX=[0.0,0.1,0.2,0.3,0.4];Y=[0.5000,0.5398,0.5793,0.6179,0.7554];a=linspace(0,0.4,21);NUM=21;L=zeros(1,NUM);N=zeros(1,NUM);S=zeros(1,NUM);B=zeros(1,NUM);for i=1:NUMxi=a(i);L(i)=lang(X,Y,xi);% 拉格朗日插值N(i)=newton(X,Y,xi);%牛顿插值B(i)=normcdf(xi,0,1);%原函数S(i)=Threch1(X,Y,0.4,0.36,xi);%三次样条函数第一类边界条件endplot(a,B,'--r');hold on;plot(a,L,'b');hold on;plot(a,N,'r');hold on;plot(a,S,'r+');hold on;legend('原函数','拉格朗日插值','牛顿插值','三次样条插值',2);hold off7增加插值节点观察误差变化main4.mclear;clc;N=5;%4.5第一问Ini=zeros(1,1001);a=linspace(-1,1,1001);Ini=1./(1+25*a.^2);for i=1:3 %节点数量变化次数N=2*N;t=linspace(-1,1,N+1);%插值节点ft=1./(1+25*t.^2);%插值节点函数值val=linspace(-1,1,101);for j=1:101L(j)=lang(t,ft,val(j));S(j)=Threch1(t,ft,0.074,-0.074,val(j));%三样条第一类边界条件插值endplot(a,Ini,'k')%原函数图象hold onplot(val,L,'r')%拉格朗日插值函数图像hold onplot(val,S,'b')%三次样条插值函数图像str=sprintf('插值节点为%d时的插值效果',N);title(str);legend('原函数','拉格朗日插值','三次样条插值');%显示图例hold offfigureend8车门曲线main5.mclearclcX=[0,1,2,3,4,5,6,7,8,9,10];Y=[0.0,0.79,1.53,2.19,2.71,3.03,3.27,2.89,3.06,3.19,3.29]; dy0=0.8;dyn=0.2;n=length(X)-1;d=zeros(n+1,1);h=zeros(1,n-1);f1=zeros(1,n-1);f2=zeros(1,n-2);for i=1:nh(i)=X(i+1)-X(i);f1(i)=(Y(i+1)-Y(i))/h(i);endfor i=2:nf2(i)=(f1(i)-f1(i-1))/(X(i+1)-X(i-1));d(i)=6*f2(i);endd(1)=6*(f1(1)-dy0)/h(1);d(n+1)=6*(dyn-f1(n-1))/h(n-1); A=zeros(n+1,n+1);B=zeros(1,n-1);C=zeros(1,n-1);for i=1:n-1B(i)=h(i)/(h(i)+h(i+1));C(i)=1-B(i);endA(1,2)=1;A(n+1,n)=1;for i=1:n+1A(i,i)=2;endfor i=2:nA(i,i-1)=B(i-1);A(i,i+1)=C(i-1);endM=A\d;x=zeros(1,n);S=zeros(1,n);for i=1:nx(i)=X(i)+0.5;S(i)=Y(i)+(f1(i)-(M(i)/3+M(i+1)/6)*h(i))*(x(i)-X(i))+M(i)/2*(x(i)-X(i ))^2+(M(i+1)-M(i))/(6*h(i))*(x(i)-X(i))^3;endplot(X,Y,'k'); hold on;plot(x,S,'o');title('三次样条插值效果图');legend('已知插值节点','三次样条插值');hold off实验结果:4.31计算插值节点处的函数值xi=0.13时Xi=0.36时2将多种插值函数即原函数图像画在同一张图上4.5.1增加插值节点观察误差变化4.5.2 车门曲线。
三次样条插值多项式实验的目的及意义:为了取得理想结果:在不增加更多的插值条件下,能够求得一个插值多项式,既有良好的逼近效果,又有好的光滑性,引进三次样条插值 多项式。
如果已知函数y=f(x)在节点a=x0<x1<…<xn=b 处的函数值和导数值:()i i x f y =,i=0,1,2,…,n如果S(x)满足条件:1. S(x)是一个分段的三次多项式且()i i y x S =;2. S(x)在[a,b]具有二阶连续导数。
则称S(x)是三次样条插值函数。
S(x)的具体形式为:()()()()⎪⎪⎩⎪⎪⎨⎧∈∈∈=-]12,121,01,[,...............][,][,n n n x x x x s x x x x s x x x x s x s其中()x S i 在[]i i x x ,1-上是三次多项式()iiiiid x c x b x a x S +++=23由插值条件()ii y x S =,i=0,1,2,…,n ,得n+1个条件。
边界条件一:()()nn y x S y x S '',''00== 边界条件二:()()nn y x S y x S '''',''''00==数学公式:()()2211133[2]()[2()]()i i i i i i i i i i ih x x x x h x x x x H x y y h h ---+-----=++2211122()()()()i i i i i i i ix x x x x x x x m m h h -------+ 算法描述:Step1:输入未知数X 及(xi,yi),i=0,1,…,n ; Step2:计算步长H[i]; Step3:计算[][]()⎪⎪⎪⎩⎪⎪⎪⎨⎧+=-=+=-+++i i i i i i i ii i i i i x x f x x f u g u hh h ,,311111λλλStep4:根据边界条件,求解相应的方程得到m0,m1,…, mn Step5:判断X 属于[]i i x x ,1-,i=1,2,…,n 中的哪一个Step6:计算()x s y i i ≈Step7:输出y. 程序原代码如下: #include "stdio.h" #define N 5 void main() { int i,k; float X,s,y0,yn;float a[N][N+1],h[N],u[N],v[N],g[N],m[N],p[N],q[N],w[N];printf("please input X:"); //X 为未知数的大小scanf("%f",&X);printf("please input x:"); //输入x的大小for(i=0;i<N;i++)scanf("%f",&a[i][0]);printf("please input y:"); //输入y的大小for(i=0;i<N;i++)scanf("%f",&a[i][1]);for(i=1;i<N;i++)h[i]=a[i][0]-a[i-1][0]; //计算步长for(i=1;i<N;i++){v[i]=h[i+1]/(h[i]+h[i+1]);u[i]=1-v[i];g[i]=3*u[i]*(a[i+1][1]-a[i][1])/h[i+1]+3*v[i]*(a[i][1]-a[i-1][1])/h[i]; }printf("\t(1)已知边界条件1\n");printf("\t(2)已知边界条件2\n");printf("请选择边界条件序号:");scanf("%d",&k);if(k==1){printf("请输入y0和yn的一阶导:"); //输入边界条件一scanf("%f%f",&m[0],&m[N-1]);p[0]=0; //用追赶法求解m[N]q[0]=0;g[1]=g[1]-v[1]*m[0];g[N-2]=g[N-2]-u[N-2]*m[N-1];for(i=1;i<N;i++){w[i]=2-u[i]*p[i-1];p[i]=v[i]/w[i];q[i]=(g[i]-u[i]*q[i-1])/w[i];}m[N-2]=q[N-2];for(i=N-3;i>0;i--)m[i]=q[i]-p[i]*m[i+1];printf("输出m[i]的值:\n");for(i=0;i<N;i++)printf("%f\n",m[i]);for(i=1;i<N;i++) //计算最终结果if(X>a[i-1][0]&&X<a[i][0])s=(h[i]+2*(X-a[i-1][0]))*(X-a[i][0])*(X-a[i][0])*a[i-1][1]/(h[i]*h[i]*h[i]) +(h[i]-2*(X-a[i][0]))*(X-a[i-1][0])*(X-a[i-1][0])*a[i][1]/(h[i]*h[i]*h[i])+ (X-a[i-1][0])*(X-a[i][0])*(X-a[i][0])*m[i-1]/(h[i]*h[i])+(X-a[i-1][0])*(X-a[i-1][0])*(X-a[i][0])*m[i]/(h[i]*h[i]);printf("s(%f)=%f\n",X,s);}if(k==2){printf("请输入y0和yn的二阶导:"); //输入边界条件二scanf("%f%f",&y0,&yn);g[0]=3*(a[1][1]-a[0][1])/h[1]-h[1]*y0/2;g[N-1]=3*(a[N-1][1]-a[N-2][1])/h[N-1]+h[N-1]*yn/2;q[0]=g[0];u[0]=1;v[N-1]=1;w[0]=2;for(i=1;i<N;i++){w[i]=2-v[i]*u[i-1]/w[i-1];q[i]=g[i]-v[i]*q[i-1]/w[i-1];}m[N-1]=q[N-1]/w[N-1];for(i=N-2;i>=0;i--)m[i]=(q[i]-u[i]*m[i+1])/w[i];printf("输出m[i]的值:\n");for(i=0;i<N;i++)printf("%f\n",m[i]);for(i=1;i<N;i++)if(X>=a[i-1][0]&&X<=a[i][0])s=(h[i]+2*(X-a[i-1][0]))*(X-a[i][0])*(X-a[i][0])*a[i-1][1]/(h[i]*h[i]*h[i]) +(h[i]-2*(X-a[i][0]))*(X-a[i-1][0])*(X-a[i-1][0])*a[i][1]/(h[i]*h[i]*h[i])+ (X-a[i-1][0])*(X-a[i][0])*(X-a[i][0])*m[i-1]/(h[i]*h[i])+(X-a[i-1][0])*(X-a[i-1][0])*(X-a[i][0])*m[i]/(h[i]*h[i]);printf("s(%f)=%f\n",X,s);}}数值计算:已知y=f(x)的如下数值求三次样条插值函数S(x),满足条件1.s’(0)=0,s’(4)=482.s’’(0)=0,s’’(4)=24Please input X:2.5Please input x:0 1 2 3 4 Please input y:-8 -7 0 19 56(1)已知边界条件1(2)已知边界条件2请选择边界条件的序号:1请输入y0和yn的一阶导:0 48 0.0000003.00000012.00000027.00000048.000000s(2.500000)=7.625000press any key tocontinue请选择边界条件的序号:2请输入y0和yn的二阶导:0 24 -0.0000003.00000011.99999927.00000248.0000007.625000press any key tocontinue s(2.500000)=7.625000 对计算结果进行评价分析:()()443845h M x S x f ≤-三次样条插值函数与三次Hermite 插值函数相比,不仅光滑度有提高,而且要求求解时还不需要增加内节点处的导数值,因此比较实用。
实验报告:牛顿差值多项式&三次样条问题:在区间[-1,1]上分别取n=10、20用两组等距节点对龙格函数21()25f x x作多项式插值及三次样条插值,对每个n 值,分别画出插值函数及()f x 的图形。
实验目的:通过编程实现牛顿插值方法和三次样条方法,加深对多项式插值的理解。
应用所编程序解决实际算例。
实验要求:1. 认真分析问题,深刻理解相关理论知识并能熟练应用;2. 编写相关程序并进行实验;3. 调试程序,得到最终结果;4. 分析解释实验结果;5. 按照要求完成实验报告。
实验原理:详见《数值分析 第5版》第二章相关内容。
实验内容:(1)牛顿插值多项式1.1 当n=10时:在Matlab 下编写代码完成计算和画图。
结果如下:代码:clear allclcx1=-1:0.2:1;y1=1./(1+25.*x1.^2);n=length(x1);f=y1(:);for j=2:nfor i=n:-1:jf(i)=(f(i)-f(i-1))/(x1(i)-x1(i-j+1));endendsyms F x p ;F(1)=1;p(1)=y1(1);for i=2:nF(i)=F(i-1)*(x-x1(i-1));p(i)=f(i)*F(i);endsyms PP=sum(p);P10=vpa(expand(P),5);x0=-1:0.001:1;y0=subs(P,x,x0);y2=subs(1/(1+25*x^2),x,x0);plot(x0,y0,x0,y2)grid onxlabel('x')ylabel('y')P10即我们所求的牛顿插值多项式,其结果为:P10(x)=-220.94*x^10+494.91*x^8-9.5065e-14*x^7-381.43*x^6-8.504e-14*x^5+123.36*x^4+2.0 202e-14*x^3-16.855*x^2-6.6594e-16*x+1.0并且这里也能得到该牛顿插值多项式的在[-1,1]上的图形,并和原函数进行对比(见Fig.1)。
Fig.1 牛顿插值多项式(n=10)函数和原函数图形从图形中我们可以明显的观察出插值函数在两端点处发生了剧烈的波动,产生了极大的误差,即龙格现象,当n=20时,这一现象将更加明显。
1.2 当n=20时:对n=10的代码进行修改就可以得到n=20时的代码。
将“x1=-1:0.2:1;”改为“x1=-1:0.1:1;”即可。
运行程序,我们得到n=20时的牛顿插值多项式,结果为:P20(x)= 260188.0*x^20 - 1.0121e6*x^18 + 2.6193e-12*x^17 + 1.6392e6*x^16 + 2.248e-11*x^15 - 1.4429e6*x^14 - 4.6331e-11*x^13 + 757299.0*x^12 + 1.7687e-11*x^11 - 245255.0*x^10 + 2.1019e-11*x^9 + 49318.0*x^8 + 3.5903e-12*x^7 - 6119.2*x^6 - 1.5935e-12*x^5 + 470.85*x^4 + 1.3597e-14*x^3 - 24.143*x^2 - 1.738e-14*x + 1.0同样的,这里得到了该牛顿插值多项式的在[-1,1]上的图形,并和原函数进行对比(见Fig.2)。
Fig.2牛顿插值多项式(n=20)函数和原函数图形当n=20时,端点处发生了更加剧烈的震荡。
表明随着分段不断增加,牛顿插值多项式与原函数的误差不但没有减少,反而变得更大了。
(2)三次样条2.1 当n=10时:在Matlab下编写代码完成计算和画图。
结果如下:代码:clear allclcx1=-1:0.2:1;y1=1./(1+25.*x1.^2);syms xm1=subs(diff(1/(1+25*x^2)),-1);m2=subs(diff(1/(1+25*x^2)),1);n=length(x1);syms a b h f dfor i=1:n-1h(i)=x1(i+1)-x1(i);f(i)=(y1(i+1)-y1(i))/(x1(i+1)-x1(i));enda(n)=1;b(1)=1;for i=2:n-1a(i)=h(i-1)/(h(i-1)+h(i));b(i)=h(i)/(h(i-1)+h(i));endd(1)=6/h(1)*(f(1)-m1);d(n)=6/h(n-1)*(m2-f(n-1));for i=2:n-1d(i)=6*(f(i)-f(i-1))/(h(i-1)+h(i));endD=d';A=2.*eye(n);for i=1:n-1A(i,i+1)=b(i);A(i+1,i)=a(i+1);endM=A^-1*D;for i=1:n-1s(i)=M(i)*(x1(i+1)-x)^3/h(i)/6+M(i+1)*(x-x1(i))^3/h(i)/6+(y1(i)-M(i)*h(i)^2/6)*(x1( i+1)-x)/h(i)+(y1(i+1)-M(i+1)*h(i)^2/6)*(x-x1(i))/h(i);endS=vpa(expand(s.'),5);for i=1:n-1x0=-1-(2/(n-1))+(2/(n-1))*i:0.001:-1+(2/(n-1))*i;y0=subs(s(i),x,x0);plot(x0,y0)hold onendy2=subs(1/(1+25*x^2),x,-1:0.001:1);plot(-1:0.001:1,y2,'r')grid onxlabel('x')ylabel('y')S即为我们所求的三次样条,其结果为:S10(x) =0.08225*x^3+0.36953*x^2+0.56627*x+0.31745[-1,-0.8]0.96279*x^3+2.4828*x^2+2.2569*x+0.76829[-0.8,-0.6]0.81773*x^3+2.2217*x^2+2.1002*x+0.73696[-0.6,-0.4]13.413*x^3+17.336*x^2+8.1461*x+1.5431[-0.4,-0.2]-54.471*x^3-23.394*x^2-1.8741e-17*x+1.0[-0.2,0]54.471*x^3-23.394*x^2+1.9683e-17*x+1.0[0,0.2]-13.413*x^3+17.336*x^2-8.1461*x+1.5431[0.2,0.4]-0.81773*x^3+2.2217*x^2-2.1002*x+0.73696[0.4,0.6]-0.96279*x^3+2.4828*x^2-2.2569*x+0.76829[0.6,0.8]-0.08225*x^3+0.36953*x^2-0.56627*x+0.31745[0.8,1]并且这里可以得到该三次样条的在[-1,1]上的图形,并和原函数进行对比(见Fig.3)。
Fig.3 三次样条(n=10)函数和原函数图形从图形我们可以看出,三次样条图形和原函数图形非常接近,误差相对较小。
2.2 当n=20时:同样的,将上面代码中的“x1=-1:0.2:1;”改为“x1=-1:0.1:1;”即可。
运行程序,我们可以得到n=20时的三次样条,其结果为S20(x) =0.16461*x^3+0.59746*x^2+0.77505*x+0.38066[-1,-0.9]0.27191*x^3+0.88717*x^2+1.0358*x+0.45888[-0.9,-0.8]0.46379*x^3+1.3477*x^2+1.4042*x+0.55712[-0.8,-0.7]0.86962*x^3+2.1999*x^2+2.0008*x+0.69632[-0.7,-0.6]1.5804*x^3+3.4792*x^2+2.7683*x+0.84984[-0.6,-0.5]3.5442*x^3+6.425*x^2+4.2412*x+1.0953[-0.5,-0.4]5.7284*x^3+9.046*x^2+5.2896*x+1.2351[-0.4,-0.3]12.534*x^3+15.171*x^2+7.1273*x+1.4189[-0.3,-0.2]-32.789*x^3-12.023*x^2+1.6884*x+1.0563[-0.2,-0.1]-89.07*x^3-28.907*x^2+5.1067e-17*x+1.0[-0.1,0]89.07*x^3-28.907*x^2-3.9148e-17*x+1.0[0,0.1]32.789*x^3-12.023*x^2-1.6884*x+1.0563[0.1,0.2]-12.534*x^3+15.171*x^2-7.1273*x+1.4189[0.2,0.3]-5.7284*x^3+9.046*x^2-5.2896*x+1.2351[0.3,0.4]-3.5442*x^3+6.425*x^2-4.2412*x+1.0953[0.4,0.5]-1.5804*x^3+3.4792*x^2-2.7683*x+0.84984[0.5,0.6]-0.86962*x^3+2.1999*x^2-2.0008*x+0.69632[0.6,0.7]-0.46379*x^3+1.3477*x^2-1.4042*x+0.55712[0.7,0.8]-0.27191*x^3+0.88717*x^2-1.0358*x+0.45888[0.8,0.9]-0.16461*x^3+0.59746*x^2-0.77505*x+0.38066[0.9,1.0]并且这里也能得到该三次样条的在[-1,1]上的图形,并和原函数进行对比(见Fig.4)。
Fig.4 三次样条(n=20)函数和原函数图形当分段数达到20时,三次样条的图像与原函数基本重合,表明随着分段数的增加,三次样条的误差也不断减少。
友情提示:本资料代表个人观点,如有帮助请下载,谢谢您的浏览!。