分段三次hermite插值C程序
分段线性插值和分段Hermit插值课程设计
一、前言本文建立在数值分析的理论基础上，能够在Matlab 环境中运行，给出了理论分析、具体实例、程序清单以及程序运行结果，对设计任务中的函数进行了分段线性插值和分段三次Hermit 插值，分别画出分段线性插值多项式和分段三次Hermit 插值多项式的图，最后对着两种不同类型的多项式进行比较和误差分析，找出这两种不同插值方法各自的优劣。
发现分段三次Hermit 插值比分段线性插值的效果要好，在步长越小时分段三次Hermit 插值与插值函数逼近效果更明显，相应的误差越小，而分段线性插值在步长越小时在个别点会出现较大的误差，但总体效果还是可以的，三次Hermit 插值总体上比分段线性插值更光滑，这也符合理论。
二、具体理论知识点(1)分段线性插值近似一条曲线的最简单的方法是过曲线上若干点作一条折线，这就是分段线性插值问题，它的确切提法是：设)(x f 在区间 ],[b a 上的差值数据为)(i i x f y = ，i i h x x n i -+=-≤≤1max 10，求一个函数)(x h φ满足：（1）[]b a C x h ,)(∈φ;（2）在每个子区间[]1,+i i x x )1,,1,0(-=n i 上1)(P x h ∈φ；（3） .,,1,0,)(n i y x i i h ==φ。
我们可以用Lagrange 插值的思想来构造分段线性插值函数)(x h φ，设满足上述条件(1)和(2)的所有函数构成的线性空间为h Φ。
先找线性空间 h Φ 的基函数)(,x l i n ),,1,0(n i =,使得:n j i x l ij j i n ,1,0,,)(,=∂=不难得出,)(,x l i n 的表达式为：[][]⎪⎩⎪⎨⎧∈∈--=,,,0,,,)(1101010,n n x x x x x x x x x x x l[][][]⎪⎪⎪⎩⎪⎪⎪⎨⎧∉∈--∈--=+-+++---,,,0,,,,,,)(11111111,i i i i i i i i i i i i i n x x x x x x x x x x x x x x x x x x l 1,2,1-=n i ,[][]⎪⎩⎪⎨⎧∈∈--=---,,,0,,,)(101101,n n n n n n x x x x x x x x x x x l所以有了基函数)(,x l i n ，满足条件(1)-(3)的分段线性插值函数)(x h φ就可写为:()∑=⋅=ni i n i h x l y x 0,)(φ(2)分段三次Hermit 插值分段三次Hermit 插值问题提法为:给定],[b a 上的1+n 个插值节点的插值数，仍记ii n i h x x -+-≤≤=110max ，并设被插值函数)(x f 在这些节点出的函数值()i i x f y = 和导数值()i i x f m '=都已知，要求)(x f 的一个插值函数)(x H h ，使之满足：（1）对n i ,,1,0 =，有i i h y x H =)(，i i h m x H =')(；（2）在每个子区间 []1,+i i x x 上，)(x H h 是不超过三次的多项式。
计算方法-插值方法实验
实验一插值方法一. 实验目的（1）熟悉数值插值方法的基本思想，解决某些实际插值问题，加深对数值插值方法的理解。
（2）熟悉Matlab 编程环境，利用Matlab 实现具体的插值算法，并进行可视化显示。
二. 实验要求用Matlab 软件实现Lagrange 插值、分段线性插值、三次Hermite 插值、Aitken 逐步插值算法，并用实例在计算机上计算和作图。
三. 实验内容1. 实验题目 （1）已知概率积分dxe y xx ⎰-=22π的数据表构造适合该数据表的一次、二次和三次Lagrange 插值公式，输出公式及其图形，并计算x =0.472时的积分值。
答：①一次插值公式：输入下面内容就可以得到一次插值结果 >> X=[0.47,0.48];Y=[0.4937452,0.5027498]; >> x=0.472;>> (x-X(2))/(X(1)-X(2))*Y(1)+(x-X(1))/(X(2)-X(1))*Y(2)ans =0.495546120000000>>②两次插值公式为：输入下面内容就可以得到两次插值结果>> X=[0.46,0.47,0.48];Y=[0.4846555,0.4937452,0.5027498]; >> x=0.472;>>(x-X(2))*(x-X(3))/((X(1)-X(2))*(X(1)-X(3)))*Y(1)+(x-X(1))*(x-X(3))/((X(2)-X(1))*(X(2)-X(3)))*Y(2)+(x-X(2))*(x-X(1))/((X(3)-X(2))*(X(3)-X(1)))*Y(3)i 0123x 0.46 047 0.48 0.49 y0.4846555 0.4937452 0.5027498 0.5116683ans =0.495552928000000>>③三次插值公式为：输入下面内容就可以得到三次插值结果>> X=[0.46,0.47,0.48,0.49];Y=[0.4846555,0.4937452,0.5027498,0.5116683];>> x=0.472;>>(x-X(2))*(x-X(3))*(x-X(4))/((X(1)-X(4))*(X(1)-X(2))*(X(1)-X(3)))*Y(1)+(x-X(4))*( x-X(1))*(x-X(3))/((X(2)-X(4))*(X(2)-X(1))*(X(2)-X(3)))*Y(2)+(x-X(4))*(x-X(2))*( x-X(1))/((X(3)-X(4))*(X(3)-X(2))*(X(3)-X(1)))*Y(3)+(x-X(3))*(x-X(2))*(x-X(1))/(( X(4)-X(1))*(X(4)-X(2))*(X(4)-X(3)))*Y(4)ans =0.495552960000000输入下面内容，绘出三点插值的图：>> X=[0.46,0.47,0.48,0.49];Y=[0.4846555,0.4937452,0.5027498,0.5116683];>> x=linspace(0.46,0.49);>>y=(x-X(2)).*(x-X(3)).*(x-X(4))/((X(1)-X(4))*(X(1)-X(2))*(X(1)-X(3)))*Y(1)+(x-X(4) ).*(x-X(1)).*(x-X(3))/((X(2)-X(4))*(X(2)-X(1))*(X(2)-X(3)))*Y(2)+(x-X(4)).*(x-X(2) ).*(x-X(1))/((X(3)-X(4))*(X(3)-X(2))*(X(3)-X(1)))*Y(3)+(x-X(3)).*(x-X(2)).*(x-X(1) )/((X(4)-X(1))*(X(4)-X(2))*(X(4)-X(3)))*Y(4);>>plot(x,y)（注意上面的“.*”不能用“*”替代）；（2）将区间[-5,5]分为10等份，求作211)(x x f +=的分段线性插值函数，输出函数表达式及其图形，并计算x =3.3152时的函数值。
埃尔米特(Hermite)插值
实验二埃尔米特（Hermite）插值一、实验目的：1.掌握埃尔米特插值算法原理；2.使用C语言编程实现埃尔米特插值算法。
二、实验准备：阅读《数值分析》2.4节二、实验要求：某人从甲地开车去乙地，每隔一段时间对行车距离和速率进行一次采样，得到在n+1 个采样时刻点t i 的里程s i和速率v i（i=0, 1, ..., n）。
要求编程构造埃尔米特插值多项式H2n+1(t)，满足H2n+1(t i)=s i，H'2n+1(t i)=v i，对所有i=0, 1, ..., n成立，并据此计算m个给定时刻的里程和速率。
函数接口定义：void Hermite_Interpolation( int N, double t[], double s[], double v[], int m, double ht[], double hs[], double hv[] );其中N为采样点个数（注意这个N不是公式中的最大下标n，而是等于n+1），采样时刻点t i、里程s i、速率v i分别通过t、s、v传入；m是需要估算的给定时刻的个数，ht传入给定的时刻点，相应计算出的里程和速率应分别存储在hs和hv中。
裁判程序如下：裁判输入数据：20.0 1.00.0 1.00.0 0.050.0 0.2 0.5 0.8 1.030.0 0.5 1.0100.0 170.0 200.030.0 150.0 0.050.0 0.25 0.5 0.75 1.050.0 1.0 2.0 3.0 4.00.0 60.0 160.0 260.0 300.05.0 70.0 100.0 120.0 20.0100.5 1.0 1.5 2.0 2.5 3.0 3.5 3.8 3.95 4.0标准输出数据：0.0000 0.1040 0.5000 0.8960 1.00000.0000 0.9600 1.5000 0.9600 0.0000100.0000 127.9297 170.0000 195.9766 200.000030.0000 165.4688 150.0000 52.9688 0.000030.2222 60.0000 105.9303 160.0000 206.3438 260.0000 307.9764 305.7687 299.9796 300.000062.6024 70.0000 109.0488 100.0000 92.9745 120.0000 41.2374 -44.8421 -16.2783 20.0000#include<stdio.h>#define MAXN 5 /* 最大采样点个数 */#define MAXM 10 /* 最大估算点个数 */void Hermite_Interpolation( int N, double t[], double s[], double v[], int m, double ht[], double hs[], double hv[] ){double l[10],p[10],h1[10],h2[10],x,ll[10],pp[10];int kk;for(kk=0;kk<m;kk++){x=ht[kk];hs[kk]=0;hv[kk]=0;int i;for(i=0;i<N;i++){l[i]=1;ll[i]=1;int j;for(j=0;j<N;j++){if(i!=j){l[i]=l[i]*(x-t[j])/(t[i]-t[j]);}}p[i]=0;pp[i]=0;int k;for(k=0;k<N;k++){if(i!=k){p[i]=p[i]+l[i]/(x-t[k]);pp[i]=pp[i]+ll[i]/(t[i]-t[k]);}}h1[i]=(1-2*pp[i]*(x-t[i]))*l[i]*l[i];h2[i]=(x-t[i])*l[i]*l[i];hs[kk]=hs[kk]+s[i]*h1[i]+v[i]*h2[i];int kkk;for(kkk=0;kkk<N;kkk++){if(x==t[kkk])break;}if(x==t[kkk])hv[kk]=v[kkk];elsehv[kk]=hv[kk]+s[i]*(2*p[i]*l[i]-4*l[i]*p[i]*(x-t[i])*pp[i]-2*pp[i]*l[ i]*l[i])+v[i]*(l[i]*l[i]+2*l[i]*p[i]*(x-t[i]));}}}int main(){int N, m;double t[MAXN], s[MAXN], v[MAXN]; /* 用于构造的数据 */double ht[MAXM], hs[MAXM], hv[MAXM]; /* 用估算的数据 */int i;while ( scanf("%d", &N) != EOF ) {for ( i=0; i<N; i++ )scanf("%lf", &t[i]);for ( i=0; i<N; i++ )scanf("%lf", &s[i]);for ( i=0; i<N; i++ )scanf("%lf", &v[i]);scanf("%d", &m);for ( i=0; i<m; i++ )scanf("%lf", &ht[i]);Hermite_Interpolation( N, t, s, v, m, ht, hs, hv );for ( i=0; i<m; i++ )printf("%.4lf ", hs[i]);printf("\n");for ( i=0; i<m; i++ )printf("%.4lf ", hv[i]);printf("\n\n");}return 0; }。
hermite插值报告
Hermite 插值实验的目的及意义：分段线性插值多项式S(x)在差值区间[a,b]上只能保证连续性，而不光滑。
要想得到在插值区间上光滑的分段线性插值多项式，采用Hermite 插值。
（带有导数的插值多项式）。
如果已知函数y=f(x)在节点a=x0<x1<…<xn=b 处的函数值和导数值：()()n i x f y x f y i i i i ,...,2,1,0,'',===则在小区间],[1i i x x -上有四个插值条件：()11--=i i x f y ， ()i i x f y =()11''--=i i x f y ，()i i x f y ''=，故能构造一个三次多项式()x H i ，并称为三次Hermite 插值多项式。
这时在整个[a,b]上可以用分段三次Hermite 插值多项式来逼近f(x).()()()()10121,21,[,],[]......,[,]n n n H x x x x H x x x x H x H x x x x -∈⎧⎪∈⎪=⎨⎪⎪∈⎩ 数学公式：()()2211133[2]()[2()]()i i i i i i i i i iih x x x x h x x x x H x y y h h ---+-----=++2211122()()()()''i i i i i i iix x x x x x x x y y h h -------+算法描述：Step1: 输入未知数X 及ii i y y x ',,其中i=0,1,…,n ；Step2: For i=0,1,…,n 对于指定X,判断X 是否满足条件1i i x X x -〈〈； Step3:如果满足计算()()2211133[2]()[2()]()i i i i i i i i i iih x x x x h x x x x H x y y h h ---+-----=++2211122()()()()''i i i i i i iix x x x x x x x y y h h -------+，如果不满足不执行循环。
分段线性插值
摘要用函数来表示变量间的数量关系广泛应用于各学科领域，但是在实际问题中，往往是通过实验、观测以及计算等方法，得到的是函数在一些点上的函数值。
如何通过这些离散数据找到函数的一个满足精度要求且便于使用的近似表达式，是经常遇到的问题。
对于这类问题我们解决的方法为插值法，而最常用也最简单的插值方法就是多项式插值。
当然用插值法得到的近似表达式必须满足插值条件即假设给定了n+1个点的自变量的值以及函数值，近似函数必须要过这n+1（x）通个点。
多项式插值，从几何角度看，就是寻求n次代数曲线y=Pn过n+1个点作为f（x）的近似。
但是随着插值节点个数的增加，高次插值多项式的近似效果并不理想。
根据大量实验得出，在进行高次多项式插值时，会出现龙格现象。
龙格（Runge）现象即当n趋于无穷大时，x在某一邻域内,f(x)收敛，而在这个区域外f(x)发散。
因此，为了解决这样的一个问题，我们可以通过缩小插值区间的办法达到减小误差的目的，所以本实验将针对低次分段插值多项式来做具体的讨论和学习。
关键词：龙格现象分段差值1、实验目的1)通过对分段线性插值算法程序的编写，提高自己编写程序的能力2)体会分段线性插值是如何消除龙格现象的。
3)用实验报告的形式展现，提高自己在写论文方面的能力2、算法理论设在节点处的函数值为,i=0,1,,n。
为了提高近似程度，可以考虑用分段线性插值来逼近原函数，这时的插值函数为分段函数：在区间上的线性函数为误差为：易见，是平面上以点为节点的折线，有如下的特点：1.在上为次数不超过一次的多项式；2.；3.；如果,由线性插值的误差公式得到令,则有关于整体误差：可以按如下方式考虑，若记则对任一都有于是，当时，说明分段线性插值收敛于。
3、数值算例已知点坐标如下表所示：xiyi用分段线性插值法，求解当x为时，对应y的值解：具体程序如下所示：#include ""float Fdline(float x[],float y[],float x1,int len){int i=0;float s=0;for(i=0;i<len-1;i++){if(x1>=x[i] && x1<x[i+1])break;}s=(x1-x[i])/(x[i-1]-x[i])*y[i-1]+(x1-x[i-1])/(x[i]-x[i-1])* y[i];return s;}float Fdline(float x[],float y[],float x1,int len);void main(){float x[]={,,,,};float y[]={,,,,};int len=sizeof(x)/sizeof(x[0]);float x1=0;float s=0;printf("请输入要求解的x1的值:\n");scanf("%f",&x1);s=Fdline(x,y,x1,len);printf("经过分段三次Hermite插值的结果为:\n");printf("%f\n",s);}运行结果：5、对结果进行分析根据分段线性插值的原理，可以看出分段线性插值虽然有很好的收敛性质，但却不是光滑的，所以线性插值的结果和实际的结果差距较大。
数值分析实验六（分段三次Hermite插值）
数值分析实验六（分段三次Hermite插值）《数值分析》实验报告实验编号：实验六课题名称：分段三次Hermite插值一、算法介绍给定的函数为f(x)=1/(25*x*x+1)，将给定区间分成10分，得到11个节点：x[0],x[1],...,x[10]，构造插值函数的基函数。
当x在(x[0],x[1])区间上时，H[0] = (x-x[0])*[((x-x[1])/(x[0]-x[1]))^2]。
其余的区间为H[0]=0。
h[0]= [1+2*(x-x[0])/(x[1]-x[0])]*[((x-x[1])/(x[0]-x[1]))^2]。
当x在[x[i-1],x[i]] (i=1,2,3, (9)区间上时，H[i]=(x-x[i])*[((x-x[i-1])/(x[i]-x[i-1]))^2]，h[i]=[1+2*(x-x[i])/(x[i-1]-x[i])]*[((x-x[i-1])/(x[i]-x[i-1]))^2)。
当x在(x[i],x[i+1]](i=1,2,3,…,10)区间上其余的区间为H[i]=(x-x[i])[((x-x[i+1])/(x[i]-x[i+1]))^2],h[i]=[1+2*(x-x[i])/(x[i+1]-x[i])]*[((x-x[i+1 ])/(x[i]-x[i+1]))^2]。
其余区间上均为H[i]=0,h[i]=0(i=1,2,…,10)。
当x在(x[9],x[10])区间上时，H[10] = (x-x[9])(((x-x[10])/(x[9]-x[10]))^2).其余的区间为H[10]=0.h[10]= (1+2*((x-x[9])/(x[10]-x[9])))(((x-x[10])/(x[9]-x[10]))^2).其余区间h[10]=0。
构造函数H(x) =∑(y[i]*h[i]+y'[i]*H[i],(i=0,1,…,10)。
二、程序代码// testV iew.cpp : implementation of the CT estV iew class//#include "stdafx.h"#include "test.h"#include "testDoc.h"#include "testView.h"#ifdef _DEBUG#define new DEBUG_NEW#undef THIS_FILEstatic char THIS_FILE[] = __FILE__;#endif/////////////////////////////////////////////////////////////////// //////////// CTestV iewIMPLEMENT_DYNCREA TE(CTestView, CView)BEGIN_MESSAGE_MAP(CTestView, CView)//{{AFX_MSG_MAP(CTestView)// NOTE - the ClassWizard will add and remove mapping macros here.// DO NOT EDIT what you see in these blocks of generated code!//}}AFX_MSG_MAP// Standard printing commandsON_COMMAND(ID_FILE_PRINT, CView::OnFilePrint)ON_COMMAND(ID_FILE_PRINT_DIRECT, CV iew::OnFilePrint) ON_COMMAND(ID_FILE_PRINT_PREVIEW,CView::OnFilePrintPreview)END_MESSAGE_MAP()/////////////////////////////////////////////////////////////////// //////////// CTestV iew construction/destructionCTestView::CTestV iew(){// TODO: add construction code here}CTestView::~CT estView(){}BOOL CTestView::PreCreateWindow(CREA TESTRUCT& cs){// TODO: Modify the Window class or styles here by modifying // the CREA TESTRUCT csreturn CV iew::PreCreateWindow(cs);}/////////////////////////////////////////////////////////////////// //////////// CTestV iew drawingvoid CTestView::OnDraw(CDC* pDC){CTestDoc* pDoc = GetDocument();ASSERT_V ALID(pDoc);// TODO: add draw code for native data hereint i,j,k;double x,y,p_x,p_y,l,xx[100],f[100],F[100],sum,p_sum;CPen MyPen,*OldPen;pDC->SetViewportOrg(400,400); //定义坐标原点for(i=-500;i<500;i++){pDC->SetPixel(0,i,RGB(0,0,0));pDC->SetPixel(i,0,RGB(0,0,0)); //画出坐标}pDC->TextOut(-210,5,"-1");pDC->TextOut(196,5,"1");//原函数MyPen.CreatePen(PS_SOLID,1,RGB(255,0,0));//定义画笔颜色OldPen=pDC->SelectObject(&MyPen);x=-1.0,y=1/(1+25*x*x);p_x=x*200;p_y=-y*200;pDC->MoveTo(p_x,p_y);for (x=-1.0;x<=1.0;x+=0.0001){y=1/(1+25*x*x);p_x=x*200;p_y=-y*200;pDC->LineT o(p_x,p_y);}pDC->SelectObject(OldPen);MyPen.DeleteObject();//分段三次Hermite插值MyPen.CreatePen(PS_SOLID,1,RGB(0,0,0)); OldPen=pDC->SelectObject(&MyPen); x=-1.0,y=1.0/(1+25*x*x);p_x=x*200;p_y=-y*200;pDC->MoveTo(p_x,p_y);x=-1.0;for(i=0;i<=10;i++){f[i]=1/(1+25*x*x);xx[i]=x;F[i]=-50*x/(1+25*x*x)/(1+25*x*x); //导数x+=0.2;}x=-1.0;for(j=0;j<=1000;j++){sum=0;for(i=0;i<=10;i++){if(x==xx[i]){sum=f[i];p_x=x*200,p_y=-sum*200;pDC->LineT o(p_x,p_y);break;}if(xxx[i]){y=(1+2*(x-xx[i])/(xx[i+1]-xx[i]))*(x-xx[i+1])*(x-xx[i+1])/(xx[i]-xx[i+1])/(xx[i]-xx[i+1]);sum+=f[i]*y;y=(1+2*(x-xx[i+1])/(xx[i]-xx[i+1]))*(x-xx[i])*(x-xx[i])/(xx[i+1]-xx[i])/(xx[i+1]-xx[i]);sum+=f[i+1]*y;y=(x-xx[i])*(x-xx[i+1])*(x-xx[i+1])/(xx[i]-xx[i+1])/(xx[i]-xx[i+1]);sum+=F[i]*y;y=(x-xx[i+1])*(x-xx[i])*(x-xx[i])/(xx[i+1]-xx[i])/(xx[i+1]-xx[i]);sum+=F[i+1]*y;p_x=x*200;p_y=-sum*200;pDC->LineT o(p_x,p_y);break;}}x+=0.002;}pDC->SelectObject(OldPen);MyPen.DeleteObject();/////////////////////////////////////////////////////////////////// //////////// CTestV iew printingBOOL CTestView::OnPreparePrinting(CPrintInfo* pInfo){// default preparationreturn DoPreparePrinting(pInfo);}void CTestView::OnBeginPrinting(CDC* /*pDC*/, CPrintInfo* /*pInfo*/){// TODO: add extra initialization before printing}void CTestView::OnEndPrinting(CDC* /*pDC*/, CPrintInfo* /*pInfo*/){// TODO: add cleanup after printing}/////////////////////////////////////////////////////////////////// //////////// CTestV iew diagnostics#ifdef _DEBUGvoid CTestView::AssertV alid() const{CView::AssertV alid();}void CTestView::Dump(CDumpContext& dc) const{CView::Dump(dc);CTestDoc* CT estV iew::GetDocument() // non-debug version is inline{ASSERT(m_pDocument->IsKindOf(RUNTIME_CLASS(CT estD oc)));return (CT estDoc*)m_pDocument;}#endif //_DEBUG/////////////////////////////////////////////////////////////////// //////////// CTestV iew message handlers三、运算结果截屏红色的曲线为原函数图像，黑色曲线为分段三次Hermite插值曲线四、算法分析上述图像中黑色的曲线为分段分段三次Hermite插值多项式所对应的图像，由图像可看出黑色的分段三次Hermited插值函数图像和拉格朗日、分段线性插值相比与红色被逼近函数的重合度最好，说明分段三次Hermite插值在函数的各节点两边插值函数的导数是相等的，保证了在各节点的平滑性，且不会出现Runge现象。
5.4 三次Hermite插值
由于0 ( x0 ) 0 ( x1 ) 0( x1 ) 0, 故 0 ( x) 含有 ( x x0 )( x x1 ) 因子。可设
0 ( x) c( x x0 )( x x1 )2
2
其中c为待定系数。
1 . 由 0( x0 ) 1, 可得 c 2 ( x0 x1 )
称 H3 ( x) 为三次Hermite插值多项式，上述条件称 为此问题的插值条件。
采用基函数的方法来构造 H3 ( x) 。
将 H3 ( x) 表示为：
H3 ( x) y00 ( x) y11( x) m00 ( x) m11( x)
其中 0 ( x),1 ( x), 0 ( x), 1 ( x)为插值基函数，且均为次数 不超过3的多项式。为满足插值条件，它们应满足
并给出余项公式。
( x1 ) m1, (i 0,1,2), H3 ( xi ) yi , H3 H3 ( x) y00 ( x) y11( x) y22 ( x) m11( x)
函数值
导数值
x0
x1
0
1
x2
0 0
1
x1
0 0 0
1
0 ( x)
1( x)
2 ( x )
作为多项式插值,三次已是较高的次数，次数再高就有 可能发生Runge现象,因此，对有n+1节点的插值问题， 我们可以使用分段两点三次Hermite插值
例1 对给定数表
x
y
x0 y0
x1 y1 m1
x2 y2
y
求一个三次Hermite插值多项式 H3 ( x) 满足条件
( x1 ) m1, (i 0,1,2), H3 ( xi ) yi , H3
三次Hermite插值
04 实例分析
CHAPTER
实例一：已知数据点的插值
总结词
利用已知数据点进行插值，可三次Hermite插值方法，利用已知的数据点来估计未知点的值。这 种方法能够更好地处理数据点的变化，并提高插值的精度。
CHAPTER
插值多项式的构造
定义
Hermite插值法是一种通过已知的离散数据点来构造一个多 项式，使其能够准确地经过这些数据点，并尽可能地平滑地 连接这些点的方法。
构造方法
Hermite插值多项式由两个部分组成，一个是线性函数，另 一个是二次函数。线性函数部分用于确保插值多项式能够准 确地经过数据点，而二次函数部分则用于保证插值多项式的 平滑性。
实例二：未知数据点的插值
总结词
在未知数据点的情况下，可以通过三次 Hermite插值方法，预测并估计未知点的值。
详细描述
在数据点未知的情况下，可以利用三次 Hermite插值方法，根据已知的数据点来预 测和估计未知点的值。这种方法能够为后续 的数据分析和处理提供重要的参考依据。
实例三：复杂函数的插值
三次Hermite插值能够提供高精度的插值结果，特别是在处理
复杂函数时。
稳定性好
02
该方法在处理大数据集时表现出良好的稳定性，不易受到噪声
和异常值的影响。
易于实现
03
三次Hermite插值的算法相对简单，易于在计算机上实现和优
化。
三次Hermite插值的局限性
对初始数据敏感
三次Hermite插值的结果对初始数据的选择 较为敏感，不同的初始数据可能导致不同的 插值结果。
实验六 分段三次Hermite插值画函数图像
实验六: 分段三次Hermite插值画函数图像学号: 姓名:指导老师：马季骕班级:计算机科学与技术（非师范）1、算法说明：分段三次Hermit插值的做法是在每一个小区间上作三次Hermit插值，因此在每一个插值节点上都需要构造两个插值基函数，然后再作它们的线性组合。
分段三次Hermit插值基函数如下：H(x)=Σ(yihi(x)+y’iHi(x))给定的函数为f(x)=1/(25*x*x+1)，将给定区间分成10分，得到11个节点：x[0],x[1],...,x[10]，构造插值函数的基函数。
当x在(x[0],x[1])区间上时，H[0] = (x-x[0])*[((x-x[1])/(x[0]-x[1]))^2]。
其余的区间为H[0]=0。
h[0]= [1+2*(x-x[0])/(x[1]-x[0])]*[((x-x[1])/(x[0]-x[1]))^2]。
当x在[x[i-1],x[i]] (i=1,2,3,...,9)区间上时，H[i]=(x-x[i])*[((x-x[i-1])/(x[i]-x[i-1]))^2]，h[i]=[1+2*(x-x[i])/(x[i-1]-x[i])]*[((x-x[i-1])/(x[i]-x[i-1]))^2)。
当x在(x[i],x[i+1]](i=1,2,3, (10)区间上其余的区间为H[i]=(x-x[i])[((x-x[i+1])/(x[i]-x[i+1]))^2],h[i]=[1+2*(x-x[i])/(x[i+1]-x[i])]*[((x-x[i+1])/(x[i]-x[i+ 1]))^2]。
其余区间上均为H[i]=0,h[i]=0(i=1,2,…,10)。
当x在(x[9],x[10])区间上时，H[10] = (x-x[9])(((x-x[10])/(x[9]-x[10]))^2).其余的区间为H[10]=0.h[10]= (1+2*((x-x[9])/(x[10]-x[9])))(((x-x[10])/(x[9]-x[10]))^2).其余区间h[10]=0。
三次Hermite插值曲线的能量优化
下 面 讨 论 当 自 由 变 量 α1 α 2 β1 β 2 取 何 值 时 ， 曲线
p(u) 的能量函数 f (α1 α 2 β1 β 2) 最小。
之满足如下条件：
ì p(u i ) = p i , p(u i + 1/3) = α1, p(u i + 2/3) = α 2 , p(u i + 1) = p i + 1 í î p′(u i ) = T i , p′(u i + 1/3) = β1, p′(u i + 2/3) = β 2 , p′(u i + 1) = T i + 1 α1 、 α 2 、β1 和 β 2 是四个自由变量。 其中，
光顺性是一个在 CAGD 中应用很普遍又很重要的
概念， 国内外许多学者对此作了大量研究， 提出了很多 光顺方法， 如 Kjellander 法、 量法及最小二乘法等 [13]。其 中， 能量法是一种整体优化方法， 其光顺效果好， 为人们 普遍采用的一种曲线光顺方法。光顺法的关键是： 能量 函数的确定， 优化问题的求解。 对于曲线 p(u) ， 一般选用 ||p(u)″||2 du 和 p ‴(u)du 作 为曲线的能量函数。其中，p″(u) 为 p(u) 的二阶导数， 体现了曲线的曲率因素。 p ‴(u) 为 p(u) 的三阶导数， 体 现了曲线的挠率因素。由于上述两个能量函数不依赖 曲线的参数化， 能取得较好的光顺效果， 且计算量较小， 便于在计算机上实现， 故在曲线光顺优化中得到了普遍 的应用。
f (α1 α 2 β1 β 2) = u ||p″(u)||2 du = u
i
ui + 1
ui + 1 3
i
||p″(u)||2 du + ||p″(u)||2 ห้องสมุดไป่ตู้u
数值分析实验报告Hermite插值法、Runge现象,比较Language插值、分段线性插值、分段三次Hermie插值
山东师范大学数学科学学院实验报告x 0.1 0.5 1 1.5 2 2.5 3y 0.95 0.84 0.86 1.06 1.5 0.72 1.9y' 1 1.5 2 2.5 3 3.5 4求质点在时刻1.8时的速度，并画出插值多项式的图像。
1)运用Hermite插值法画出图像，如图4-1，并求质点在时刻1.8时的速度。
>>clear>>clc>>X=[0.1 0.5 1 1.5 2 2.5 3;0.95 0.84 0.86 1.06 1.5 0.72 1.9;1 1.5 2 2.5 3 3.5 4];>> x=0.1:0.01:3;>> H=Hermite1(X,x);>> plot(x,H)>> hold on>> plot(X(1,:),X(2,:),'r*')>> H1_8=Hermite(X,1.8);>> plot(1.8,H1_8,'go')>> legend('插值图像','原始点','目标点');图4-1二、验证高次插值的Runge现象问题分析和算法设计（一）Language插值代码function [Ln] =Lagrange(X,x)%请输入2*n+1矩阵X，X中第一行每个元素都是插值节点，X中第二行每个元素都是插值节点对应的函数值；%第二章P24例一拉格朗日插值n=size(X,2);d=0;for m=1:1:nif x==X(1,m);d=m;breakendend运行结果和总结 运行结果 例：给定函数55,11)(2≤≤-+=x xx f ； （1） 验证表2-10的误差结果（高次插值的Runge 现象）；（2） 以0.1为步长分别进行Language 插值、分段线性插值、分段三次Hermite插值，画出三种插值函数以及f(x)的图像，比较三种插值结果。
分段Hermite插值
6.6 分段埃尔米特插值及其MATLAB 程序6.6.2 分段埃尔米特插值的MATLAB 程序调用格式一：YI=interp1(X,Y,XI,'pchip')调用格式二：YI=interp1(X,XI,'pchip')例6.6.5 试用MA TLAB 程序计算例6.6.1中在各小区间中点处分段三次埃尔米特插值)(2/1+i n x H 及其相对误差.解 在MATLAB 工作窗口输入程序>> h=0.2;x0=-1:h:1;y0=1./(1+25.*x0.^2); xi=-0.9:h:0.9; fi=1./(1+25.*xi.^2); yi=interp1(x0,y0,xi,'pchip'); Ri=abs((fi-yi)./fi); xi,fi,yi,Ri,i=[xi',fi',yi',Ri'] 运行后屏幕显示各小区间中点x i 处的函数值f i ，插值s i ，相对误差值R i (略).6.6.3 作有关分段埃尔米特插值图形的MATLAB 程序作有关分段埃尔米特插值图形的MATLAB 主程序function H=hermitetx(x0,y0,xi,x,y)H= interp1(x0,y0,xi,'pchip');Hn= interp1(x0,y0,x,'pchip');plot(x0,y0,'o',x,Hn,'-',xi,H,'*',x,y,'-.')legend('节点(xi,yi)', '分段埃尔米特插值函数','插值点(x,H)','被插值函数y')我们也可以直接在在MATLAB 工作窗口编程序，例如，>> x0 =-6:6; y0 =sin(x0); xi = -6:.25:6;yi = interp1(x0,y0,xi,'pchip');x=-6:0.001:6; y=sin(x); plot(x0,y0,'o',xi,yi,x,y,':'), legend('节点(xi,yi)','分段埃尔米特插值函数','被插值函数y=sinx') title(' y=sinx 及其分段埃尔米特插值函数和节点的图形')>> x0 =-6:6; y0 =cos(x0);xi = -6:.25:6;yi = interp1(x0,y0,xi,'pchip');x=-6:0.001:6; y=cos(x); plot(x0,y0,'o',xi,yi,x,y,':'), legend('节点(xi,yi)','分段埃尔米特插值函数','被插值函数y=cosx') title(' y=cosx 及其分段埃尔米特插值函数和节点的图形')例6.6.6 设函数211)(x x f +=定义在区间]5,5[-上，节点(X (i ),f (X (i )))的横坐标向量X 的元素是首项a =-5，末项b =5，公差h =1.5的等差数列，构造三次分段埃尔米特插值函数)(3,x H n .把区间]5,5[-分成20等份，构成20个小区间，用MA TLAB 程序计算各小区间中点i x 处)(3,x H n 的值，并作出节点，插值点，)(x f 和)(3,x H n 的图形.解 在MATLAB 工作窗口输入程序>>x0=-5:1.5:5;y0=1./(1+x0.^2); x1=-4.75:0.5:4.75;x=-5:0.001:5; y=1./(1+x.^2); H= hermitetx(x0,y0,x1,x,y) title('函数y=1/(1+x^2)及其分段埃尔米特插值函数，插值，节点(xi,yi) 的图形')运行后屏幕显示各小区间中点i x 处)(3,x H n 的值，出现节点，插值点，)(x f 和)(3,x H n 的图形(略).例6.6.7 设函数x x x f cos 5.0)(-=定义在区间],[ππ-上，取7=n ，按等距节点构造分段埃尔米特插值函数)(3,7x H ，用MA TLAB 程序计算各小区间中点i x 处)(3,7x H 的值，作出节点，插值点，)(x f 和)(3,7x H 的图形.解 记节点的横坐标7,,2,1,0,7/2, =π=+π-=i h ih x i 插值点)(21121+++=i i i x x x，6,,2,1,0 =i .在MA TLAB 工作窗口输入程序>> h=2*pi/7; x0=-pi:h:pi;y0=0.5.*x0-cos(x0); xi=-pi+h/2:h:pi-h/2;b=max(x0); a=min(x0); x=a:0.001:b;y=0.5.*x-cos(x); H= hermitetx(x0,y0,xi,x,y)title('函数y=0.5x-cos(x)及其分段埃尔米特插值函数，插值，节点(xi,yi) 的图形')运行后屏幕显示各小区间中点i x 处)(3,7x H 的值，出现节点，插值点，)(x f 和)(3,7x H 的图形(略).6.6.4 用MATLAB 计算有关分段埃尔米特插值的误差例6.6.8 设函数22511)(x x f +=定义在区间]1,1[-上，取10=n ，按等距节点构造分段埃尔米特插值函数)(3,x H n ，用MA TLAB 程序在]1,1[-上计算)(max )4(11x fx ≤≤-和)(3,x H n 的误差公式和误差限. 解 在MATLAB 工作窗口输入程序>> syms h,x=-1:0.0001:1;yxxxx=150000000./(1+25.*x.^2).^5.*x.^4-4500000./(1+25.*x.^2).^4.*x.^2+15000./(1+25.*x.^2).^3;myxxxx=max(yxxxx), R=(h^4)* abs(myxxxx/384)运行后输出)(x f 的4阶导数在区间]1,1[-上绝对值的最大值myxxxx 和)(3,x H n 在区间]1,1[-上的误差公式myxxxx 为myxxxx = R =15000 625/16*h^4（4）在MATLAB 工作窗口输入程序>> h=0.2; R =625/16*h^4运行后输出误差限为R =0.06250000000000例6.6.9 设函数))432sin 3tan(cos()(2x x x f ++=定义在区间],[ππ-上，取9=n ，按等距节点构造分段埃尔米特插值函数)(3,x H n .(1)用MA TLAB 程序计算各小区间中点i x 处)(3,x H n 的值，作出节点，插值点，)(x f 和)(3,x H n 的图形； (2) 并用MA TLAB 程序计算各小区间中点处)(3,x H n 的值及其相对误差；(3) 用MA TLAB 程序求)(max )4(x f x ππ≤≤-和)(3,x H n 在区间],[ππ-上的误差公式和各插值的误差限.解 （1）记节点的横坐标9,,2,1,0,9/2, =π=+π-=i h ih x i ，插值点)(21121+++=i i i x x x，8,,2,1,0 =i .在MATLAB 工作窗口输入程序>>h=2*pi/9; x0=-pi:h:pi;y0=tan(cos((sqrt(3)+sin(2*x0))./(3+4*x0.^2)));xi=-pi+h/2:h:pi-h/2;fi=tan(cos((3^(1/2)+sin(2*xi))./(3+4*xi.^2)));b=max(x0); a=min(x0); x=a:0.001:b;y=tan(cos((3^(1/2)+sin(2.*x))./(3+4*x.^2)));Hi= hermite tx (x0,y0,xi,x,y);Ri=abs((fi-Hi)./fi); xi,fi,Hi,Ri,i=[xi',fi',Hi',Ri']title('函数y=tan(cos((sqrt(3)+sin(2x))/(3+4x^2)))及其分段埃尔米特插值函数，插值，节点(xi,yi) 的图形')运行后屏幕显示各小区间中点x i 处的函数值f i ，插值H i ，相对误差值R i ，并且作出节点，插值点，)(x f 和)(3,x H n 的图形(略).（2）在MATLAB 工作窗口输入程序>> syms xy=tan(cos((3^(1/2)+sin(2*x))/(3+4*x^2)));yxxxx=diff(y,x,4),%simple(yxxxx)运行后屏幕显示函数)(x f 的4阶导数)()4(x f，然后将输出的)()4(x f 编程求)(m a x )4(x f x ππ≤≤-和)(3,x H n 及其在区间],[ππ-上的误差限的MA TLAB 程序如下>>syms h,x=-pi:0.0001:pi;yxxxx=-12.*(1.+tan(cos((3.^(1./2)+sin(2.*x))./(3.+4.*x.^2))).^2).^2.*sin((3.^(1./2)+sin(2.*x))./(3.+4.*x.^2)).^3.*(2.*cos(2.*x)./(3.+4.*x.^2)-8.*(3.^(1./2)+sin(2.*x))./(3.+4.*x.^2).^2.*x).^2.*(-4.*sin(2.*x)./(3.+4.*x.^2)-32.*cos(2.*x)./(3.+4.*x.^2).^2.*x+128.*(3.^(1./2)+sin(2.*x))./(3.+4.*x.^2).^3.*x.^2.-8.*(3.^(1./2)+sin(2.*x))./(3.+4.*x.^2).^2)+16.*(1.+tan(cos((3.^(1./2)+sin(2.*x))./(3.+4.*x.^2))).^2).^2.*sin((3.^(1./2)+sin(2.*x))./(3.+4.*x.^2)).^4.*(2.*cos(2.*x)./(3.+4.*x.^2)-8.*(3.^(1./2)+sin(2.*x))./(3.+4.*x.^2).^2.*x).^4.*tan(cos((3.^(1./2)+sin(2.*x))./(3.+4.*x.^2)))+8.*tan(cos((3.^(1./2)+sin(2.*x))./(3.+4.*x.^2))).^3.*(1.+tan(co s((3.^(1./2)+sin(2.*x))./(3.+4.*x.^2))).^2).*sin((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).^4.*(2.*cos(2.*x)./(3+4.*x.^2)-8.*(3.^(1./2)+sin(2.*x))./(3.+4.*x.^2).^2.*x).^4-8.*tan(cos((3.^(1./2)+sin(2.*x ))./(3.+4.*x.^2))).*(1.+tan(cos((3.^(1./2)+sin(2.*x))./(3.+4*x.^2))).^2).*sin((3.^(1./2)+sin(2.*x))./(3.+4.*x.^2)).^2.*(2.*cos(2.*x)./(3.+4.*x.^2)-8*(3.^(1./2)+sin(2.*x))./(3.+4*x.^2).^2.*x).^4+6.*(1+tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).^2).*sin((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).*(2.*cos(2.*x)./(3+4.*x.^2)-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2.*x).^2.*(-4.*sin(2.*x)./(3+4.*x .^2)-32.*cos(2.*x)./(3+4.*x.^2).^2.*x+128.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^3.*x.^2-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2)+(1+tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).^2).*cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).*(2.*cos(2.*x)./(3+4.*x.^2)-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2.*x).^4-3.*(1+tan(cos((3.^(1./2)+si n(2.*x))./(3+4.*x.^2))).^2).*cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).*(-4.*sin(2.*x)./(3+4.*x.^2)-32.*cos(2.*x)./(3+4.*x.^2).^2.*x +128.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^3.*x.^2-8.*(3.^(1./2)+s in(2.*x))./(3+4.*x.^2).^2).^2-4.*(1+tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).^2).*cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).*(2.*cos(2.*x)./(3+4.*x.^2)-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2.*x).*(-8.*cos(2.*x)./(3+4.*x.^2)+96.*sin(2.*x)./(3+4.*x.^2).^2.*x +768.*cos(2.*x)./(3+4.*x.^2).^3.*x.^2-48.*cos(2.*x)./(3+4.*x.^2).^2-3072.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^4.*x.^3+384.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^3.*x)-(1+tan(cos((3.^(1./2)+sin(2.*x ))./(3+4.*x.^2))).^2).*sin((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).*(16.*sin(2.*x)./(3+4.*x.^2)+256.*cos(2.*x)./(3+4.*x.^2).^2.*x-3072.*sin(2.*x)./(3+4.*x.^2).^3.*x.^2+192.*sin(2.*x)./(3+4.*x.^2).^2-24576.*cos(2.*x)./(3+4.*x.^2).^4.*x.^3+3072.*cos(2.*x)./(3+4.*x.^2).^3.*x+98304.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^5.*x.^4-18432.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^4.*x.^2+384.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^3)-12.*(1+tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).^2).^2*sin((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).^2.*(2.*cos(2.*x)./(3+4.*x.^2)-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2.*x).^4.*cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))-24.*tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).^2.*(1+tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).^2).*sin((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).^3.*(2.*cos(2.*x)./(3+4.*x.^2)-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2.*x).^2.*(-4.*sin(2.*x)./(3+4.*x.^2)-32.*cos(2.*x)./(3+4.*x .^2).^2.*x+128.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^3.*x.^2-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2)-24.*tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).^2.*(1+tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).^2).*sin((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).^2.*(2.*cos(2.*x)./(3+4.*x.^2)-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2.*x).^4.*cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))+36.*tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).*(1+tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).^2).*sin((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).*(2.*cos(2.*x)./(3+4.*x.^2)-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2.*x).^2.*cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).*(-4.*sin(2.*x)./(3+4.*x.^2)-32.*cos(2.*x)./(3+4.*x.^2).^2.*x+128.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^3.*x.^2-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2)+6.*tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).*(1+tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).^2).*cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).^2.*(2.*cos(2.*x)./(3+4.*x.^2)-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2.*x).^4+6.*tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).*(1+tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).^2).*s in((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).^2.*(-4.*sin(2.*x)./(3+4.*x.^2)-32.*cos(2.*x)./(3+4.*x.^2).^2.*x+128.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^3.*x.^2-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2).^2+8.*tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).*(1+tan(cos((3.^(1./2)+sin(2.*x))./(3+4.*x.^2))).^2).*sin((3.^(1./2)+sin(2.*x))./(3+4.*x.^2)).^2.*(2.*cos(2.*x)./(3+4.*x.^2)-8.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^2.*x).*(-8.*cos(2.*x)./(3+4.*x.^2)+96.*sin(2.*x)./(3+4.*x.^2).^2.*x+768.*cos(2.*x)./(3+4.*x.^2).^3.*x.^2-48.*c os(2.*x)./(3+4.*x.^2).^2-3072.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^4.*x.^3+384.*(3.^(1./2)+sin(2.*x))./(3+4.*x.^2).^3.*x)myxx=max(yxxxx), R=(h.^4).* abs(myxx./384)将其保存为名为myxx.m 的M 文件，然后在MATLAB 工作窗口输入该文件名>> myxx运行后屏幕显示myxx =)(max )4(x f x π≤≤π-和)(3,x H n 在区间],[ππ-上的误差公式)(max384)4(4x f h R x π≤≤π-=如下myxx = R =73.94706841647552 1734520780029061/9007199254740992*h^4最后在MATLAB 工作窗口输入>>h=2*pi/9; R =1734520780029061/9007199254740992*h^4运行后屏幕显示)(3,x H n 在区间],[ππ-上的误差限R =0.04574453029948。
数值分析实验,用程序实现Hermite插值法
《数值分析》实验报告实验序号：实验六 实验名称: Hermite 插值法1. 实验目的：学会Hermite 插值法,并应用该算法于实际问题.2. 实验内容：求一个函数ϕ(x )用来近似函数f (x )，用分段三次Hermit 插值的方法来求解近似函数ϕ(x )并画出近似函数图像及原函数图像。
设在区间[a,b]上，给定n+1个插值节点b x x x x a n =<<<<=...210和相应的函数值n y y y ,...,,10以及一阶导数值''1'0,...,,ny y y ，求一个插值函数)(x H ，满足以下条件： （1）),...,2,1,0()(,)(''n i y x H y x H i i i i === （2） )(x H 在每一个小区间[1,+j j x x ]上是三次多项式。
对于给定函数11-,2511)(2≤≤+=x x x f 。
在区间[]11-，上画出f (x )和分段三次Hermit 插值函数)(x H 的函数图像。
3. 实验分析：算法分析：1. 分段三次Hermit 插值的算法思想：分段三次Hermit 插值的做法是在每一个小区间上作三次Hermit 插值，因此在每一个插值节点上都需要构造两个插值基函数)(),(x H x h i i ，然后再作它们的线性组合。
分段三次Hermit 插值基函数如下：⎪⎩⎪⎨⎧≤≤----+=其它 0 ))(21()(1021010100x x x x x x x x x x x x h ⎪⎩⎪⎨⎧≤≤---=其它 0 ))(()(10210100x x x x x x x x x x H1,...,2,1 0 ))(21( ))(21()(1211112111-=⎪⎪⎪⎩⎪⎪⎪⎨⎧≤<----+≤≤----+=++++---n i x x x x x x x x x x x x x x x x x x x x x x x h i i i i i i i i i i-i i i i i i i 其它1,...,2,1 0 ))(( ))(()(12111211-=⎪⎪⎪⎩⎪⎪⎪⎨⎧≤<---≤≤---=+++--n i x x x x x x x x x x x x x x x x x x x H i i i i i i i i-i i i i i 其它⎪⎩⎪⎨⎧≤<----+=---其它 0 ))(21()(1-n 2111n n n n n n n n x x x x x x x x x x x x h ⎪⎩⎪⎨⎧≤<---=--其它 0 ))(()(1-n 211n n n n n n x x x x x x x x x x H 分段三次Hermit 插值函数是：∑=+=n i i i i i x H y x h y x H 0'))()(()( 4. 实验代码：// LDlg.cpp : implementation file//#include "stdafx.h"#include "L.h"#include "LDlg.h"#ifdef _DEBUG#define new DEBUG_NEW#undef THIS_FILEstatic char THIS_FILE[] = __FILE__;#endif/////////////////////////////////////////////////////////////////////////////// CAboutDlg dialog used for App Aboutclass CAboutDlg : public CDialog{public:CAboutDlg();// Dialog Data//{{AFX_DATA(CAboutDlg)enum { IDD = IDD_ABOUTBOX };//}}AFX_DATA// ClassWizard generated virtual function overrides//{{AFX_VIRTUAL(CAboutDlg)protected:virtual void DoDataExchange(CDataExchange* pDX); // DDX/DDV support //}}AFX_VIRTUAL// Implementationprotected://{{AFX_MSG(CAboutDlg)//}}AFX_MSGDECLARE_MESSAGE_MAP()};CAboutDlg::CAboutDlg() : CDialog(CAboutDlg::IDD){//{{AFX_DATA_INIT(CAboutDlg)//}}AFX_DATA_INIT}void CAboutDlg::DoDataExchange(CDataExchange* pDX){CDialog::DoDataExchange(pDX);//{{AFX_DATA_MAP(CAboutDlg)//}}AFX_DATA_MAP}BEGIN_MESSAGE_MAP(CAboutDlg, CDialog)//{{AFX_MSG_MAP(CAboutDlg)// No message handlers//}}AFX_MSG_MAPEND_MESSAGE_MAP()///////////////////////////////////////////////////////////////////////////// // CLDlg dialogCLDlg::CLDlg(CWnd* pParent /*=NULL*/): CDialog(CLDlg::IDD, pParent){//{{AFX_DATA_INIT(CLDlg)// NOTE: the ClassWizard will add member initialization here//}}AFX_DATA_INIT// Note that LoadIcon does not require a subsequent DestroyIcon in Win32m_hIcon = AfxGetApp()->LoadIcon(IDR_MAINFRAME);}void CLDlg::DoDataExchange(CDataExchange* pDX){CDialog::DoDataExchange(pDX);//{{AFX_DATA_MAP(CLDlg)// NOTE: the ClassWizard will add DDX and DDV calls here//}}AFX_DATA_MAP}BEGIN_MESSAGE_MAP(CLDlg, CDialog)//{{AFX_MSG_MAP(CLDlg)ON_WM_SYSCOMMAND()ON_WM_PAINT()ON_WM_QUERYDRAGICON()ON_BN_CLICKED(IDC_LARGRI, OnLargri)ON_BN_CLICKED(IDC_BUTTON2, OnButton2)ON_BN_CLICKED(IDC_HERMITE, OnHermite)//}}AFX_MSG_MAPEND_MESSAGE_MAP()///////////////////////////////////////////////////////////////////////////// // CLDlg message handlersBOOL CLDlg::OnInitDialog(){CDialog::OnInitDialog();// Add "About..." menu item to system menu.// IDM_ABOUTBOX must be in the system command range.ASSERT((IDM_ABOUTBOX & 0xFFF0) == IDM_ABOUTBOX);ASSERT(IDM_ABOUTBOX < 0xF000);CMenu* pSysMenu = GetSystemMenu(FALSE);if (pSysMenu != NULL){CString strAboutMenu;strAboutMenu.LoadString(IDS_ABOUTBOX);if (!strAboutMenu.IsEmpty()){pSysMenu->AppendMenu(MF_SEPARATOR);pSysMenu->AppendMenu(MF_STRING, IDM_ABOUTBOX, strAboutMenu);}}// Set the icon for this dialog. The framework does this automatically // when the application's main window is not a dialogSetIcon(m_hIcon, TRUE); // Set big iconSetIcon(m_hIcon, FALSE); // Set small icon// TODO: Add extra initialization herereturn TRUE; // return TRUE unless you set the focus to a control}void CLDlg::OnSysCommand(UINT nID, LPARAM lParam){if ((nID & 0xFFF0) == IDM_ABOUTBOX){CAboutDlg dlgAbout;dlgAbout.DoModal();}else{CDialog::OnSysCommand(nID, lParam);}// If you add a minimize button to your dialog, you will need the code below // to draw the icon. For MFC applications using the document/view model, // this is automatically done for you by the framework.void CLDlg::OnPaint(){if (IsIconic()){CPaintDC dc(this); // device context for paintingSendMessage(WM_ICONERASEBKGND, (WPARAM) dc.GetSafeHdc(), 0);// Center icon in client rectangleint cxIcon = GetSystemMetrics(SM_CXICON);int cyIcon = GetSystemMetrics(SM_CYICON);CRect rect;GetClientRect(&rect);int x = (rect.Width() - cxIcon + 1) / 2;int y = (rect.Height() - cyIcon + 1) / 2;// Draw the icondc.DrawIcon(x, y, m_hIcon);}else{CDialog::OnPaint();}}// The system calls this to obtain the cursor to display while the user drags // the minimized window.HCURSOR CLDlg::OnQueryDragIcon(){return (HCURSOR) m_hIcon;}void CLDlg::OnOK(){int x00=300,y00=350,i,j;double x;CDC *pDC=GetDC();pDC->SetMapMode(MM_LOMETRIC);pDC->SetViewportOrg(x00,y00);//画坐标轴与原函数for(i=-700; i<=700; i++){pDC->SetPixel(i,0,RGB(0,0,0));pDC->SetPixel(0,i,RGB(0,0,0));}for(x=-1; x<=1; x+=0.001){double j=400.0/(1+25*x*x);pDC->SetPixel(x*500,j,RGB(255,0,0));}pDC->TextOut(-30,-10,"0");pDC->TextOut(-30,430,"1");pDC->TextOut(490,-10,"1");pDC->TextOut(-490,-10,"-1");pDC->MoveTo(-10,680); //x箭头pDC->LineTo(0,700);pDC->MoveTo(0,700);pDC->LineTo(10,680);pDC->MoveTo(680,10); //y箭头pDC->LineTo(700,0);pDC->MoveTo(700,0);pDC->LineTo(680,-10);pDC->TextOut(-30,700,"y");pDC->TextOut(700,-10,"x");}void CLDlg::OnLargri(){int x00=300,y00=350,i,j;double x;CDC *pDC=GetDC();pDC->SetMapMode(MM_LOMETRIC);pDC->SetViewportOrg(x00,y00);//画坐标轴for(i=-700; i<=700; i++){pDC->SetPixel(i,0,RGB(0,0,0));pDC->SetPixel(0,i,RGB(0,0,0));}double yx[]={-1,-0.8,-0.6,-0.4,-0.2,0,0.2,0.4,0.6,0.8,1}; pDC->TextOut(-30,-10,"0");pDC->TextOut(-30,430,"1");pDC->TextOut(490,-10,"1");pDC->TextOut(-490,-10,"-1");pDC->MoveTo(-10,680); //x箭头pDC->LineTo(0,700);pDC->MoveTo(0,700);pDC->LineTo(10,680);pDC->MoveTo(680,10); //y箭头pDC->LineTo(700,0);pDC->MoveTo(700,0);pDC->LineTo(680,-10);pDC->TextOut(-30,700,"y");pDC->TextOut(700,-10,"x");// 拉格朗日差值的函数double yy[12],lx[12],ly[12];double l_fenzi[12],l_fenmu[12];double l_x,l_y;for(i=0; i<=10; i++){yy[i]=1.0/(1+25*yx[i]*yx[i]);for(i=0; i<=10; i++){l_fenmu[i]=1.0;for(j=0; j<=10; j++){if(i!=j)l_fenmu[i]=l_fenmu[i]*(yx[i]-yx[j]);}}double qq,pp;for(qq=-1; qq<=1; qq+=0.0003){for(i=0; i<=10; i++){l_fenzi[i]=1.0;for(j=0; j<=10; j++){if(i!=j)l_fenzi[i]=l_fenzi[i]*(qq-yx[j]);}}pp=0;for(i=0; i<=11; i++){pp=pp+1.0/(1+25*yx[i]*yx[i])*l_fenzi[i]/l_fenmu[i];}pDC->SetPixel(qq*500,pp*390+5,RGB(132,112,225));}void CLDlg::OnButton2(){int x00=300,y00=350,i,j;double x;CDC *pDC=GetDC();pDC->SetMapMode(MM_LOMETRIC);pDC->SetViewportOrg(x00,y00);//画坐标轴与原函数for(i=-700; i<=700; i++){pDC->SetPixel(i,0,RGB(0,0,0));pDC->SetPixel(0,i,RGB(0,0,0));}double yx[]={-1,-0.8,-0.6,-0.4,-0.2,0,0.2,0.4,0.6,0.8,1}; double yy[14];for(i=0; i<=10; i++){yy[i]=1.0/(1+25*yx[i]*yx[i]);}pDC->TextOut(-30,-10,"0");pDC->TextOut(-30,430,"1");pDC->TextOut(490,-10,"1");pDC->TextOut(-490,-10,"-1");pDC->MoveTo(-10,680); //x箭头pDC->LineTo(0,700);pDC->MoveTo(0,700);pDC->LineTo(10,680);pDC->MoveTo(680,10); //y箭头pDC->LineTo(700,0);pDC->MoveTo(700,0);pDC->LineTo(680,-10);pDC->TextOut(-30,700,"y");pDC->TextOut(700,-10,"x");// 线性分段差值的图像CPen pen;CPen*oldpen;pen.CreatePen(PS_SOLID,5,RGB(0,0,0));oldpen=pDC->SelectObject(&pen);for(i=0; i<10; i++){pDC->MoveTo(yx[i]*480,yy[i]*400);pDC->LineTo(yx[i+1]*480,yy[i+1]*400); }}void CLDlg::OnHermite(){int x00=300,y00=350,i,j;double x;CDC *pDC=GetDC();pDC->SetMapMode(MM_LOMETRIC);pDC->SetViewportOrg(x00,y00);//画坐标轴与原函数for(i=-700; i<=700; i++){pDC->SetPixel(i,0,RGB(0,0,0));pDC->SetPixel(0,i,RGB(0,0,0));}double yx[]={-1,-0.8,-0.6,-0.4,-0.2,0,0.2,0.4,0.6,0.8,1};double yy[12];for(i=0; i<=10; i++){yy[i]=1.0/(1+25*yx[i]*yx[i]);}pDC->TextOut(-30,-10,"0");pDC->TextOut(-30,430,"1");pDC->TextOut(490,-10,"1");pDC->TextOut(-490,-10,"-1");pDC->MoveTo(-10,680); //x箭头pDC->LineTo(0,700);pDC->MoveTo(0,700);pDC->LineTo(10,680);pDC->MoveTo(680,10); //y箭头pDC->LineTo(700,0);pDC->MoveTo(700,0);pDC->LineTo(680,-10);pDC->TextOut(-30,700,"y");pDC->TextOut(700,-10,"x");//分段三次Hermite差值的函数double x0,x1,yd1,yd0,y1,y0;for(i=0; i<10; i++){x0=yx[i],x1=yx[i+1];y0=1.0/(1+25*x0*x0);y1=1.0/(1+25*x1*x1);yd0=-(50*x0)*1.0/((1+25*x0*x0)*(1+25*x0*x0));yd1=-(50*x1)*1.0/((1+25*x1*x1)*(1+25*x1*x1));for(double qq=x0; qq<x1; qq+=0.005){double pp= y0*(1+2*(qq-x0)/(x1-x0)) * (qq-x1)/(x0-x1) * (qq-x1)/(x0-x1)+y1*(1+2*(qq-x1)/(x0-x1)) * (qq-x0)/(x1-x0) * (qq-x0)/(x1-x0)+yd0*(qq-x0) * (qq-x1)/(x0-x1) * (qq-x1)/(x0-x1)+yd1*(qq-x1) * (qq-x0)/(x1-x0) * (qq-x0)/(x1-x0);pDC->SetPixel(qq*500,pp*400,RGB(225,185,15));}}}5.实验截图6. 实验结果分析：通过本次实验我对分段三次Hermit插值有了更深刻更全面的掌握，它在给定了节点处的函数值和导数值以后，构造了一个整体上具有一阶连续微商的插值函数。
131张艳-埃尔米特(Hermite) 插值逼近的C语言程序
论文题目：埃尔米特（Hermite）插值逼近的C语言程序院系：数学科学学院专业：数学与应用数学*名：**学号：********指导教师：***完成时间：2007-5-15埃尔米特（Hermite ） 插值逼近的C 语言程序张艳包头师范学院数学科学学院摘要：本文主要探讨埃尔米特（Hermite ） 插值逼近的C 语言程序算法，着重分析其推导过程，并给出了其C 语言程序以及埃尔米特（Hermite ）插值逼近的简单应用.关键词：Hermite 插值多项式;插值条件; Hermite 插值基函数.一、Hermite 插值多项式定义定义： 设基点n x x x ,,,10 互异.给定)(i x f ,)('i x f ,n i ,1,0 =.要求插值多项式)(x H 满足)()()()(''i i i i x f x H x f x H == n i ,1,0 = （1）则称)(x H 为二重密切Hermite 插值多项式，简称为Hermite 插值多项式.n x x x ,,,10 称为二重插值基点.（1）式共有22+n 个条件，因此Hermite 插值多项式通常次数不超过12+n ，故可将)(x H 记为)(12x H n +.二、埃尔米特（Hermite ）插值多项式的存在唯一性定理：)(x f 关于互异基点n x x x ,,,10 满足条件（1）的二重密切12+n 次Hermite 插值多项式存在且唯一.证明：设有2n+1次多项式12121012)(++++++=n n n x a x a a x H （2）满足条件（1）即)()()()(''i i i i x f x H x f x H == n i ,1,0 = （3）由（2）式知（3）式是一个关于1210,,,+n a a a 的22+n 阶线性方程组.)(12x H n +的存在唯一性决定于（3）式为齐次线性方程系组，即当)(i x f =0，)('i x f =0（n i ,1,0 =）时，（3）式仅有平凡解)12,,1,0(0+==n i a i .现用反证法证明：若齐次方程组有非平凡解，则表示存在一个次数不高于12+n 的多项式)(12x H n +满足0)()('1212==++i n i n x H x H n i ,1,0 =i x 为)(12x H n +的互异的二重零点，即12+n 次多项式)(12x H n +有22+n 个零点（包括重数），这和代数基本定理相矛盾. 三、埃尔米特（Hermite ）插值多项式的构造由定理知)(12x H n +存在且唯一，我们用类似于拉格朗日插值多项式的构造方法来构造Hermite 插值多项式设)(x i α，)(x i β，n i ,1,0 =，分别满足插值条件⎪⎩⎪⎨⎧==0)()('j i ijj i x x αδα n j ,1,0 =（4） ⎪⎩⎪⎨⎧==ijj i j i x x δββ)(0)('n j ,1,0 =（5） (其中ij δ表示克罗内克（Kronecker ）符号.当j i =， 时ij δ=1；当j i ≠，时ij δ=0.)的12+n 次多项式，于是次数不超过12+n 次多项式)()()()()(0'012x x f x x f x H i ni ＝i n i i i n βα∑∑+==+ 能够满足插值条件⎩⎨⎧==++)()()()(''1212i i n i i n x f x Hx f x Hn i ,1,0 =.因而就是所要求的12+n 次多项式.因此只要构造出满足条件（4）（5）的)(x i α和)(x i β即可.我们把满足插值条件（4）（5）的12+n 次多项式)(x i α和)(x i β（n i ,1,0 =）称为Hermite 插值基函数.下面构造)(x i α和)(x i β 由于关于基点n x x x ,,,10 的拉格朗日基函数)(x l i 满足⎪⎩⎪⎨⎧===0)()(2)]([0)(''22j i j i j i j i x l x l x l x l （ i j ≠，n j ,1,0 =）且)(2x l i 是n 2次多项式，结合插值条件（4）设)(x i α为)()()(2x l B x A x i i i i +=α，则依据条件（4）要求应有⎪⎩⎪⎨⎧=++=++=+==0)()()(2)()]()[()()()()()('2'22'2j i j i i j i j i i j i i j i j i i ji j i i j i ij j i x l x l B x A x l A x l B x A x l A x x l B x A x αδαn j i ,1,0, =当j i ≠时，由于0)(=j i x l 故0)(=j i x α 0)('=j i x α当j i =时应有 ⎩⎨⎧=++==+=0)()(2)(1)(''i i i i i i i i i i i i i x l B x A A x B x A x αα 从中解出i A i B 可得∑∑≠=≠=-+=-=--=-=nik k ki i i i i nik k ki i ii x x x x A B x x x l A 00'121112)(2从而得到 ∑≠=--+=-+=nik k i ki i ii ii i x l x x x x x l x l x x x 022')(]1)(21[)()]()(21[)(α设)(x i β为 )()()(2x l x x C x i i i i -=β则 )()()(2)()('2'x l x l x x C x C l C x i i i i i i i i -+=β 依据条件（5）要求应有⎪⎩⎪⎨⎧=-+==-=ijj i j i i j i j i i j i j i i j i j i x l x l x x C x l C x x l x x C x δββ)()()(2)()(0)()()('2'2n j i ,1,0, = 当j i ≠时，由于0)(=j i x l 故0)(=j i x β 0)('=j i x β当j i =时，故也应有1)()()(2)()('2'=-+=i i i i i i i i i i i i i x l x l x x C x C l C x β 而1)(=i i x l故1=i C 即)()()(2x l x x x i i i -=β 因此，我们得到埃尔米特插值函数的基函数为⎪⎪⎩⎪⎪⎨⎧-=--+=-+=∑≠=)()()()(]1)(21[)()]()(21[)(2022'x l x x x x l x x x x x l x l x x x i i ini k k i k i i i i i i i βα n i ,1,0 = 根据插值条件，利用二重密切的Hermite 插值基函数的性质，Hermite 插值多项式可简单地表示为)()()()(]1)(21[)()()()()()(2'2000'012x l x x x f x l x x x x x f x x f x x f x H i i ni i i n ik k k i i ni i ni i i i n i i n -+--+=+=∑∑∑∑∑=≠====+βα 四、埃尔米特（Hermite ）插值多项式误差在求解某些数学问题时，用有限的过程代替无限过程所产生的误差称为截断误差（或方法误差）.定理：a ：设)(x f 的导数)()12(x f n +于[a,b]连续，)()22(x fn +于(a,b)内存在，],[b a x i ∈（n i ,1,0 =）互异；b:)(12x H n +为Hermite 插值多项式；则 210)22(1212)]())([()!22()()()()(n n n n x x x x x x n f x H x f x R ---+=-=+++ η. 其中],[b a ∈η与x 有关。
实习：Matlab作业hermite插值
题目：利用Matlab实现数据的Hermite插值和分段三次Hermite插值小组成员：王晓波（38）蔡明宇（20）一、程序实现意义：一般的，从各种试验得来的数据总有一定的数量，而利用插值技术能够从有限的数据中获取整体的状态。
而Hermite插值不仅保证了插值函数与原函数在给定数据点处得拟合，同时保证了在相应点处导数的相同，从而在很大程度上保证了曲线的“光滑性”。
因此，通过Matlab实现Hermite插值具有很普遍的意义。
二、实现过程：1、Hermite插值由于并不是所有的Matlab版本都提供现有的Hermite插值函数包，故我们首先编写了实现给定五个观测点的Hermite插值的M程序，代码如下：function [f,f0] = Hermite1(x,y,y_1)syms t;f = ;!if(length(x) == length(y))if(length(y) == length(y_1))n = length(x);elsedisp('y和y的导数的维数不相等');return;endelsedisp('x和y的维数不相等！ ');return;end*for i=1:nh = ;a = ;for j=1:nif( j ~= i)h = h*(t-x(j))^2/((x(i)-x(j))^2);a = a + 1/(x(i)-x(j));endendf = f + h*((x(i)-t)*(2*a*y(i)-y_1(i))+y(i));<endf0 = subs(f,'t');其中x为给定点横坐标数组，y为给定点纵坐标数组，y_1为原函数在给定点处的导数数组。
测试证明该程序可以实现，例如输入如下数组：x=1::;y_1=[ ];y=[1 ];>> [f,f0] = Hermite1(x,y,y_1)；运行结果如下：f =$(390625*((3972231*t)/35 - 28321/0000)*(t - 1)^2*(t - 7/5)^2*(t - 8/5)^2*(t - 9/5)^2)/36 - (390625*(t - 1)^2*(t - 6/5)^2*(t - 7/5)^2*(t - 9/5)^2*((28557*t)/28 - 23/2000))/36 + (390625*((64*t)/3 - 61/3)*(t - 6/5)^2*(t - 7/5)^2*(t - 8/5)^2*(t - 9/5)^2)/576 + (390625*((763*t)/1984 + 043/6240000)*(t - 1)^2*(t - 6/5)^2*(t - 8/5)^2*(t - 9/5)^2)/16 - (390625*((77623*t)/28 - 931/60000)*(t - 1)^2*(t - 6/5)^2*(t - 7/5)^2*(t - 8/5)^2)/576f0 =．利用matlab绘制图像：2、程序的窗口化：利用Matlab提供的GUIDE工具以及callback函数实现相应函数的窗口化，GUI代码如下：function varargout = untitled(varargin)?% UNTITLED M-file for% UNTITLED, by itself, creates a new UNTITLED or raises the existing% singleton*.%% H = UNTITLED returns the handle to a new UNTITLED or the handle to% the existing singleton*.%% UNTITLED('CALLBACK',hObject,eventData,handles,...) calls the local% function named CALLBACK in with the given input arguments.%% UNTITLED('Property','Value',...) creates a new UNTITLED or raises the，% existing singleton*. Starting from the left, property value pairs are% applied to the GUI before untitled_OpeningFcn gets called. An% unrecognized property name or invalid value makes property application% stop. All inputs are passed to untitled_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 untitled%% Last Modified by GUIDE 15-Sep-2011 22:24:48% Begin initialization code - DO NOT EDITgui_Singleton = 1;gui_State = struct('gui_Name', mfilename, ...'gui_Singleton', gui_Singleton, ...'gui_OpeningFcn', @untitled_OpeningFcn, ...'gui_OutputFcn', @untitled_OutputFcn, ...'gui_LayoutFcn', [] , ...'gui_Callback', []);】if nargin && ischar(varargin{1})= str2func(varargin{1});endif nargout[varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:}); elsegui_mainfcn(gui_State, varargin{:});end% End initialization code - DO NOT EDIT<% --- Executes just before untitled is made visible.function untitled_OpeningFcn(hObject, eventdata, handles, varargin) % This function has no output args, see OutputFcn.% hObject handle to figure% eventdata reserved - to be defined in a future version of MATLAB % handles structure with handles and user data (see GUIDATA)% varargin command line arguments to untitled (see VARARGIN)% Choose default command line output for untitled= hObject;>% Update handles structureguidata(hObject, handles);% UIWAIT makes untitled wait for user response (see UIRESUME)% uiwait;% --- Outputs from this function are returned to the command line.function varargout = untitled_OutputFcn(hObject, eventdata, handles)% varargout cell array for returning output args (see VARARGOUT);…% hObject handle to figure% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA)% Get default command line output from handles structurevarargout{1} = ;function edit1_Callback(hObject, eventdata, handles)% hObject handle to edit1 (see GCBO)…% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA)% Hints: get(hObject,'String') returns contents of edit1 as text% str2double(get(hObject,'String')) returns contents of edit1 as a double guidata(hObject, handles);% --- Executes during object creation, after setting all properties.function edit1_CreateFcn(hObject, eventdata, handles)% hObject handle to edit1 (see GCBO)（% eventdata reserved - to be defined in a future version of MATLAB% handles empty - handles not created until after all CreateFcns called% Hint: edit controls usually have a white background on Windows.% See ISPC and COMPUTER.if ispc && isequal(get(hObject,'BackgroundColor'),get(0,'defaultUicontrolBackgroundColor'))set(hObject,'BackgroundColor','white');end￥function edit2_Callback(hObject, eventdata, handles)% hObject handle to edit2 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA)% Hints: get(hObject,'String') returns contents of edit2 as text% str2double(get(hObject,'String')) returns contents of edit2 as a double guidata(hObject, handles);% --- Executes during object creation, after setting all properties.、function edit2_CreateFcn(hObject, eventdata, handles)% hObject handle to edit2 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles empty - handles not created until after all CreateFcns called% Hint: edit controls usually have a white background on Windows.% See ISPC and COMPUTER.if ispc && isequal(get(hObject,'BackgroundColor'),get(0,'defaultUicontrolBackgroundColor'))set(hObject,'BackgroundColor','white');end—function edit3_Callback(hObject, eventdata, handles)% hObject handle to edit3 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA)% Hints: get(hObject,'String') returns contents of edit3 as text% str2double(get(hObject,'String')) returns contents of edit3 as a double guidata(hObject, handles);·% --- Executes during object creation, after setting all properties.function edit3_CreateFcn(hObject, eventdata, handles)% hObject handle to edit3 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles empty - handles not created until after all CreateFcns called% Hint: edit controls usually have a white background on Windows.% See ISPC and COMPUTER.if ispc && isequal(get(hObject,'BackgroundColor'),get(0,'defaultUicontrolBackgroundColor'))set(hObject,'BackgroundColor','white');·endfunction edit4_Callback(hObject, eventdata, handles)% hObject handle to edit4 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA)% Hints: get(hObject,'String') returns contents of edit4 as text% str2double(get(hObject,'String')) returns contents of edit4 as a double '% --- Executes during object creation, after setting all properties.function edit4_CreateFcn(hObject, eventdata, handles)% hObject handle to edit4 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles empty - handles not created until after all CreateFcns called% Hint: edit controls usually have a white background on Windows.% See ISPC and COMPUTER.if ispc && isequal(get(hObject,'BackgroundColor'),get(0,'defaultUicontrolBackgroundColor'))；set(hObject,'BackgroundColor','white');end% --- Executes on button press in pushbutton1.function pushbutton1_Callback(hObject, eventdata, handles)% hObject handle to pushbutton1 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA)x=str2num(get,'string'));y=str2num(get,'string'));<y_1=str2num(get,'string'));x0=str2num(get,'string'));syms t;f = ;if(length(x) == length(y))if(length(y) == length(y_1))n = length(x);elsedisp('yºÍyµÄµ¼ÊýµÄάÊý²»ÏàµÈ');return;end—elsedisp('xºÍyµÄάÊý²»ÏàµÈ£¡ ');return;endfor i=1:nh = ;a = ;for j=1:nif( j ~= i)h = h*(t-x(j))^2/((x(i)-x(j))^2);a = a + 1/(x(i)-x(j));、endendf = f + h*((x(i)-t)*(2*a*y(i)-y_1(i))+y(i));endf0 = subs(f,'t',x0);plot,x,y,'*');^function edit5_Callback(hObject, eventdata, handles)% hObject handle to edit5 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA)% Hints: get(hObject,'String') returns contents of edit5 as text% str2double(get(hObject,'String')) returns contents of edit5 as a doubleguidata(hObject, handles);{% --- Executes during object creation, after setting all properties.function edit5_CreateFcn(hObject, eventdata, handles)% hObject handle to edit5 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles empty - handles not created until after all CreateFcns called% Hint: edit controls usually have a white background on Windows.% See ISPC and COMPUTER.if ispc && isequal(get(hObject,'BackgroundColor'),get(0,'defaultUicontrolBackgroundColor'))~set(hObject,'BackgroundColor','white');end程序运行结果：其中左上方纵列的三个对话框从上到下分别输入给定点的横坐标x，纵坐标y以及导数值y_1，右侧空白框输入维数，下方坐标图显示插值函数图像，例如仍插入上面所给定的点列，得出结果：从图上看拟合程度还是比较不错的。
分段三次 hermite 插值多项式的数学表达
分段三次 hermite 插值多项式的数学表达Hermite插值多项式是一种用于在给定的点集上进行插值的数学工具。
与其他插值方法不同的是，Hermite插值多项式不仅考虑了函数在各个插值点上的函数值，还考虑了函数在该点上的导数值。
这使得Hermite插值多项式能够更准确地拟合函数的曲线特征，特别是在存在函数奇点或不连续点的情况下。
要理解Hermite插值多项式的具体数学表达，首先需要了解插值点和插值条件的概念。
设给定的插值点集为{(x0, f0, f'0), (x1, f1, f'1), ..., (xn, fn, f'n)}，其中xi为插值点的横坐标，fi为插值点的纵坐标，f'i为插值点的导数值。
我们的目标是构造一个多项式P(x)，满足以下条件：1.在每个插值点(xi, fi)处，多项式P(x)的函数值等于fi：P(xi) = fi；2.在每个插值点(xi, fi)处，多项式P(x)的导数值等于f'i：P'(xi) = f'i。
根据这些插值条件，我们可以得到Hermite插值多项式的数学表达式。
首先，我们需要定义一个Lagrange插值基函数Lk(x)，用于描述在插值点xi处的多项式P(x)的函数值。
Lagrange插值基函数可以通过以下公式计算得到：Lk(x) = Π(j ≠ k) [(x - xj) / (xk - xj)]其中Π是乘积符号，j和k分别表示插值点的索引。
然后，我们可以构造Hermite插值多项式Hk(x)，它的数学表达式可以通过Lagrange插值基函数和插值点的函数值、导数值得到：Hk(x) = [1 - 2(x - xi)L'i(xi)]Li(x)^2 + (x - xi)Li(x)^2 其中Li(x)表示第i个插值点处的Lagrange插值基函数，L'i(xi)表示第i个插值点处的Lagrange插值基函数的导数值。
数值分析各算法流程图
01,,n1,,n1,,)n x及数值分析各算法流程图一、插值1、 拉格朗日插值流程图：( 相应程序：lagrintp(x,y,xx))2,,n ,,j n 1,2,,n 1,,)n 2、 牛顿插值流程图(1)产生差商表的算法流程图(相应程序：divdiff(x,y))注：1、另一程序divdiff1(x,y)，输出的矩阵包含了节点向量。
而divdiff(x,y)不含节点向量。
2、另一程序tableofdd(x,y,m)，输出的是表格形式，添加了表头。
1,,),,n m 及1,,m (2)非等距节点的牛顿插值流程图(相应程序：newtint11(x,y,xx,m)) 、注：1、虽然程序newtint11(x,y,xx,m)考虑了多种情形，看上去很复杂，但基本流程结构还是如上图所示。
2、程序中调用的子程序是divdiff 。
若调用的子程序是divdiff1的话，流程图中的第三，第四，第五步要相应的改一下数字。
2,3,,1m +1,,j1,2,,n=1,2,,)n m 及(3)求差分表的流程图(相应程序：difference(y,m))注：1、difference 输出的是矩阵D 。
而另一程序tableofd(y,m)，输出的是带有表头的差分表。
n x m1,,),,1,,m注：1、程序newtforward1(x,y,xx,m))的结构与上述流程图一致，xx可以是数组。
2、另一程序newtforward(x,y,xx,m))先求出插值多项式，再求插值多项式在插值点的函数值。
基本结构还是和上面的流程图一样。
n x m1,,),,-x x1,,m注：1、程序newtbackward1(x,y,xx,m))的结构与上述流程图一致，xx可以是数组。
2、另一程序newtbackward(x,y,xx,m))先求出插值多项式，再求插值多项式在插值点的函数值。
基本结构还是和上面的流程图一样。
1,2,,n1,2,,n ,2,,)n x及3、Hermite 插值流程图(1) 已知条件中一阶导数的个数与插值节点的个数相等时的Hermite 插值流程图。
hermite插值以及两种MATLAB程序
给定矢量P 0, P 1, R 0, R 1,称满足下列条件的参数三次多项式曲线P (t ),t ∈[0,1]为Hermite 曲线：H (x 0)=y 0,H (x 1)=y 1,H ′(x 0)=m 0,H ′(x 1)=m 1,即Hermite 曲线两个端点为P 0, P 1,在两端点的切矢量分别R 0, R 1。
记几何矩阵和基矩阵分别为G H , M H , G H , M H 是未知的.取G H =[P 0,P 1,R 0,R 1],则只要M H 就可以了。
一般的曲线经过多项式分解, 得到参数多项式曲线的矩阵表示:P (t )=G ∙M ∙T将（1）式代入（2）得到：G H ∙M H ∙T H |t=0=G H ∙M H ∙(1,0,0,0)T =P 0, G H ∙M H ∙T H |t=1=G H ∙M H ∙(1,1,1,1)T =P 1, G H ∙M H ∙T H |t=0=G H ∙M H ∙(0,1,0,0)T =R 0, G H ∙M H ∙T H |t=0=G H ∙M H ∙(0,1,2,3)T =R 1, 将上面四个式子合并如下形式：G H ∙M H ∙[101001111020103]=[P 0,P 1,R 0,R 1]=G H 上面方程的解不唯一，不妨取M H =[101001111020103]−1=[1000−320−3−21−2100−11] 从而得到三次Hermite 曲线的方程：P(t)=G H∙M H∙T其中M H∙T确定了一组Hermite基函数G0(t)，G1(t)，H0(t)，H1(t)，即M H∙T=[10−320−3−21−2100−11][1tt2t3]=[1−3t2+2t33t2−2t3t−2t2+t3−t2+t3]附：MATLAB程序function yy=hermite(x,y,dy,xx)% 输入X——左右两个端点的X轴坐标Y——左右两个端点的Y轴坐标dy——左右两个端点的切矢xx——中间插值的点X轴坐标%输出yy——中间插值的点Y轴坐标function yy=hermite(x,y,dy,xx)k=length(xx);z=zeros(1,k);for i=1:k;s=0;xaix=xx(i);a=1-3.*(xaix)^2+2.*(xaix)^3;b=2.*(xaix)^2-2.*(xaix)^3;c=xaix-2.*(xaix)^2+(xaix)^3;d=-2.*(xaix)^2+(xaix)^3;s=y(1)*a+y(2)*b+dy(1)*c+dy(2)*d;z(i)=s;endyy=z;function yy=hermite(x,y,dy,xx) % 输入X——左右两个端点的X轴坐标Y——左右两个端点的Y轴坐标dy——左右两个端点的切矢xx——中间插值的点X轴坐标%输出yy——中间插值的点Y轴坐标m=length(x);n=length(y);l=length(dy);k=length(xx);if m~=n,error('向量长度不一样'); end;if n~=l,error('向量长度不一样'); end;z=zeros(1,k);for i=1:k;s=0;a=xx(i)-x(1);b=x(1)-x(2);c=xx(i)-x(2);a1=(1-2*a/b)*(c/b)^2;aa=xx(i)-x(2);a2=(1+2*aa/b)*(a/b)^2;b1=a*(c/b)^2;b2=c*(a/b)^2;s=y(1)*a1+y(2)*a2+dy(1)*b1+dy(2)*b2; z(i)=s;endyy=z;。