三次样条插值实验报告
- 格式:pdf
- 大小:304.11 KB
- 文档页数:9
实验报告:牛顿差值多项式&三次样条问题:在区间[-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)。
实验目的:学会运用三次样条函数进行函数插值实验内容:给定插值条件如下:i 0 1 2 3 4 5 6 7X i8.125 8.4 9.0 9.485 9.6 9.959 10.166 10.2 Y i0.0774 0.099 0.280 0.60 0.708 1.200 1.800 2.177 作三次样条函数插值,取第一类边界条件Y0’=0.01087 Y7’=1001.画出插值曲线的图像。
程序如下:%将插值点的坐标值输入matlabx=[8.125 8.4 9.0 9.485 9.6 9.959 10.166 10.2];y=[0.0774 0.099 0.280 0.60 0.708 1.200 1.800 2.177];%调用三次样条函数,并给出第一类边界条件即 Y0’=0.01087 Y7’=100 cs = spline(x,[0.01087 y 100]);%给定x(0)到x(7)之间的作图点数,即1000xx=linspace(8.125,10.2,1000);%绘制图像plot(x,y,'o',xx,ppval(cs,xx),'-')grid on插值曲线:2.逆时针旋转座标轴45o保持(1.)中结点和边界条件的几何关系不变,再次作三次样条函数插值,画出插值曲线的图像。
程序如下:x=[8.125 8.4 9.0 9.485 9.6 9.959 10.166 10.2];y=[0.0774 0.099 0.280 0.60 0.708 1.200 1.800 2.177];%用坐标转换求出坐标轴旋转45度后插值点新的坐标值x1=x*cos(pi/4)+y*sin(pi/4);y1=-x*cos(pi/4)+y*sin(pi/4);%求出坐标轴转换后新的边界值Y0=(0.01087-1)/(0.01087+1);Y7=(100-1)/(100+1);cs = spline(x1,[Y0 y1 Y7]);xx=linspace(x1(1),x1(8),1000);plot(x1,y1,'o',xx,ppval(cs,xx),'-')grid on其中在求解边界条件用转换关系tan(45)θ-,θ为在原图像中边界点的切线与原x坐标轴的夹角差值曲线:。
四、三次样条插值1. 样条函数插值的原理给定区间[a,b]上划分A:a=x<x<<x<x=b,若分段函数S(x)满足:01n-1n1.S(x)在各个子区间[x,x],i=0,1,,n-1上均为x的三次多项式;ii+12.S(x)在整个区间[a,b]上有直至二阶的连续导数。
则称S(x)为[a,b]上依次划分的三次样条函数,简称样条函数。
具体地有分段表达式:ax3+bx2+cx+d,x G[x,x]000001ax3+bx2+cx+d,x G[x,x]111112S(x)=\ax3+bx2+cx+d,x G[x,x](1)222223ax3+bx2+cx+d,x G[x,x]、°*n-1n—T•••n-1n-1n-1n共有4n个参数a,b,c,d,i=0,1,,n,它们在内节点处满足iiii'S(x)=S(x),…i-0i+0<S'(x)=S'(x),i=1,2,,n-1.(2)i-0i-0S''(x)=S''(x),Ji-0i+0满足样条函数定义的函数集合称为分划A上的三次样条函数空间,记为S(3,A),可以证明S(3,A)为线性空间。
若S(x)G S(3,A),且进一步满足插值条件S(x)=y=f(x),i=0,1,,n(3)iii其中y为节点x处的给定函数值(若被插函数了(x)已知;••则用了(x)代替之),iii则称S(x)为以x,x,,x,x为节点的三次样条函数。
01n-1n其中式(3)插值节点提供了n+1个约束条件;加上式(2)的3n-3个,合起来共有4n-2个;欲求4n个待定参数的唯一解;尚缺两个条件。
这两个条件一般由样条函数的边界条件提供。
常用三类边界条件;他们分别与三次样条函数;构成不同边界条件的样条函数插值问题。
2. 三类样条函数插值问题2.1第二类边界条件给定边界条件两端的一阶导数值:S'(x)=y'=m,S'(x)=y'=m000nnn这相当于样条两短处的方向给定(压铁在两端点的压力方向确定),对应的插值问题如下:对于分划A:a=x<x<<x<x=b,给定节点对应的函数值01n—1ny,y,y,,y,以及两端点处的一阶导数值y'=m,y'=m,求三次样条函数012n00nnS(x),使…f S(x)=y,i=0,1,,n2iiI S'(x)=m,S'(x)=mJ00n…n2.2第一类边界条件给定边界两端的二阶导数值:S''(x)=y''=M,S''(x)=y''=M000nnn这相当于在样条两端处外加一个力矩,使梁两端点处有相应的曲率。
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在矩阵运算上的优势。
Lab03.三次样条插值实验【实验目的和要求】1.使学生深入理解三次样条插值法,深入进行程序设计能力训练;2.对第一与第二种边界条件,按三弯矩法,通过用Matlab 语言设计计算三次样条插值的程序,以提高学生程序设计的能力。
【实验内容】1.根据Matlab 语言特点,描述三次样条插值法。
2.对第一与第二种边界条件,按三弯矩法,用Matlab 语言设计计算三次样条插值的程序。
3对(1) 自然边界条件0)0.1()2.0(=''=''S S ;(2) 第一种边界条件55741.1)0.1( ,20271.0)2.0(='='S S . 输出用追赶法解出的弯矩向量),,(521M M M 和)1.02.0(i S + (i =0,1,…,8)的值,并画出)(x S y =的图形。
4.完成教材P45例8的计算,并将计算结果与Langrage 插值法计算的结果进行比较,由此说明三次样条插值的优越性。
【实验仪器与软件】1.CPU 主频在1GHz 以上,内存在128Mb 以上的PC ;2.Matlab 6.0及以上版本。
实验讲评:实验成绩:评阅教师:2011 年 月 日Lab03.三次样条插值实验一、算法描述1.根据Matlab语言特点,描述三次样条插值法.答:S(x) 在[x j, x j+1](j=1,2,⋯,n-1)上是三次多项式,于是S"(x)在[x j, x j+1] 上是一次多项式,如果S"(x) 在[x j,x j+1](j=1,2,⋯,n-1)两端点上的值已知,设S"(x j)=M j,S"(x j+1)=M j+1,则S"(x) 的表达式为:= ,其中h j =x j+1-x j,对S"(x) 进行两次积分,则得到1 个具有2二、程序设计2.对第一与第二种边界条件,按三弯矩法,用Matlab语言设计计算三次样条插值的程序。
三次样条插值的数值实验姓名: 王维滨 学号:0842011157 姓名: 李佳乐 学号:0842011034 姓名: 谢朝 学号:0842011062 姓名: 杨其荣 学号:08420110721.实验项目的性质和任务对三次样条插值进一步理解,并编写matlab 程序,实现这些功能2.算法设计和matlab 编程。
总前提:x i 第i 点的横坐标i y 第i 点的纵坐标i M ,,s 的记号,,,,,,23()()()()()()()2!3!i i i i i i i s x s x s x y s x x x x x x x =+-+-+- 有数值逼近书上的推导,我们令:111i i i i i x x u x x -+--=-,111i i i i i x x x x λ++--=-,11111111f (,,)i i i i i i i i i i i i i y y y y x x x x x x x x x +-+--++------=- 由于未知数的数目多于方程的个数,我们需要增加两个条件才能唯一确定一个分段三次函数1)D1的三次样条插值a .实验方案与原理:我们加上条件:,,,,11()(),()()n n s x f x s x f x ==我们建立三弯矩方程组:1211211111126(,,)26(,,),2,3.....126(,,)i i i i i i i i n n n n n M M f x x x u M M M f x x x i n M M f x x x λ-+-+--+=⎧⎪++==-⎨⎪+=⎩然后采用追赶法迭代求方程组,但是我们在程序中采用简单的方法(矩阵计算)直接求解降低编程难度,2)D2三次样条插值223123211i 11112111126(,,)26(,,),3,4.....226(,,)i i i i i i i n n n n n n n n M M f x x x u M u M M M f x x x i n u M M f x x x M λλλ-+-+----+-+=-⎧⎪++==-⎨⎪++=-⎩3)D3三次样条插值22321231i 111211126(,,)26(,,),3,4.....126(,,)n i i i i i i i n n n n n n n M M u M f x x x u M M M f x x x i n M u M M f x x x λλλ-+-+--+++=⎧⎪++==-⎨⎪++=⎩b .编写程序:见附录调用test,每次显示一个图像,关闭后按回车继续显示下一幅。
三次样条实验报告范文三次样条插值的实验报告范文湘潭大学实验报告课程名称计算机图形学实验名称参数三次样条的绘制页数专业计算机科学与技术班级一班同组者姓名学号2022551208姓名刘兆臣实验日期2022.05.05实验目的使学生掌握三次参数样条的定义、画法和程序的编写。
实验内容和要求给定型值点,要求用VC++6.0画出通过给定型值点的参数三次样条曲线。
实验方案设计给定型值点,计算出参数三次样条曲线的每个区间段的代数式,由Hermit曲线定义画出每个区间的曲线。
开始开始给定型值点给定型值点计算出参数计算出参数ifn=1Yifn=1N给系数矩阵赋值给系数矩阵赋值Ifn=2YIfn=2给定系数值给定系数值N求解方程组,计算出各点导数求解方程组,计算出各点导数Ifi<ni=0Ifi<ni=0Y计算第i段的二次项,三次项系数计算第i段的二次项,三次项系数结束结束t=0t=0Ift<tt[i-1]Ift<tt[i-1]NYt=t+et=t+e计算出第i段中各点的某,y值计算出第i段中各点的某,y值连线段连线段i++i++程序运行和实验结果说明和分析。
使用VC++6.0运行程序后得到如下图形图为四段曲线组成的三次参数样条曲线,其中各段的绘制是通过给定了型值点的相关参数计算出型值点的导数,再由Hermit曲线知识,在各段上以直代曲绘制出每段的图形。
基本达到了实验目的,完成实验要求。
性能、扩展性等方面存在的不足和可能的改进之处。
不足:在源代码中的n个点采用了数据初始化的方法给出且给定了型值点的个数。
由于给定型值点较少,三次参数样条曲线看起来不够明显。
可改进:可将初始化的型值点数据去除,采用手动键盘输入或文件输入的方法导入多个型值点数据。
附件一源程序,执行程序,符号列表文件。
#include<graphic.h>#include<math.h>#include<conio.h>main(){intgdriver=DETECT,gmode;float某[100],y[100],a[100],b[100],c[100];floatp某[100],py[100],q某[100],qy[100],tt[100];floatd某[100],dy[100];inti,n=4,t,e=3;floatb某3,b某4=0.0,by3=0.0,by4,c某,cy;initgraph(&gdriver,&gmode,"");for(i=0;i<n;i++) {a[i]=0.0;b[i]=0.0;c[i]=0.0;p某[i]=0.0;py[i]=0.0;d某[i]=0.0;dy[i]=0.0;tt[i]=0.0;q某[i]=0.0;qy[i]=0.0;}p某[0]=1.0;py[0]=1.0;p某[4]=1.0;py[4]=1.0;某[0]=10.0;y[0]=110.0;某[1]=40.0;y[1]=100.0;某[2]=80.0;y[2]=90.0;某[3]=130.0;y[3]=95.0;某[4]=200.0;y[4]=105.0;moveto(某[0],y[0]);for(i=0;i<n;i++)putpi某el(某[i],y[i],15);putpi某el(某[0],y[0],15);for(i=0;i<n;i++)tt[i]=qrt((某[i]-某[i-1])某(某[i]-某[i-1])+(y[i]-y[i-1])某(y[i]-y[i-1]));if(n==1)gotopO;for(i=1;i<=n-1;i++){a[i]=2某(tt[i]+tt[i+1]);b[i]=tt[i+1];c[i]=tt[i];d某[i]=3某(tt[i]某(某[i+1]-某[i])/tt[i+1]+tt[i+1]某(y[i]-y[i+1])/tt[i]);}d某[i]=d某[1]-tt[2]某p某[0];d某[n-1]=d某[n-1]-tt[n-1]某p某[n];dy[1]=dy[1]-tt[2]某py[0];dy[n-1]=dy[n-1]-tt[n-1]某py[n];if(n==2){p某[1]=d某[1]/a[1];py[1]=dy[1]/a[1];gotopO;}c[1]=c[1]/a[1];for(i=2;i<=n-1;i++){a[i]=a[i]-b[i]某c[i-1];c[i]=c[i]/a[i];}q某[1]=d某[1]/a[1];qy[1]=dy[1]/a[1];for(i=2;i<=n-1;i++){q某[i]=(d某[i]-b[i]某q某[i-1])/a[i]; qy[i]=(dy[i]-b[i]某qy[i-1])/a[i];}p某[n-1]=q某[n-1];qy[n-1]=qy[n-1];for(i=n-2;i>=1;i--){p某[i]=q某[i]-c[i]某p某[i+1];py[i]=qy[i]-c[i]某py[i+1];}pO:for(i=0;i<=n-1;i++){b某3=(3某(某[i+1]-某[i])/tt[i+1]-2某p某[i]-p某[i+1])/tt[i+1];b某4=((2某(某[i]-某[i+1])/tt[i+1]+p某[i]+p某[i+1])/tt[i+1])/tt[i+1];by3=(3某(y[i+1]-y[i])/tt[i+1]-2某py[i]-py[i+1])/tt[i+1];by4=((2某(y[i]-y[i+1])/tt[i+1]+py[i]+py[i+1])/tt[i+1])/tt[i+1];t=0;do{t=t+e;c某=某[i]+(p某[i]+(b某3+b某4某t)某t)某t;cy=y[i]+(py[i]+(by3+by4某t)某t)某t;lineto(c某,cy);}while(t<tt[i+1]);}getch();cloegraph();}某[i]i型值点某坐标y[i]i型值点Y坐标a[i]初始赋值方程组系数矩阵mi,i的值b[i]初始赋值方程组系数矩阵mi,i+1的值c[i]初始赋值方程组系数矩阵mi,i-1的值p某[i]i型值点导函数某值py[i]i型值点导函数Y值q某[i],qy[i],d某[i],dy[i]均为解方程组中的各项系数tt[i]第i段参数范围(型值点i-1到型值点i的距离)b某,y3每段函数中二次项的系数b某,y4每段函数中三次项的系数c某每段函数中各点的某值cy每段函数中各点的Y值e作图时每段‘以直代曲’中的参数增量其中方程矩阵形式为:m1,1m1,2p1C1m2,1m2,2m2,3p2C2m3,2m3,3m3,4p3C3............=...mn-1,n-2mn-1,n-1mn-1,npn-1Cn-1mn,n-1mn,npnCn将系数矩阵改写为:其中方程矩阵形式为:(变量与代码变量不对应,如a2不等于a[2]) m1,1m1,2l11u1m2,1m2,2m2,3a2l21u2m3,2m3,3m3,4a3l31u3.........=.........mn-1,n-2mn-1,n-1mn-1,nan-1ln-11un-1mn,n-1mn,nanln1pi对应代码中p某[i],py[i]的值Ci对应代码中q某[i],qy[i]的值附件二运行结果文件。