分段线性插值函数的编程实现
用MATLAB实现拉格朗日插值和分段线性插值

用MATLAB真止推格朗日插值战分段线性插值之阳早格格创做1、真验真质:用MATLAB真止推格朗日插值战分段线性插值.2、真验手段:1)教会使用MATLAB硬件;2)会使用MATLAB硬件举止推格朗日插值算法战分段线性好值算法;3、真验本理:利用推格朗日插值要领举止多项式插值,并将图形隐式出去.4、真验步调及运止截止(1)真止lagrange插值1)定义函数:f = 1/(x^2+1) 将其保存正在f.m 文献中,简直步调如下:function y = f1(x)y = 1./(x.^2+1);2)定义推格朗日插值函数:将其保存正在lagrange.m 文献中,简直真止步调编程如下:function y = lagrange(x0,y0,x)m = length(x); /区间少度/n = length(x0);for i = 1:nl(i) = 1;endfor i = 1:mfor j = 1:nfor k = 1:nif j == kcontinue;endl(j) = ( x(i) -x0(k))/( x0(j) - x0(k) )*l(j);endendendy = 0;for i = 1:ny = y0(i) * l(i) + y;end3)修坐尝试步调,保存正在text.m文献中,真止绘图:x=-5:0.001:5;y=(1+x.^2).^-1;p=polyfit(x,y,n);py=vpa(poly2sym(p),10)plot_x=-5:0.001:5;f1=polyval(p,plot_x);figureplot(x,y,‘r',plot_x,f1)输进n=6,出现底下的图形:通过上图不妨瞅到当n=6是不很佳的模拟.于是沉新运止text.M并采用n=11由此可睹n=11时的图像是不妨很佳的真止模拟(2)分段线性插值:修坐div_linear.m文献.简直编程如下/*分段线性插值函数:div_linear.m 文献*/function y = div_linear(x0,y0,x,n)%for j = 1:length(x)for i = 1:n-1if (x >= x0(i)) && (x <= x0(i+1))y = (x - x0(i+1))/(x0(i) - x0(i+1))*y0(i) + ( x - x0(i))/(x0(i+1) - x0(i))*y0(i+1);elsecontinue;endend%end尝试步调(text2.m):n = input(‘输进n =:’);x0 = linspace( -5,5,n);for x = -5:0.01:5y = div_linear(x0,f(x0),x,n);hold on;plot(x,y,'r');plot(x,f(x),'b');end2)运止尝试步调,那是会出现:输进n=:2)输进n=6,并按Enter键,出现:4)闭掉图形界里后,沉新运止步调,输进n=11,并按enter键后出现:5)再次闭掉图形界里,输进n=100,并按enter键,出现:此时.图形将于本函数图形基础符合,证明分隔区间越多,图像交近真正在的图像.(3)用lagrange插值瞅察y = |si n(k*π*x)|的缺点分解:1)编写函数文献,保存正在f2.m 中x=0:0.01:1;k= input('输进k:')n= input('输进n:');y=abs(sin(k*pi*x));p=polyfit(x,y,n-1);py=vpa(poly2sym(p),8);plot_x=0:0.01:1;f1=polyval(p,plot_x);plot(x,y,plot_x,f1);2)运止该步调:输进k=:1输进n=:2出现如下图形界里:闭掉图形界里后沉新运止f2.m,输进k=:1,n=:3出现如下界里:再次闭掉图形界里,输进k=:1,n=:6 后出现:此时图形基础符合.类推,输进k=2,n=3后出现:k =2, n =11,出现如下图形:k =2,n =15,那时出现:k =2,n =19,出现:当k=2,n=21时,图形如下:此时基础符合.5、真验归纳:通过本次课程安排,尔发端掌握了MATLAB使用,加深了对付于百般线性插值的明白;培植了独力处事本领战创制力;概括使用博业及前提知识,办理本质数教问题的本领;正在本次课程安排中,正在教授的粗心指挥下,支益匪浅.共时对付数教的钻研有了更深进的认识.。
插值法——线性分段插值

插值法——线性分段插值 1.插值函数%%分段线性插值function PLI = Piecewise_linear_interpolation(X,f,precision)[m,n] = size(X);a = min(X);b = max(X);X = sort(X);F = subs(f,X);for k = 1:n-1B = Basic_fun(X,k);I = B(1)*F(k)+B(2)*F(k+1);PLI{1,k} = [X(k),X(k+1)];PLI{2,k} = I;t{k} = X(k):(X(k+1)-X(k))/precision:X(k+1);T{k} = subs(I,t{k});Y_real{k} = subs(f,t{k});endfor k = 1:n-1t_((precision+1)*(k-1)+1:(precision+1)*k) = t{k};T_((precision+1)*(k-1)+1:(precision+1)*k) = T{k};Y_real_((precision+1)*(k-1)+1:(precision+1)*k)= Y_real{k};endh = figure;set(h,'color','w');plot(X,F,'r*',t_,T_,'g',t_,Y_real_,'b');xlabel('x shaft');ylabel('y shaft');legend('F:节点对应函数值','T:分段线性插值函数图像','Y_real:真实函数图像');title('分段线性插值');grid onend 2.基函数%%基函数,max(X)>k>0function BF = Basic_fun(X,k)X = sort(X);syms x;BF(1) = (x-X(k+1))/(X(k)-X(k+1));BF(2) = (x-X(k))/(X(k+1)-X(k));end 3.拟合值函数%%线性插值拟合值function LIV = Linear_interpolation_value(X,f,precision,x_value)[m,n] = size(X);a = min(X);b = max(X);X = sort(X);Answer = Piecewise_linear_interpolation(X,f,precision);for i = 1:n-1if x_value >= X(i) && x_value <= X(i+1)s = i;endendLIV{1,1} = '线性插值拟合值';LIV{2,1} = vpa(subs(Answer{2,s},x_value),6);LIV{1,2} = '真实值';LIV{2,2} = vpa(subs(f,x_value),6);LIV{1,3} = '误差';LIV{2,3} = abs(LIV{2,1}-LIV{2,2});end 4.例⼦clear allclcX = -5:1:5;syms x;f = - 0.08858*x^8 + 3.694*x^7 - 64.7*x^6 + 617.8*x^5 - 3490.0*x^4 + 11820.0*x^3 - 23150.0*x^2 + 23580.0*x - 9319.0; precision = 200;%%分段线性插值disp('分段线性插值');Piecewise_linear_interpolation(X,f,precision) 结果分段线性插值S =2×10 cell 数组列 1 ⾄ 4{1×2 double} {1×2 double} {1×2 double} {1×2 double} {1×1 sym } {1×1 sym } {1×1 sym } {1×1 sym }列 5 ⾄ 8{1×2 double} {1×2 double} {1×2 double} {1×2 double} {1×1 sym } {1×1 sym } {1×1 sym } {1×1 sym }列 9 ⾄ 10{1×2 double} {1×2 double}{1×1 sym } {1×1 sym }>> S{2,:}ans =(227077586881*x)/50000 + 37695704689/2500ans =(3983468847*x)/2000 + 60987657739/12500ans =(7723057429*x)/10000 + 30518164433/25000ans =(2518396259*x)/10000 + 4494858583/25000ans =(3136314129*x)/50000 - 9319ans =(465835271*x)/50000 - 9319ans =(422501*x)/10000 - 1113617/25000ans =4111433/25000 - (622509*x)/10000ans =- (271*x)/80 - 151661/12500ans =2072089/2500 - (10681481*x)/50000 图像如下。
lagrange插值分段线性插值matlab代码

Lagrange插值:x=0:3;y=[-5,-6,-1,16];n=length(x);syms q;for k=1:nfenmu=1;p=1;for j=1:nif(j~=k)fenmu=fenmu*(x(k)-x(j))p=conv(p,poly(x(j)))endendc(k,:)=p*y(k)/fenmuenda=zeros(1,n);for i=1:nfor j=1:na(i)=a(i)+c(j,i)endend输出结果:fenmu =-1p =1 -1fenmu =2p =1 -3 2fenmu =-6p =1 -6 11 -6c =0.8333 -5.0000 9.1667 -5.0000 fenmu =1p =1 0fenmu =-1p =1 -2 0fenmu =2p =1 -5 6 0c =0.8333 -5.0000 9.1667 -5.0000-3.0000 15.0000 -18.0000 0 fenmu =2p =1 0fenmu =2p =1 -1 0fenmu =-2p =1 -4 3 0c =0.8333 -5.0000 9.1667 -5.0000-3.0000 15.0000 -18.0000 00.5000 -2.0000 1.5000 0 fenmu =3p =1 0fenmu =6p =1 -1 0fenmu =6p =1 -32 0c =0.8333 -5.0000 9.1667 -5.0000-3.0000 15.0000 -18.0000 00.5000 -2.0000 1.5000 02.6667 -8.0000 5.3333 0a =0.8333 0 0 0a =-2.1667 0 0 0 a =-1.6667 0 0 0a =1 0 0 0a =1 -5 0 0a =1 10 0 0a =1 8 0 0a =1 0 0 0a =1.0000 0 9.1667 0a =1.0000 0 -8.8333 0a =1.0000 0 -7.3333 0a =1.0000 0 -2.0000 0a =1.0000 0 -2.0000 -5.0000a =1.0000 0 -2.0000 -5.0000a =1.0000 0 -2.0000 -5.0000a =1.0000 0 -2.0000 -5.0000 分段线性插值:先保存M文件:x=1:6;y=[7 16 8 25 12 24];u=5.3;delta=diff(y)./diff(x);n=length(x);for j=2:(n-1)if x(j)<uk=j;endend在command window中输入:s=u-x(k);v=y(k)+s.*delta(k)输出结果:v =15.6000解:第一种做法,用spline,共55个点,其中,54个有效首先保存你一个M文件:figure('position',get(0,'screensize'))axes('position',[0 0 1 1])[x,y] = ginput;然后在command window里输入以下内容:n = length(x);s = (1:n)';t = (1:.05:n)';u = spline(s,x,t);v = spline(s,y,t);clf resetplot(x,y,'.',u,v,'-');对应的x、y值:0.3572917 0.25361450.3572917 0.29096390.3503472 0.34036140.3461806 0.42590360.3427083 0.52710840.3253472 0.61626510.3065972 0.68734940.290625 0.75240960.2892361 0.79337350.2954861 0.7969880.3225694 0.75481930.340625 0.68493980.3690972 0.6150602 0.3864583 0.6126506 0.3899306 0.7259036 0.3927083 0.8066265 0.3920139 0.8993976 0.4024306 0.9295181 0.4239583 0.8933735 0.4239583 0.8078313 0.4295139 0.7343373 0.4315972 0.6451807 0.4440972 0.6439759 0.4565972 0.7439759 0.4704861 0.8451807 0.4767361 0.9054217 0.4961806 0.9463855 0.5086806 0.876506 0.5045139 0.8186747 0.5010417 0.7524096 0.4892361 0.6403614 0.503125 0.6295181 0.5052083 0.6271084 0.5322917 0.7090361 0.5510417 0.763253 0.5739583 0.8355422 0.5961806 0.8572289 0.5947917 0.7837349 0.5753472 0.7090361 0.5579861 0.6391566 0.5357639 0.5668675 0.5322917 0.5283133 0.5350694 0.4789157 0.565625 0.536747 0.5947917 0.5933735 0.6253472 0.610241 0.6322917 0.5728916 0.615625 0.5331325 0.6003472 0.4993976 0.5788194 0.4415663 0.559375 0.3716867 0.5295139 0.2957831 0.4975694 0.2403614 0.4711806 0.2018072 0.6607639 0.3090361第二种做法,用pchip,共52个点,全部有效首先保存一个M文件:figure('position',get(0,'screensize'))axes('position',[0 0 1 1])[x,y] = ginput;然后在command window里输入以下内容:n = length(x);s = (1:n)';t = (1:.05:n)';u = pchip (s,x,t);v = pchip (s,y,t);clf resetplot(x,y,'.',u,v,'-');对应的x、y值:0.5190972 0.84879520.5052083 0.75120480.4947917 0.67891570.5100694 0.66927710.5399306 0.73554220.5753472 0.81746990.596875 0.86204820.6190972 0.87771080.6149306 0.81385540.5878472 0.74277110.5878472 0.74277110.5635417 0.67168670.5350694 0.6030120.528125 0.5632530.528125 0.52590360.565625 0.58012050.6052083 0.62710840.634375 0.61867470.6190972 0.57168670.5878472 0.5234940.5364583 0.41265060.4961806 0.32108430.459375 0.2753012我更喜欢第一种,用spline的,这个能将之间画出弧度,而pchip更像是直接用线段将点依次连接得到的。
线性分段插值函数

数值分析实验报告任课教师:马季骕班级:11级计算机科学与技术1实验目的及要求2程序的源代码3实验操作4实验结果及分析1实验目的及要求学会使用分段线性插值的方法求原函数的逼近函数。
随着插值节点的增加,插值多项式的插值多项式的次数也增加,而对于高次的插值容易带来剧烈的震荡,带来数值的不稳定(Runge现象)。
为了既要增加插值的节点,减小插值的区间,以便更好的逼近插值函数,又要不增加插值多项式的次数以减少误差,可采用分段线性插值.当给出了n+1个节点上f(x)的一张函数表后,用分段线性插值法求一个函数φ(x),并满足φ(x)是一个不超过n次的多项式;当插值节点取的足够多时逼近函数φ(x)能够很好的逼近被逼近函数f(x)。
而插值函数φ(x)的次数就会相应地升高,高次的插值多项式收敛到相应的被逼近函数。
所取的n值越大,逼近的越准确。
2程序源代码3// 数值分析Dlg.cpp : implementation file4//56#include "stdafx.h"7#include "数值分析.h"8#include "数值分析Dlg.h"910#ifdef _DEBUG11#define new DEBUG_NEW12#undef THIS_FILE13static char THIS_FILE[] = __FILE__;14#endif1516/////////////////////////////////////////////////////////////////////////// //17// CAboutDlg dialog used for App About1819class CAboutDlg : public CDialog20{21public:22CAboutDlg();2324// Dialog Data25//{{AFX_DATA(CAboutDlg)26enum { IDD = IDD_ABOUTBOX };27//}}AFX_DATA2829// ClassWizard generated virtual function overrides30//{{AFX_VIRTUAL(CAboutDlg)31protected:32virtual void DoDataExchange(CDataExchange* pDX); // DDX/DDV support33//}}AFX_VIRTUAL3435// Implementation36protected:37//{{AFX_MSG(CAboutDlg)38//}}AFX_MSG39DECLARE_MESSAGE_MAP()40};4142CAboutDlg::CAboutDlg() : CDialog(CAboutDlg::IDD)43{44//{{AFX_DATA_INIT(CAboutDlg)45//}}AFX_DATA_INIT46}4748void CAboutDlg::DoDataExchange(CDataExchange* pDX)49{50CDialog::DoDataExchange(pDX);51//{{AFX_DATA_MAP(CAboutDlg)52//}}AFX_DATA_MAP53}5455BEGIN_MESSAGE_MAP(CAboutDlg, CDialog)56//{{AFX_MSG_MAP(CAboutDlg)57// No message handlers58//}}AFX_MSG_MAP59END_MESSAGE_MAP()6061/////////////////////////////////////////////////////////////////////////////62// CMyDlg dialog6364CMyDlg::CMyDlg(CWnd* pParent /*=NULL*/)65: CDialog(CMyDlg::IDD, pParent)66{67//{{AFX_DATA_INIT(CMyDlg)68// NOTE: the ClassWizard will add member initialization here69//}}AFX_DATA_INIT70// Note that LoadIcon does not require a subsequent DestroyIcon in Win32 71m_hIcon = AfxGetApp()->LoadIcon(IDR_MAINFRAME);72}7374void CMyDlg::DoDataExchange(CDataExchange* pDX)75{76CDialog::DoDataExchange(pDX);77//{{AFX_DATA_MAP(CMyDlg)78// NOTE: the ClassWizard will add DDX and DDV calls here79//}}AFX_DATA_MAP80}8182BEGIN_MESSAGE_MAP(CMyDlg, CDialog)83//{{AFX_MSG_MAP(CMyDlg)84ON_WM_SYSCOMMAND()85ON_WM_PAINT()86ON_WM_QUERYDRAGICON()87ON_BN_CLICKED(IDC_NEW, OnNew)88ON_BN_CLICKED(IDC_LAG, OnLag)89ON_BN_CLICKED(IDC_LINE, OnLine)90ON_BN_CLICKED(IDC_HER, OnHer)91//}}AFX_MSG_MAP92END_MESSAGE_MAP()9394/////////////////////////////////////////////////////////////////////////// //95// CMyDlg message handlers9697BOOL CMyDlg::OnInitDialog()98{99CDialog::OnInitDialog();100101// Add "About..." menu item to system menu.102103// IDM_ABOUTBOX must be in the system command range.104ASSERT((IDM_ABOUTBOX & 0xFFF0) == IDM_ABOUTBOX);105ASSERT(IDM_ABOUTBOX < 0xF000);106107CMenu* pSysMenu = GetSystemMenu(FALSE);108if (pSysMenu != NULL)109{110CString strAboutMenu;111strAboutMenu.LoadString(IDS_ABOUTBOX);112if (!strAboutMenu.IsEmpty())113{114pSysMenu->AppendMenu(MF_SEPARATOR);115pSysMenu->AppendMenu(MF_STRING, IDM_ABOUTBOX, strAboutMenu);116}117}118119// Set the icon for this dialog. The framework does this automatically120// when the application's main window is not a dialog121SetIcon(m_hIcon, TRUE); // Set big icon122SetIcon(m_hIcon, FALSE); // Set small icon123124// TODO: Add extra initialization here125126return TRUE; // return TRUE unless you set the focus to a control127}128129void CMyDlg::OnSysCommand(UINT nID, LPARAM lParam)130{131if ((nID & 0xFFF0) == IDM_ABOUTBOX)132{133CAboutDlg dlgAbout;134dlgAbout.DoModal();135}136else137{138CDialog::OnSysCommand(nID, lParam);139}140}141142// If you add a minimize button to your dialog, you will need the code below 143// to draw the icon. For MFC applications using the document/view model, 144// this is automatically done for you by the framework.145146void CMyDlg::OnPaint()147{148if (IsIconic())149{150CPaintDC dc(this); // device context for painting151152SendMessage(WM_ICONERASEBKGND, (WPARAM) dc.GetSafeHdc(), 0);153154// Center icon in client rectangle155int cxIcon = GetSystemMetrics(SM_CXICON);156int cyIcon = GetSystemMetrics(SM_CYICON);157CRect rect;158GetClientRect(&rect);159int x = (rect.Width() - cxIcon + 1) / 2;160int y = (rect.Height() - cyIcon + 1) / 2;161162// Draw the icon163dc.DrawIcon(x, y, m_hIcon);164}165else166{167CDialog::OnPaint();168}169}170171// The system calls this to obtain the cursor to display while the user drags 172// the minimized window.173HCURSOR CMyDlg::OnQueryDragIcon()174{175return (HCURSOR) m_hIcon;176}177178void CMyDlg::OnNew()179{180// TODO: Add your control notification handler code here181 int x00=225,y00=225,i,j;182double x;183184CDC *pDC=GetDC();185pDC->SetMapMode(MM_LOMETRIC);186187pDC->SetViewportOrg(x00,y00);188189//画坐标轴与原函数190for(i=-650; i<=650; i++)191{192pDC->SetPixel(i,0,RGB(0,0,0));193pDC->SetPixel(0,i,RGB(0,0,0));194}195196for(x=-1; x<=1; x+=0.001)197{198double j=400.0/(1+25*x*x);199pDC->SetPixel(x*500,j,RGB(255,0,0));200}201202}203204void CMyDlg::OnLag()205{206// TODO: Add your control notification handler code here 207int x00=225,y00=225,i,j;208double x;209210CDC *pDC=GetDC();211pDC->SetMapMode(MM_LOMETRIC);212213pDC->SetViewportOrg(x00,y00);214215//画坐标轴216for(i=-650; i<=650; i++)217{218pDC->SetPixel(i,0,RGB(0,0,0));219pDC->SetPixel(0,i,RGB(0,0,0));220}221222223double yx[]={-1,-0.8,-0.6,-0.4,-0.2,0,0.2,0.4,0.6,0.8,1}; 224225// 拉格朗日差值的函数226double yy[12],lx[12],ly[12];227double l_fenzi[12],l_fenmu[12];228 double l_x,l_y;229for(i=0; i<=10; i++)230{231yy[i]=1.0/(1+25*yx[i]*yx[i]);232}233234for(i=0; i<=10; i++)235{236l_fenmu[i]=1.0;237for(j=0; j<=10; j++)238{239if(i!=j)240l_fenmu[i]=l_fenmu[i]*(yx[i]-yx[j]);241}242}243 double qq,pp;244for(qq=-1; qq<=1; qq+=0.0003)245{246for(i=0; i<=10; i++)247{248l_fenzi[i]=1.0;249for(j=0; j<=10; j++)250{251if(i!=j)252 l_fenzi[i]=l_fenzi[i]*(qq-yx[j]);253}254}255pp=0;256for(i=0; i<=11; i++)257{258 pp=pp+1.0/(1+25*yx[i]*yx[i])*l_fenzi[i]/l_fenmu[i]; 259}260pDC->SetPixel(qq*500,pp*390+5,RGB(150,255,0));261}262263}264265void CMyDlg::OnLine()266{267// TODO: Add your control notification handler code here 268int x00=225,y00=225,i,j;269double x;270271CDC *pDC=GetDC();272pDC->SetMapMode(MM_LOMETRIC);273274pDC->SetViewportOrg(x00,y00);275276//画坐标轴与原函数277for(i=-650; i<=650; i++)278{279pDC->SetPixel(i,0,RGB(0,0,0));280pDC->SetPixel(0,i,RGB(0,0,0));281}282283double yx[]={-1,-0.8,-0.6,-0.4,-0.2,0,0.2,0.4,0.6,0.8,1}; 284double yy[14];285for(i=0; i<=10; i++)286{287yy[i]=1.0/(1+25*yx[i]*yx[i]);288}289290// 线性分段差值的图像291CPen pen;292293CPen*oldpen;294295pen.CreatePen(PS_SOLID,5,RGB(0,0,0));296297oldpen=pDC->SelectObject(&pen);298for(i=0; i<10; i++)299{300 pDC->MoveTo(yx[i]*480,yy[i]*400);301pDC->LineTo(yx[i+1]*480,yy[i+1]*400);302 }303}304305306void CMyDlg::OnHer()307{308// TODO: Add your control notification handler code here 309int x00=225,y00=225,i,j;310double x;311312CDC *pDC=GetDC();313pDC->SetMapMode(MM_LOMETRIC);314315pDC->SetViewportOrg(x00,y00);316317//画坐标轴与原函数318for(i=-650; i<=650; i++)319{320pDC->SetPixel(i,0,RGB(0,0,0));321pDC->SetPixel(0,i,RGB(0,0,0));322}323324 double yx[]={-1,-0.8,-0.6,-0.4,-0.2,0,0.2,0.4,0.6,0.8,1};325326 double yy[12];327 for(i=0; i<=10; i++)328 {329 yy[i]=1.0/(1+25*yx[i]*yx[i]);330 }331332 //分段三次Hermite 差值的函数333334 double x0,x1,yd1,yd0,y1,y0;335 for(i=0; i<10; i++)336 {337338 x0=yx[i],x1=yx[i+1];339 y0=1.0/(1+25*x0*x0);340 y1=1.0/(1+25*x1*x1);341 yd0=-(50*x0)*1.0/((1+25*x0*x0)*(1+25*x0*x0));342 yd1=-(50*x1)*1.0/((1+25*x1*x1)*(1+25*x1*x1));343344 for(double qq=x0; qq<x1; qq+=0.005)345 {346 double pp= y0*(1+2*(qq-x0)/(x1-x0)) * (qq-x1)/(x0-x1) *(qq-x1)/(x0-x1)347 +y1*(1+2*(qq-x1)/(x0-x1)) * (qq-x0)/(x1-x0) * (qq-x0)/(x1-x0) 348 +yd0*(qq-x0) * (qq-x1)/(x0-x1) * (qq-x1)/(x0-x1)349 +yd1*(qq-x1) * (qq-x0)/(x1-x0) * (qq-x0)/(x1-x0);350 pDC->SetPixel(qq*500,pp*400,RGB(0,255,255));351 }352 }353354 }程序操作算法分析:1. 分段线性插值的算法思想:分段线性插值需要在每个插值节点上构造分段线性插值基函数)(x l j ,然后再作它们的线性组合。
分段插值法

// 分段插值法.cpp : 定义控制台应用程序的入口点。
//分段线性插值和分段抛物差值 ;#include"stdafx.h"#include<iostream>#include<iomanip>#include<math.h>using namespace std;double Paixu (double a[][40],int m){double t;int k;for (int i=0;i<m-1;i++){k=i;for (int j=i+1;j<m;j++)if (a[0][k]>a[0][j])k=j;if (k!=i){t=a[0][i];a[0][i]=a[0][k];a[0][k]=t;t=a[1][i];a[1][i]=a[1][k];a[1][k]=t;}}return 0;}double Xianxing (double a[][40], int m,double w){double s;int t;for (int i=1;i<m;i++)if (w<=a[0][i])t=i;if (w>a[0][m-1])t=m-1;s=((w-a[0][t])/(a[0][t-1]-a[0][t]))*a[1][t-1]+(w-a[0][t-1])/(a[0][t]-a[0][t-1])*a[1][t];cout <<"插值结果为:( "<<w<<","<< s<<" )"<<endl;return 0;}double Paowu (double a[][40],int m,double w){double s;int t=0;if (w<a[0][1])t=1;for (int i=2;i<m-1;i++){if (w<=a[0][i]){if (fabs(w-a[0][i])>=fabs(w-a[0][i-1]))t=i-1;elset=i;break;}}if (w>a[0][m-2])t=m-2;s=((w-a[0][t])*(w-a[0][t+1])/((a[0][t-1]-a[0][t])*(a[0][t-1]-a[0][t+1])))*a[1][t-1]+ ((w-a[0][t-1])*(w-a[0][t+1])/((a[0][t]-a[0][t-1])*(a[0][t]-a[0][t+1])))*a[1][t]+ ((w-a[0][t-1])*(w-a[0][t])/((a[0][t+1]-a[0][t-1])*(a[0][t+1]-a[0][t])))*a[1][t+1];cout <<"插值结果为:( "<<w<<","<<s<<" )"<<endl;return 0;}int main( ){int t,m;double a[2][40]={0},w;cout <<"请输入已知点的个数:"<<endl;cout <<"m=";cin >>m;cout <<"请输入w的值:"<< endl;cout <<"w=";cin >>w;cout <<"请输入函数矩阵:"<< endl;for (int i=0;i<2;i++)for(int j=0;j<m;j++)cin >>a[i][j];Paixu(a,m);cout <<"经过排序后矩阵为:A="<<endl;for (int i=0;i<2;i++){for(int j=0;j<m;j++)cout <<setw(5)<<a[i][j];cout <<endl;}for (int i=0;;i++){cout <<"如果选择分段线性插值法请输入“”,如果选择抛物插值法请输入“”:";cin >>t;if (t==1){Xianxing(a,m,w);break;}if (t==2){Paowu(a,m,w);break;}if (t!=1 && t!=2)cout <<"输入有误!!!!请重新输入!"<<endl;}cin >>t;return 0;}。
用MATLAB实现拉格朗日插值和分段线性插值.

用MATLAB实现拉格朗日插值和分段线性插值1、实验内容:用MATLAB实现拉格朗日插值和分段线性插值。
2、实验目的:1学会使用MATLAB软件;2会使用MATLAB软件进行拉格朗日插值算法和分段线性差值算法;3、实验原理:利用拉格朗日插值方法进行多项式插值,并将图形显式出来。
4、实验步骤及运行结果(1实现lagrange插值1定义函数:f = 1/(x^2+1 将其保存在 f.m 文件中,具体程序如下:function y = f1(xy = 1./(x.^2+1;2定义拉格朗日插值函数:将其保存在lagrange.m 文件中,具体实现程序编程如下:function y = lagrange(x0,y0,xm = length(x; /区间长度/n = length(x0;for i = 1:nl(i = 1;endfor i = 1:mfor j = 1:nfor k = 1:nif j == kcontinue;endl(j = ( x(i -x0(k/( x0(j - x0(k *l(j;endendendy = 0;for i = 1:ny = y0(i * l(i + y;end3建立测试程序,保存在text.m文件中,实现画图:x=-5:0.001:5;y=(1+x.^2.^-1;p=polyfit(x,y,n;py=vpa(poly2sym(p,10plot_x=-5:0.001:5;f1=polyval(p,plot_x;figureplo t(x,y,‘r',plot_x,f1输入n=6,出现下面的图形:通过上图可以看到当n=6是没有很好的模拟。
于是重新运行text.M并选择n=11由此可见n=11时的图像是可以很好的实现模拟(2分段线性插值:建立div_linear.m文件。
具体编程如下/*分段线性插值函数:div_linear.m 文件*/ function y = div_linear(x0,y0,x,n%for j = 1:length(xfor i = 1:n-1if (x >= x0(i && (x <= x0(i+1y = (x - x0(i+1/(x0(i - x0(i+1*y0(i + ( x - x0(i/(x0(i+1 - x0(i*y0(i+1;elsecontinue;endend%end测试程序(text2.m:输入n =:’;n = input(‘x0 = linspace( -5,5,n;for x = -5:0.01:5y = div_linear(x0,f(x0,x,n;hold on;plot(x,y,'r';plot(x,f(x,'b';end2运行测试程序,这是会出现:输入n=:2输入n=6,并按Enter键,出现:4关掉图形界面后,重新运行程序,输入n=11,并按enter键后出现:5再次关掉图形界面,输入n=100,并按enter键,出现:此时。
Matlab实验报告六(三次样条与分段线性插值)

实验名称插值与拟合
所属课程数学软件与实验
实验类型综合型实验
专业信息与计算科学
班级
学号
姓名
指导教师
一、实验概述
【实验目的】
学会在matlab环境下使用几种不同的插值法和拟合两种方法构造函数依据已经知道的某些特殊点来推测实际问题中需要知道但又不便于测量出来的量。
【实验原理】
1.z=interp2(x0,y0,z0,x,y,’method’): 要求x0,y0单调;x, y可取为矩阵, 或x取行向量, y取为列向量, x,y的值分别不能超出x0,y0的范围。
2.分段线性插值与计算量与n无关;n越大, 误差越小.
3.三次样条插值比分段线性插值更光滑。
4.‘linear’ : 分段线性插值;‘spline’ : 三次样条
二、实验内容
问题1 对函数, x([-5,5], 分别用分段线性插值和三次样条插值作插值(其中插值节点不少于20), 并分别作出每种插值方法的误差曲线.
1180 1320 1450 1420 1400 1300 700 900];
mesh(x,y,z)
xi=0:20:2800;
yi=0:20:2400;
zi=interp2(x,y,z,xi',yi,'cubic');
mesh(xi,yi,zi)
3.结果
4.结论及分析
通过实验,结果正确,分析无误。
三、实验小结
1270 1500 1200 1100 1350 1450 1200 1150
1230 1390 1500 1500 1400 900 1100 1060
1180 1320 1450 1420 1400 1300 700 900
lagrange插值分段线性插值matlab代码

Lagrange插值:x=0:3;y=[-5,-6,-1,16];n=length(x);syms q;for k=1:nfenmu=1;p=1;for j=1:nif(j~=k)fenmu=fenmu*(x(k)-x(j)) p=conv(p,poly(x(j)))endendc(k,:)=p*y(k)/fenmuenda=zeros(1,n);for i=1:nfor j=1:na(i)=a(i)+c(j,i)endend输出结果:fenmu =-1p =1 -1fenmu =2p =1 -3 2fenmu =-6p =1 -6 11 -6c =0.8333 -5.0000 9.1667 -5.0000 fenmu =1p =1 0fenmu =-1p =1 -2 0fenmu =2p =1 -5 6 0c =0.8333 -5.0000 9.1667 -5.0000-3.0000 15.0000 -18.0000 0 fenmu =2p =1 0fenmu =2p =1 -1 0fenmu =-2p =1 -4 3 0c =0.8333 -5.0000 9.1667 -5.0000-3.0000 15.0000 -18.0000 00.5000 -2.0000 1.5000 0 fenmu =3p =1 0fenmu =6p =1 -1 0fenmu =6p =1 -32 0c =0.8333 -5.0000 9.1667 -5.0000-3.0000 15.0000 -18.0000 00.5000 -2.0000 1.5000 02.6667 -8.0000 5.3333 0a =0.8333 0 0 0a =-2.1667 0 0 0 a =-1.6667 0 0 0a =1 0 0 0a =1 -5 0 0a =1 10 0 0a =1 8 0 0a =1 0 0 0a =1.0000 0 9.1667 0a =1.0000 0 -8.8333 0a =1.0000 0 -7.3333 0a =1.0000 0 -2.0000 0a =1.0000 0 -2.0000 -5.0000a =1.0000 0 -2.0000 -5.0000a =1.0000 0 -2.0000 -5.0000a =1.0000 0 -2.0000 -5.0000 分段线性插值:先保存M文件:x=1:6;y=[7 16 8 25 12 24];u=5.3;delta=diff(y)./diff(x);n=length(x);for j=2:(n-1)if x(j)<uk=j;endend在command window中输入:s=u-x(k);v=y(k)+s.*delta(k)输出结果:v =15.6000解:第一种做法,用spline,共55个点,其中,54个有效首先保存你一个M文件:figure('position',get(0,'screensize'))axes('position',[0 0 1 1])[x,y] = ginput;然后在command window里输入以下容:n = length(x);s = (1:n)';t = (1:.05:n)';u = spline(s,x,t);v = spline(s,y,t);clf resetplot(x,y,'.',u,v,'-');对应的x、y值:0.3572917 0.25361450.3572917 0.29096390.3503472 0.34036140.3461806 0.42590360.3427083 0.52710840.3253472 0.61626510.3065972 0.68734940.290625 0.75240960.2892361 0.79337350.2954861 0.7969880.3225694 0.75481930.340625 0.6849398 0.3690972 0.6150602 0.3864583 0.6126506 0.3899306 0.7259036 0.3927083 0.8066265 0.3920139 0.8993976 0.4024306 0.9295181 0.4239583 0.8933735 0.4239583 0.8078313 0.4295139 0.7343373 0.4315972 0.6451807 0.4440972 0.6439759 0.4565972 0.7439759 0.4704861 0.8451807 0.4767361 0.9054217 0.4961806 0.9463855 0.5086806 0.876506 0.5045139 0.8186747 0.5010417 0.7524096 0.4892361 0.6403614 0.503125 0.6295181 0.5052083 0.6271084 0.5322917 0.7090361 0.5510417 0.763253 0.5739583 0.8355422 0.5961806 0.8572289 0.5947917 0.7837349 0.5753472 0.7090361 0.5579861 0.6391566 0.5357639 0.5668675 0.5322917 0.5283133 0.5350694 0.4789157 0.565625 0.536747 0.5947917 0.5933735 0.6253472 0.610241 0.6322917 0.5728916 0.615625 0.5331325 0.6003472 0.4993976 0.5788194 0.4415663 0.559375 0.3716867 0.5295139 0.2957831 0.4975694 0.2403614 0.4711806 0.2018072 0.6607639 0.3090361第二种做法,用pchip,共52个点,全部有效首先保存一个M文件:figure('position',get(0,'screensize'))axes('position',[0 0 1 1])[x,y] = ginput;然后在command window里输入以下容:n = length(x);s = (1:n)';t = (1:.05:n)';u = pchip (s,x,t);v = pchip (s,y,t);clf resetplot(x,y,'.',u,v,'-');对应的x、y值:0.5190972 0.84879520.5052083 0.75120480.4947917 0.67891570.5100694 0.66927710.5399306 0.73554220.5753472 0.81746990.596875 0.86204820.6190972 0.87771080.6149306 0.81385540.5878472 0.74277110.5878472 0.74277110.5635417 0.67168670.5350694 0.6030120.528125 0.5632530.528125 0.52590360.565625 0.58012050.6052083 0.62710840.634375 0.61867470.6190972 0.57168670.5878472 0.5234940.5364583 0.41265060.4961806 0.32108430.459375 0.2753012我更喜欢第一种,用spline的,这个能将之间画出弧度,而pchip更像是直接用线段将点依次连接得到的。
分段线性插值

1、实验目的:设在区间[a,b]上,给定n+1个插值节点a=X0<X1<X2<。
<Xn=b和相应的函数值y0,y1,y2,。
,yn ,求一个插值函数φ(x),具有下面性质:(1)φ(x)=yj(j=0,1,2,…,n)(2)φ(x)在每一个小区间[xj,yj+1]上是线性函数。
2、算法分析:●分段线性插值的算法思想:分段线性插值需要在每个插值节点上构造分段线性插值基函数,然后再作它们的线性组合。
分段线性插值基函数的特点是在对应的插值节点上函数值取1,其它节点上函数值取0。
插值基函数如下:有了基函数就可以直接写出分段线性插值函数的表达式。
●具体程序设计:for(i=0;i<=20;i++) //选取节点{ax[i]=-1+((2/20.0)*i); //选取上的21个对称的节点ay[i]=1.0/(1+25*ax[i]*ax[i]);}x1=-1;while(x1<=1){m=0;x1=x1+0.00001;for(i=0;i<=20;i++){if(i==0&&x1>=-1&&x1<=-0.9){n=ay[0]*((x1-ax[1])/(ax[0]-ax[1])); //计算并和相应的函数值组合m=n+m;}else{if(x1>=ax[i-1]&&x1<=ax[i]) //计算并和相应的函数值组合{n=ay[i]*((x1-ax[i-1])/(ax[i]-ax[i-1]));m=n+m;}else{if(x1>ax[i]&&x1<=ax[i+1]){n=ay[i]*((x1-ax[i+1])/(ax[i]-ax[i+1]));m=m+n;}else{if(i==19&&x1>ax[19]&&x1<=ax[20]){n=ay[20]*((x1-ax[19])-(ax[20]-ax[19]));//计算并和相应的函数值组合m=n+m;}}}}}}3、实验结果截图:在[-1,1]区间上选取了21个等分节点的分段线性插值函数的图像如下:4、程序代码// SHIYAN456View.cpp : implementation of the CSHIYAN456View class //#include "stdafx.h"#include "SHIYAN456.h"#include "SHIYAN456Doc.h"#include "SHIYAN456View.h"#ifdef _DEBUG#define new DEBUG_NEW#undef THIS_FILEstatic char THIS_FILE[] = __FILE__;#endif/////////////////////////////////////////////////////////////////////////////// CSHIYAN456ViewIMPLEMENT_DYNCREATE(CSHIYAN456View, CView)BEGIN_MESSAGE_MAP(CSHIYAN456View, CView)//{{AFX_MSG_MAP(CSHIYAN456View)ON_COMMAND(ID_FFunction, OnFFunction)ON_COMMAND(ID_Lagrange, OnLagrange)ON_COMMAND(ID_Subsection, OnSubsection)ON_COMMAND(ID_Hermite, OnHermite)//}}AFX_MSG_MAP// Standard printing commandsON_COMMAND(ID_FILE_PRINT, CView::OnFilePrint)ON_COMMAND(ID_FILE_PRINT_DIRECT, CView::OnFilePrint)ON_COMMAND(ID_FILE_PRINT_PREVIEW, CView::OnFilePrintPreview) END_MESSAGE_MAP()/////////////////////////////////////////////////////////////////////////////// CSHIYAN456View construction/destructionCSHIYAN456View::CSHIYAN456View(){// TODO: add construction code here}CSHIYAN456View::~CSHIYAN456View(){}BOOL CSHIYAN456View::PreCreateWindow(CREATESTRUCT& cs)// TODO: Modify the Window class or styles here by modifying // the CREATESTRUCT csreturn CView::PreCreateWindow(cs);}/////////////////////////////////////////////////////////////////////////////// CSHIYAN456View drawingvoid CSHIYAN456View::OnDraw(CDC* pDC){CSHIYAN456Doc* pDoc = GetDocument();ASSERT_VALID(pDoc);// TODO: add draw code for native data hereAfxGetMainWnd()->SetWindowText("实验四五六函数图像");for(int k=650;k<=690;k++){pDC->SetPixel(k,55,RGB(255,0,0));pDC->SetPixel(k,85,RGB(0,255,0));pDC->SetPixel(k,115,RGB(0,0,255));pDC->SetPixel(k,145,RGB(0,255,255));}pDC->TextOut(700,50,"原函数图像");pDC->TextOut(700,80,"Lagrange插值函数图像");pDC->TextOut(700,110,"分段线性插值函数图像");pDC->TextOut(700,140,"Hermite插值函数图像");for(int i=6;i<=600;i++)pDC->SetPixel(400,i,RGB(0,0,0));pDC->TextOut(395,4,"y");for(int j=100;j<=700;j++)pDC->SetPixel(j,500,RGB(0,0,0));pDC->TextOut(700,500,"x");//pDC->MoveTo(400,500);}/////////////////////////////////////////////////////////////////////////////// CSHIYAN456View printingBOOL CSHIYAN456View::OnPreparePrinting(CPrintInfo* pInfo){// default preparationreturn DoPreparePrinting(pInfo);}void CSHIYAN456View::OnBeginPrinting(CDC* /*pDC*/, CPrintInfo* /*pInfo*/) {// TODO: add extra initialization before printing}void CSHIYAN456View::OnEndPrinting(CDC* /*pDC*/, CPrintInfo* /*pInfo*/) {// TODO: add cleanup after printing}/////////////////////////////////////////////////////////////////////////////// CSHIYAN456View diagnostics#ifdef _DEBUGvoid CSHIYAN456View::AssertValid() const{CView::AssertValid();}void CSHIYAN456View::Dump(CDumpContext& dc) const{CView::Dump(dc);}CSHIYAN456Doc* CSHIYAN456View::GetDocument() // non-debug version is inline {ASSERT(m_pDocument->IsKindOf(RUNTIME_CLASS(CSHIYAN456Doc)));return (CSHIYAN456Doc*)m_pDocument;}#endif //_DEBUG/////////////////////////////////////////////////////////////////////////////// CSHIYAN456View message handlersvoid CSHIYAN456View::OnFFunction(){// TODO: Add your command handler code hereCClientDC dr(this);COLORREF rgb=RGB(255,0,0);double x1,y1,x,y;x1=-1.0;y1=1/(1+25*x1*x1);x=x1*200+400;y=-y1*200+500;while(x1<=1){dr.MoveTo(int(x),int(y));x1=x1+0.00001;y1=1/(1+25*x1*x1);x=x1*200+400;y=-y1*200+500;dr.SetPixel(int(x),int(y),rgb);}}void CSHIYAN456View::OnLagrange(){// TODO: Add your command handler code hereCClientDC dr(this);COLORREF rgb=RGB(0,255,0);int i,j;double x1=0,y1=0,x=0,y=0,m=0,n=0,ax[100],ay[100];for(i=0;i<=10;i++){ax[i]=-1+((2/10.0)*i);ay[i]=1.0/(1+25*ax[i]*ax[i]);}x1=-1;y1=1/(1+25*x1*x1);x=x1*200+400;y=-y1*200+500;dr.MoveTo(int(x),int(y));while(x1<=1){m=0;for(i=0;i<=10;i++){n=1;for(j=0;j<=10;j++){if(i!=j)n=((x1-ax[j])/(ax[i]-ax[j]))*n;}m=ay[i]*n+m;}x=x1*200+400;y=-m*200+500;dr.SetPixel(int(x),int(y),rgb);x1=x1+0.00001;}}void CSHIYAN456View::OnSubsection(){// TODO: Add your command handler code hereCClientDC dr(this);COLORREF rgb=RGB(0,0,255);int i;double x1=0,y1=0,x=0,y=0,m=0,n=0,ax[100],ay[100];for(i=0;i<=20;i++){ax[i]=-1+((2/20.0)*i);ay[i]=1.0/(1+25*ax[i]*ax[i]);}x1=-1;y1=1/(1+25*x1*x1);x=x1*200+400;y=-y1*200+500;dr.MoveTo(int(x),int(y));while(x1<=1){m=0;x1=x1+0.00001;for(i=0;i<=20;i++){if(i==0&&x1>=-1&&x1<=-0.9){n=ay[0]*((x1-ax[1])/(ax[0]-ax[1]));m=n+m;}else{if(x1>=ax[i-1]&&x1<=ax[i]){n=ay[i]*((x1-ax[i-1])/(ax[i]-ax[i-1]));m=n+m;}else{if(x1>ax[i]&&x1<=ax[i+1]){n=ay[i]*((x1-ax[i+1])/(ax[i]-ax[i+1]));m=m+n;}else{if(i==19&&x1>ax[19]&&x1<=ax[20]){n=ay[20]*((x1-ax[19])-(ax[20]-ax[19]));m=n+m;}}}}}x=x1*200+400;y=-m*200+500;dr.SetPixel(int(x),int(y),rgb);}}void CSHIYAN456View::OnHermite(){// TODO: Add your command handler code hereCClientDC dr(this);COLORREF rgb=RGB(0,255,255);int i;double x1=0,y1=0,x=0,y=0,m=0,n=0,tt=0,tt1=0,h=0,ax[100],ay[100],a[100];for(i=0;i<=10;i++){ax[i]=-1+((2/10.0)*i);ay[i]=1.0/(1+25*ax[i]*ax[i]);a[i]=(-50*ax[i])/((1+25*ax[i]*ax[i])*(1+25*ax[i]*ax[i]));}x1=-1;y1=1/(1+25*x1*x1);x=x1*200+400;y=-y1*200+500;dr.MoveTo(int(x),int(y));while(x1<=1){m=0;x1=x1+0.00001;for(i=0;i<=10;i++){if(i==0&&x1>=ax[0]&&x1<=ax[1]){n=(1+2*((x1-ax[0])/(ax[1]-ax[0])))*((x1-ax[1])/(ax[0]-ax[1]))*((x1-ax[1])/(ax[0]-ax[1]));tt=ay[0]*n;h=(x1-ax[0])*((x1-ax[1])/(ax[0]-ax[1]))*((x1-ax[1])/(ax[0]-ax[1]));tt1=a[0]*h;m=tt+tt1+m;}else{if(x1>=ax[i-1]&&x1<=ax[i]){n=(1+2*((x1-ax[i])/(ax[i-1]-ax[i])))*((x1-ax[i-1])/(ax[i]-ax[i-1]))*((x1-ax[i-1])/(ax[i]-ax[i-1]));tt=ay[i]*n;h=(x1-ax[i])*((x1-ax[i-1])/(ax[i]-ax[i-1]))*((x1-ax[i-1])/(ax[i]-ax[i-1]));tt1=a[i]*h;m=tt+tt1+m;}else{if(x1>ax[i]&&x1<=ax[i+1]){n=(1+2*((x1-ax[i])/(ax[i+1]-ax[i])))*((x1-ax[i+1])/(ax[i]-ax[i+1]))*((x1-ax[i+1])/(ax[i]-ax[i+1]));tt=ay[i]*n;h=(x1-ax[i])*((x1-ax[i+1])/(ax[i]-ax[i+1]))*((x1-ax[i+1])/(ax[i]-ax[i+1]));tt1=a[i]*h;m=tt+tt1+m;}else{if(i==9&&x1>ax[9]&&x1<=ax[10]){n=(1+2*((x1-ax[10])/(ax[9]-ax[10])))*((x1-ax[9])/(ax[10]-ax[9]))*((x1-ax[9])/(ax[10]-ax[9]));tt=ay[i]*n;h=(x1-ax[10])*((x1-ax[9])/(ax[10]-ax[9]))*((x1-ax[9])/(ax[10]-ax[9]));tt1=a[i]*h;m=tt+tt1+m;}}}}}x=x1*200+400;y=-m*200+500;dr.SetPixel(int(x),int(y),rgb);}}5、总结体会分段线性插值的方法克服了Lagrange插值法当节点不断加密时,构造的插值多项式的次数不断升高,高次多项式虽然是连续的,但是不一定都收敛到相应的被插函数而产生Runge现象。
(matlab)分段线性插值代码

disp('给出插值区间的等分数n')
n=input('n = ');
%作被插函数图象
u=-5:(10/200):5;
v=1./(1+u.^2);
plot(u,v,'r');
hold on%固化图形屏幕
%给出插值条件
x=-5:(10/n):5;
y=1./(1+x.^2);
%取需用分段线性插值函数计算其函数值的点z
break
else
i=i+1;
end
end
end
plot(z,Ih,'b')
hold off%释放固化的图形屏幕
z=-5:10/(2*n):5;
m=max(size(z));
%计算分段线性插值函数在z处的值,并作图
for k=1:m
i=1;
while i<=n
if z(k)>=x(i) & z(k)<=x(i+1)
Ih(k)=y(i)*(z(k)-x(i+1))/(x(i)-x(i+1))(i+1)-x(i));
文档之家的所有文档均为用户上传分享文档之家仅负责分类整理如有任何问题可通过上方投诉通道反馈
(matlab)分段线性插值代码
% y=1/(1+x^2) 分段线性xi=-5+10*i/n(i=1,1,...,n)。比较发现,随着n的增大,两者吻合得越来越好,龙格现象并未发生
clear all %清除命令空间中所有变量
分段线性插值函数的编程实现

1 问题的提出对21()1f x x=+在(-5,5)上进行分段线性插值,取不同节点个数n ,得到不同分段线性插值函数.虽然MATLAB 里有直接分段线性插值的函数,但为了对分段插值算法有更明确的理解,编写该程序是有必要的.需要解决的问题:1、 由已知数据节点编写程序,实现分段线性插值函数,从而能由所编函数得到非节点的函数值.2、 比较用不同节点数所得插值函数与真实函数的误差,从而得出节点数与插值效果的关系.2 理论基础分段线性插值适用于计算简单、光滑性要求不高的插值问题,且其整体逼近)(x f 的效果较好.从几何意义上看,分段线性插值就是用折线近似代替曲线错误!未找到引用源。
.设在区间[a,b]上取n+1个点.a 110b x x x x n n =<<<<=-函数)(x f 在上述节点处的函数值为)(y i i x f = ,,2,1,0)(n i =于是得到n+1个点).,(,),,(,,1100n n y x y x y x )(连接相邻两点错误!未找到引用源。
和),(11++i i y x ,,2,1,0)(n i =,得一折线函数)(x ϕ,若满足:(1))(x ϕ在[a,b]上连续;(2)错误!未找到引用源。
),,2,1,0(n i =;(3))(x ϕ在每个小区间错误!未找到引用源。
上是线性函数, 则称折线函数)(x ϕ为分段线性插值函数. 模型一:由分段线性插值函数的定义可知,)(x ϕ在每个小区间错误!未找到引用源。
上可表为,)(1111++++--+--=i ii ii i i i y x x x x y x x x x x ϕ错误!未找到引用源。
)1-,,2,1,0(n i =.)(x ϕ是一分段函数,若用基函数表示,只需对1,,2,1,0-=n i 令⎪⎪⎪⎩⎪⎪⎪⎨⎧=≤≤--=≤≤--=+++---,其他略去略去0)(,)0(,)(111111n i x x x x x x x i x x x x x x x x l i i i i i i i i i i i 显然,()i l x 是分段的线性连续函数,且满足⎩⎨⎧≠==k i ki x l k i ,0,1)( 于是∑==ni i i x l y x 0),()(ϕ b x a ≤≤模型二:首先确定间隔序列i ,使得:1+≤≤i i x x x第二个量是局部变量s ,其定义为 :i x x s -=最后一个量是一阶均差ii ii i x x y y --=++11δ则插值基函数可表示为i i ii ii k i s y x x y y x x y x L δ+=---+=++11)()(截断误差为:))((21)(1''1i i x x x x f x R --=-)(ξ []i i x x x ,1-∈ 若,)(2''M x f ≤错误!未找到引用源。
matlab实现lagrange插值和分段线性插值

数值分析作业姓名:虞驰程题目函数:在[-5,5]上,取n=10,对其进行分段线性插值和拉格朗日插值,在Matlab中实现且绘图。
Matlab 实现:首先定义函数f,在Matlab中用fun ctio n.m 文件编写,具体代码如图1所示:图1 f(x) 函数定义分段线性插值的基本函数,用function.m 文件编写,具体代码如图2所示:图2分段线性插值基本函数定义拉格朗日插值的基本函数,用function.m 文件编写,具体代码如图3所示:图3拉格朗日插值的基本函数进行分段线性插值并绘图和原函数进行对比的Matlab实现代码如图4所示:图4分段线性插值函数绘制其结果如图5所示,其中红色代表分段线性插值结果,蓝色代表原函数:图5分段线性插值和原函数对比同理可以进行拉格朗日插值并绘图,其Matlab实现代码如图6所示,其结果如图Ffe Edit Tesrt Go C«ll TooU D«bug Oesktop Wind Help *1 D 已■ | $ ■•勺◊Qi 勺-it • * F;iL ・昌电T・• ■ ■ ■ | MadgJBmr ■血*S肓匸1必」•* _U■「鬲迢[0 _L 1- IWEU I1『inpirt ftji')2 - ■Elinus:""卜氐£11*1〕3 - \-\t^£ 15* —yslftErflnfeiwQpiftrtl,!!) ■5 —told An:ft - plfft <45- r* * E X).7 - jlflrt(i p fCx)/lb J kE — d.isp(-x):9 - -end2 u"含m al=n' taufwi图6拉格朗日插值函数绘制图7拉格朗日插值和原函数对比7所示:最后我们可以将分段线性插值、拉格朗日插值和原函数进行对比,其实现代码如图8所示,最终结果如图9所示(黑色代表原函数,蓝色是分段线性插值,红色是拉格朗日插值)图8两种插值方法和原函数对比实现图9两种插值方法和原函数对比。
第二章分段低次插值

1. 分段线性Lagrange插值的算法设计
程序:lagrange1.m
2. 分段二次Lagrange插值的算法设计
程序:lagrange2.m
16
-4 -3 -2 -1 0 1 2 3 4 5
2. 分段线性插值的误差估计
由第二节定理1可知,n次Lagrange插值多项式的余项为
Rn (x) f (x) Pn (x)
f (n1) ( )
(n 1)!
n 1
(
x)
那么分段线性插值 L1(x)的余项为 R1(x) f (x) L1(x) f ( x) L(1k )( x)
yk
1
(
(x xk 1
xk xk
)( )(
x xk 1 ) xk 1 xk 1
)
yk
( x xk 1 )( x xk 1 ) ( xk xk 1 )( xk xk 1 )
yk 1
(x ( xk 1
xk 1 )( x xk ) xk 1 )( xk 1 xk )
k 1,2,,n 1
0.81343
f (0.98) L(24)(0.98) 1.09784
f (1.1) L(24)(1.1) 1.25513
15
分段低次Lagrange插值的特点 计算较容易 可以解决Runge现象 但插值多项式分段 插值曲线在节点处会出现尖点 插值多项式在节点处不可导
三、分段低次插值的算法设计
yk 1
x * xk xk 1 xk
内插
若x* x0
外插
取
y* L1(x*)
L(10)( x*)
y0
x * x1 x0 x1
y1
x * x0 x1 x0
数值分析分段线性插值样条插值Runge函数Newton-Lagrange-Chebyshev插值多项式Runge现象matlab源程序代码

题目1:对Runge 函数R(x ) =1在区间[-1,1]作下列插值逼近,并和1 + 25x 2R(x)的图像进行比较,并对结果进行分析。
= -1 + ih,h= 0.1,0 ≤ i≤ 20 绘出它的20 次Newton 插值(1)用等距节点xi多项式的图像。
分别画出在[-1,1]区间,[-0.7,0.7]区间和[-0.5,0.5]区间上的 Newton 插值多项式和Runge 函数的图像从图中可以看出: 1)在[-0.5,0.5]区间 Newton 插值多项式和 Runge 函数的图像偏差较小,这说 明 Newton 插值多项式在此区间上可以较好的逼近 Runge 函数; 2)在边界(自变量 x=-1 和 x=1)附近,Newton 插值多项式和 Runge 函数的图像 偏差很大,Newton 插值多项式出现了剧烈的震荡。
(Runge 现象) (2)用节点 x = cos(2i + 1π)(, i = 0,1,2,⋅ ⋅ ⋅ ,20),绘出它的 20 次 Lagrangei 42 插值多项式的图像。
画出在[-1,1]区间上的 Lagrange 插值多项式和 Runge 函数的图像从图中可以看出:使用 Chebyshev 多项式零点构造的 Lagrange 插值多项式和 Runge 函数的图 像偏差较小,没有出现 Runge 现象。
(3)用等距节点 x i 的图像。
= -1 + ih ,h = 0.1,0 ≤ i ≤ 20 绘出它的分段线性插值函数画出在[-1,1]区间上分段线性插值函数和 Runge 函数的图像从图中可以看出:使用分段线性插值函数和 Runge 函数的图像偏差较小,也没有出现 Runge 现象,只在自变量 x=0 处有稍许偏差。
(4)用等距节点 x i 函数的图像。
= -1 + ih ,h = 0.1,0 ≤ i ≤ 20 绘出它的三次自然样条插值画出在[-1,1]区间上三次自然样条插值函数和 Runge 函数的图像从图中可以看出:使用三次自然样条插值函数和 Runge 函数的图像偏差较小,也没有出现 Runge 现象。
分段插值法

x2
x2 )
)
y2
(0.42 ( x2
x0 )(0.42 x1 ) x0 )( x2 x1 )
0.43281
f (0.75)
L(24 ) (0.75)
y3
(0.75 ( x3
x4 )(0.75 x5 ) x4 )( x3 x5 )
y4
(0.75 ( x4
x3 )(0.75 x5 ) x5 )( x4 x5 )
并不一定会有更好的插值结果这是因为高次多项式的振荡是很厉害的
§
从上节可知,如果插值多项式的次数(cìshù)过高,可能产生 Runge现象,因此,在构造插值多项式时常采用分段
插值的方法。
设插值节点为 xi , 函数值为 yi , i 0,1,, n
hi xi1 xi , i 0,1,2,, n 1
h
max i
hi
任取两个相邻的节点 xk , xk 1 ,形成一个插值区间 [xk , xk 1 ]
构造(gòuzào)Lagrange线性插值
1
共十六页
L(1k )( x) yklk ( x) yk 1lk 1( x)
yk
x xk 1 xk xk 1
yk 1
x xk xk 1 xk
k 0,1,,n 1
0.43307
同理
f (0.75) L(13)(0.75) 0.81448
f (0.98) L(14)(0.98) 1.10051
f (1.1)
L(14 ) ( 1.1)
0.87335 1.1 1.05 0.8 1.05
1.18885 1.1 0.8 1.05 0.8
1.25195
f
(
分段线性插值函数的编程实现摘要

分段线性插值函数的编程实现摘 要在代数插值过程中,人们为了获得较好的近似效果,通常情况下是增加插值节点数.由于二次插值比线性插值近似效果好,因此容易错误地认为插值多项式次数越高越好.事实上,随着插值节点的增多,插值多项式不一定收敛到被插值函数.通过分段低次插值或样条插值可以得到较好的近似逼近函数,分段低次插值具有公式简单、运算量小、稳定性好、收敛性有保证等优点.随着子区间长度h 取得足够小,分段低次插值总能满足所要求的精度.因此分段低次插值应用十分广泛. 分段线性插值是分段低次插值中常见的方法之一,在本文中对21()1f x x =+在(-5,5)上进行分段线性插值,取不同节点个数n ,得到不同分段线性插值函数.并用MATLAB 编写分段线性插值函数,最后比较用不同节点数所得插值函数与真实函数的误差,从而得出节点数与插值效果的关系.关键字:线性插值;分段低次插值;分段线性插值函数;MATLAB 软件USING PROGRAMMING TO ACHIEVE THE PIECEWISE LINEAR INTERPORATION FUNCTIONABSTRACTIn the process of algebra interpolation, people typically increase the interpolation node number in order to get a better approximate result. As the result calculated by quadratic interpolation is better than the same aspect by linear interpolation, it is wrongly thought that the higher power of interpolation the better the result is. In fact, with the increase of the number of node interpolation method, interpolation polynomial is not convergent as interpolation function.It can get much better approximation function through the piecewise low powerof interpolation or spline interpolation. Piecewise low power of interpolation is comparably simple and the stability and the convergence can be guaranteed. As the length h of the subinterval approaches to zero, piecewise low power interpolation can meet the required accuracy. This method is used widely.Piecewise linear interpolation is one of the typical methods of piecewise low power of interpolation, in this essay, we will conduct linear interpolation from (-5,5) by different node number n and relative piecewise linear interpolation functions. Then, we will programme linear interpolation by matlab and compare the error between the result derived from linear interpolation and real function. The method is used to test the relationship between different lengths and relative results of interpolation.Key words: linear interpolation; piecewise low power of interpolation; piecewise linear interpolation function; MATLAB software目录1 问题的提出 (1)2 理论基础 (1)3 编程过程 (3)4 结果分析 (7)5 结论 (7)参考文献 (8)附录 (9)。
C#分段线性插值函数

C#分段线性插值函数由于项目需要,需要将数据采集得到的点数转化为固定点数,使用分段线性插值其实现代码如下:/// <summary>/// 分段线性插值,将一组数插值为所需点数/// </summary>/// <param name="dataIn">待插值的数据数组</param>/// <param name="n">插值点数</param>/// <returns>插值后的数据数组</returns>public static double[] Interpolation(double[] dataIn,int n){double[] dataOut = new double[n];int lenIn = dataIn.Length;double[] a = new double[lenIn];double[] divIn = new double[lenIn];double[] divOut = new double[n];divIn[0] = 0;for (int i = 1; i < lenIn; i++){divIn[i] = divIn[i - 1] + 1;}divOut[0] = 0;for (int i = 1; i < n; i++){divOut[i] = divOut[i - 1] + lenIn / Convert.ToDouble(n);}int k = 0;for (int i = k; i < n; i++){for (int j = 0; j < lenIn - 1; j++){if (divOut[i] >= divIn[j] && divOut[i] < divIn[j + 1]){dataOut[i] = (dataIn[j + 1] - dataIn[j]) * (divOut[i] - divIn[j]) / (divIn[j + 1] - divIn[j]) + dataIn[j];k = i;}}}return dataOut;}。
分段低次插值

七、分段低次插值7.1 分段线性插值7.1.1 代码PROGRAM ILIMPLICIT NONEINTEGER :: N,K,I,N2REAL*8 :: X(60),Y(60),L(60),X0(60),Y0(60)OPEN(unit=11,file='INPUT.txt')OPEN(unit=22,file='OUTPUT.txt')READ(11,*)NREAD(11,*)X(1:N)READ(11,*)Y(1:N)READ(11,*)N2READ(11,*)X0(1:N2)DO I=1,N2,1DO K=1,N-1,1IF (X0(I)>X(K)) THEN !!找寻插值点所在的区间IF (X0(I)<X(K+1)) THENY0(I)=((X0(I)-X(K+1))/(X(K)-X(K+1))*Y(K))+((X0(I)-X(K))/(X(K+1)-X(K))*Y(K+1))ENDIFENDIFENDDOENDDODO I=1,N2,1WRITE(22,*) Y0(I)ENDDOENDPROGRAM7.1.2 验证使用函数(7.1)构建了12个插值点来进行验证,()23=(7.1)f x x使用1.3, 1.45和1.55作为插值点,程序得到的插值结果为5.10, 6.33和7.23,根据公式(7.1)得到的结果分别为5.07, 6.31和7.20,结果较为接近。
输入文件:输出文件:7.2 分段三次艾尔米特插值7.2.1 代码PROGRAM IHIMPLICIT NONEINTEGER :: N,K,I,N2REAL*8 :: X(60),Y(60),DY(60),L(60),X0(60),Y0(60)OPEN(unit=11,file='INPUT.txt')OPEN(unit=22,file='OUTPUT.txt')READ(11,*)NREAD(11,*)X(1:N)READ(11,*)Y(1:N)READ(11,*)DY(1:N)READ(11,*)N2READ(11,*)X0(1:N2)DO I=1,N2,1DO K=1,N-1,1IF (X0(I)>X(K)) THENIF (X0(I)<X(K+1)) THENY0(I)=(1+2*(X0(I)-X(K))/(X(K+1)-X(K)))*(((X0(I)-X(K+1))/(X(K)-X(K+1)))**2)*Y(K)+(1+2*(X0(I)-X(K+1))/(X(K)-X(K+1)))*(((X0(I)-X(K))/(X(K+1)-X(K)))**2)*Y(K+1)+(X0(I)-X(K))*(((X0(I)-X(K+1))/(X(K)-X(K+1)))**2)*DY(K)+(X0(I)-X(K+1))*(((X0(I)-X(K))/(X(K+1)-X(K)))**2)*DY(K+1)ENDIFENDIFENDDOENDDODO I=1,N2,1WRITE(22,*) Y0(I)ENDDOENDPROGRAM7.2.2 验证使用与7.1.2节中完全相同的验证公式和数据,同样使用1.3, 1.45和1.55作为插值点,此时程序得到的插值结果为5.07, 6.31和7.21,与公式(7.1)得到的结果5.07, 6.31和7.20更为接近。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
1 问题的提出对21()1f x x=+在(-5,5)上进行分段线性插值,取不同节点个数n ,得到不同分段线性插值函数.虽然MATLAB 里有直接分段线性插值的函数,但为了对分段插值算法有更明确的理解,编写该程序是有必要的.需要解决的问题:1、 由已知数据节点编写程序,实现分段线性插值函数,从而能由所编函数得到非节点的函数值.2、 比较用不同节点数所得插值函数与真实函数的误差,从而得出节点数与插值效果的关系.2 理论基础分段线性插值适用于计算简单、光滑性要求不高的插值问题,且其整体逼近)(x f 的效果较好.从几何意义上看,分段线性插值就是用折线近似代替曲线错误!未找到引用源。
.设在区间[a,b]上取n+1个点.a 110b x x x x n n =<<<<=-函数)(x f 在上述节点处的函数值为)(y i i x f = ,,2,1,0)(n i =于是得到n+1个点).,(,),,(,,1100n n y x y x y x )(连接相邻两点错误!未找到引用源。
和),(11++i i y x ,,2,1,0)(n i =,得一折线函数)(x ϕ,若满足:(1))(x ϕ在[a,b]上连续;(2)错误!未找到引用源。
),,2,1,0(n i =;(3))(x ϕ在每个小区间错误!未找到引用源。
上是线性函数, 则称折线函数)(x ϕ为分段线性插值函数. 模型一:由分段线性插值函数的定义可知,)(x ϕ在每个小区间错误!未找到引用源。
上可表为,)(1111++++--+--=i ii ii i i i y x x x x y x x x x x ϕ错误!未找到引用源。
)1-,,2,1,0(n i =.)(x ϕ是一分段函数,若用基函数表示,只需对1,,2,1,0-=n i 令⎪⎪⎪⎩⎪⎪⎪⎨⎧=≤≤--=≤≤--=+++---,其他略去略去0)(,)0(,)(111111n i x x x x x x x i x x x x x x x x l i i i i i i i i i i i 显然,()i l x 是分段的线性连续函数,且满足⎩⎨⎧≠==k i ki x l k i ,0,1)( 于是∑==ni i i x l y x 0),()(ϕ b x a ≤≤模型二:首先确定间隔序列i ,使得:1+≤≤i i x x x第二个量是局部变量s ,其定义为 :i x x s -=最后一个量是一阶均差ii ii i x x y y --=++11δ则插值基函数可表示为i i ii ii k i s y x x y y x x y x L δ+=---+=++11)()(截断误差为:))((21)(1''1i i x x x x f x R --=-)(ξ []i i x x x ,1-∈ 若,)(2''M x f ≤错误!未找到引用源。
,[]b a x ,∈令1--=i i i x x h 当[]i i x x x ,1-∈时221141)(41))((i i i i i h x x x x x x =-≤----从而错误!未找到引用源。
2221218181)()(h M h M x P x f ≤≤- i ni h h ≤≤=1max 当0→h 时有)()(1x f x P →即)(1x P 收敛于()f x 。
3 编程过程1、模型一: 用MATLAB 分别建立m 文件:(1)原函数fd1.m(2)分段线性插值函数fd2.m(3)比较不同节点数所得分段线性插值函数的插值效果fd3.m2、选取插值节点数为偶数在MATLAB 窗口中执行:fd3 n=2的数据见附录,图像如下: n=8的图如下:-5-4-3-2-101234500.51xy原函数(实线)-插值函数(虚线)-5-4-3-2-1012345-1-0.5xR (x )误差分析-5-4-3-2-101234500.51xy原函数(实线)-插值函数(虚线)-5-4-3-2-1012345-0.4-0.200.20.4xR (x )误差分析n=20的图-5-4-3-2-101234500.51xy原函数(实线)-插值函数(虚线)-5-4-3-2-1012345-0.1-0.0500.050.1xR (x )误差分析3、模型二:用MATLAB 分别建立m 文件:(1)分段插值函数fd22(2)插值效果比较函数fd32(选取插值节点数为奇数) 源程序(参见附录).在MATLAB 窗口中执行:fd32得下图:上图为不同节点数插值函数图像与原函数图像,下图为误差图像:n=3n=5n=7n=9n=11n=3n=5n=7n=9n=113、 由上所有的图可看出,由于原函数是偶函数,等距节点所得插值函数有很强对称性,下任取节点, 编写程序fd33.m ,得图.上图为不同节点数插值函数图像与原函数图像,下图为误差图像:00.20.40.60.81-1-0.50.5n=3n=5n=7n=9n=11n=3n=5n=7n=9n=114、 比较不同节点所得插值函数与被插函数误差的平方和,程序模板为d1.m 得下图:红星由fd32得奇数节点误差平方和,绿星加圈由fd3得偶数节点误差平方和,圈由f33得随机节点误差平方和,数据见附录.051015202520406080100120140n 节点数e 误差平方和4 结果分析1、不同插值节点数所得的分段线性插值函数,在节点处与原函数的函数值一定相同.2、所得的分段线性插值函数在原函数斜率绝对值变化大的地方,与原函数的误差比较大.3、由误差平方和e ,插值节点个数越多,e 有减小的趋势,最后趋于0.只考虑奇数或偶数个节点,则随节点数增加e 严格减小.4、随机生成的节点不如等距节点使插值效果好.5 结论插值节点个数越多,分段线性插值函数与原函数误差平方和有减小趋势,插值效果越好.参考文献[1] 李庆杨,王能超.数值分析.(第4版)[M].武汉:华中科技大学出版社,2006:13-24.[2] 肖筱南,赵来军.现代数值计算方法.[M].北京:北京大学出版社,2003:146-149.[3] 吴勃英,王德明.数值分析原理.[M].北京:科学出版社,2003:132-134.[4] 刘卫国编.MATLAB程序设计教程.[M].北京:中国水利水电出版社,2001:1-180.[5] 蔺小林,蒋耀林.现代数值分析.[M].北京:国防工业出版社,2004:184-186.附录程序如下:% fd1.m线性插值原函数function y=fd1(x)y=1./(1+x.^2);% fd2.m 分段线性插值函数function yi=fd2(x,y,xi)n=length(x);m=length(y);if n~=merror('X和Y向量的长度必须相同');return;endfor k=1:n-1if abs(x(k)-x(k+1))<eps % x(k)-x(k+1) 的绝对值必须大于e error('数据有误');return;endif x(k)<=xi&xi<=x(k+1) % 保证x(k) < xi < x(k+1)temp=x(k)-x(k+1);yi=(xi-x(k+1))/temp*y(k)+(xi-x(k))/(-temp)*y(k+1)return;endend% fd3.m 比较插值效果a=-5;b=5;n=input('请输入分端节点数:');if n<=0error('你输入的数据有误!');break;endh=(b-a)/(n-1); % 求节点x=a:h:b;y=fd1(x);xx=a:0.1:b; % 用分段线性插值函数求非节点函数值yyi=fd1(xx);m1=length(xx);z=zeros(1,m1);for k1=1:m1z(k1)=fd2(x,y,xx(k1));endw=z-yyi; % 计算误差subplot(2,1,1);plot(x,y,'o',xx,yyi,'-',x,y,'k:');%插值图像xlabel('x');ylabel('y');title('原函数(实线)-插值函数(虚线)');hold onsubplot(2,1,2);plot(xx,w,'k:'); % 误差的图像xlabel('x');ylabel('R(x)');title('误差分析');hold onxx=xx';yyi=yyi';z=z';w=w';% fd22.m 分段线性插值函数function v=fd22(x,y,u)delta=diff(y)./diff(x);n=length(x);k=ones(size(u));for j=2:n-1k(x(j)<=u)=j;ends=u-x(k);v=y(k)+s.*delta(k);% fd32.m同时画不同节点的插值函数图像和误差图像clearcloset=[-5:0.01:5];a=['k' 'g' 'r' 'c' 'm'];for i=1:5n=2*i+1;x=linspace(-5,5,n); %把区间[-5 5]分为(n-1)份,算插值节点y=fd1(x);p=fd22(x,y,t);p=p'; %计算以(x,y)为插值点的插值函数在t处的各个值y1=fd1(t);y1=y1';e=p-y1; %计算误差subplot(2,1,1);plot(x,y,a(i));hold on; %画出插值函数图像及误差图像subplot(2,1,2);plot(t,e,a(i));hold on;endsubplot(2,1,1);legend('n=3','n=5','n=7','n=9','n=11')subplot(2,1,2);legend('n=3','n=5','n=7','n=9','n=11')subplot(2,1,1);fplot(@fd1,[-5 5],'k'); %画出原函数图像hold off%fd33.m 插值节点非等分区间获得closet=[-5:0.01:5];a=['k' 'g' 'r' 'c' 'm'];for i=1:5n=2*i+1;x=[-5 rand(1,n-2)*10-5 5]; %得(-5,5)上的n维随机向量x=sort(x);y=fd1(x);p=fd22(x,y,t);p=p';y1=fd1(t);y1=y1';e=p-y1;subplot(2,1,1);plot(x,y,a(i));hold on;subplot(2,1,2);plot(t,e,a(i));hold on;endsubplot(2,1,1);legend('n=3','n=5','n=7','n=9','n=11')subplot(2,1,2);legend('n=3','n=5','n=7','n=9','n=11')subplot(2,1,1);fplot(@fd1,[-5 5],'k');hold off%fd1.m 比较不同节点数误差平方和cleart=[-5:0.01:5];a=[];b=[];for i=1:10n=2*i; %n=2*i+1则是奇数节点x=linspace(-5,5,n)y=fd1(x);p=fd22(x,y,t);y1=fd1(t);e=p-y1;e=e*e';a=[a e];b=[b n];endplot(b,a,'go')xlabel('n节点数')ylabel('e误差平方和')hold onn=2X Y YI(原函数)W -5.0000 0.0385 0.0385 0-4.9000 0.0400 0.0577 -0.0177-4.8000 0.0416 0.0769 -0.0353-4.7000 0.0433 0.0962 -0.0528-4.6000 0.0451 0.1154 -0.0703-4.5000 0.0471 0.1346 -0.0876-4.4000 0.0491 0.1538 -0.1047-4.3000 0.0513 0.1731 -0.1218-4.2000 0.0536 0.1923 -0.1387-4.1000 0.0561 0.2115 -0.1554-4.0000 0.0588 0.2308 -0.1719-3.9000 0.0617 0.2500 -0.1883-3.8000 0.0648 0.2692 -0.2045-3.7000 0.0681 0.2885 -0.2204-3.6000 0.0716 0.3077 -0.2361-3.5000 0.0755 0.3269 -0.2515-3.4000 0.0796 0.3462 -0.2665-3.3000 0.0841 0.3654 -0.2813-3.2000 0.0890 0.3846 -0.2956-3.1000 0.0943 0.4038 -0.3096-3.0000 0.1000 0.4231 -0.3231-2.9000 0.1063 0.4423 -0.336-2.8000 0.1131 0.4615 -0.3484 -2.7000 0.1206 0.4808 -0.3601 -2.6000 0.1289 0.5000 -0.3711 -2.5000 0.1379 0.5192 -0.3813 -2.4000 0.1479 0.5385 -0.3905 -2.3000 0.1590 0.5577 -0.3987 -2.2000 0.1712 0.5769 -0.4057 -2.1000 0.1848 0.5962 -0.4113 -2.0000 0.2000 0.6154 -0.4154 -1.9000 0.2169 0.6346 -0.4177 -1.8000 0.2358 0.6538 -0.418 -1.7000 0.2571 0.6731 -0.416 -1.6000 0.2809 0.6923 -0.4114 -1.5000 0.3077 0.7115 -0.4038 -1.4000 0.3378 0.7308 -0.3929 -1.3000 0.3717 0.7500 -0.3783 -1.2000 0.4098 0.7692 -0.3594 -1.1000 0.4525 0.7885 -0.336 -1.0000 0.5000 0.8077 -0.3077 -0.9000 0.5525 0.8269 -0.2744 -0.8000 0.6098 0.8462 -0.2364 -0.7000 0.6711 0.8654 -0.1942 -0.6000 0.7353 0.8846 -0.1493 -0.5000 0.8000 0.9038 -0.1038 -0.4000 0.8621 0.9231 -0.061 -0.3000 0.9174 0.9423 -0.0249 -0.2000 0.9615 0.9615 0-0.1000 0.9901 0.9808 0.00930 1.0000 1.0000 0 0.1000 0.9901 0.9808 0.0093 0.2000 0.9615 0.9615 0 0.3000 0.9174 0.9423 -0.0249 0.4000 0.8621 0.9231 -0.061 0.5000 0.8000 0.9038 -0.1038 0.6000 0.7353 0.8846 -0.1493 0.7000 0.6711 0.8654 -0.1942 0.8000 0.6098 0.8462 -0.23640.9000 0.5525 0.8269 -0.27441.0000 0.5000 0.8077 -0.3077 1.1000 0.4525 0.7885 -0.336 1.2000 0.4098 0.7692 -0.3594 1.3000 0.3717 0.7500 -0.3783 1.4000 0.3378 0.7308 -0.3929 1.5000 0.3077 0.7115 -0.4038 1.6000 0.2809 0.6923 -0.4114 1.7000 0.2571 0.6731 -0.416 1.8000 0.2358 0.6538 -0.418 1.9000 0.2169 0.6346 -0.41772.0000 0.2000 0.6154 -0.41542.1000 0.1848 0.5962 -0.41132.2000 0.1712 0.5769 -0.40572.3000 0.1590 0.5577 -0.39872.4000 0.1479 0.5385 -0.39052.5000 0.1379 0.5192 -0.38132.6000 0.1289 0.5000 -0.37112.7000 0.1206 0.4808 -0.36012.8000 0.1131 0.4615 -0.34842.9000 0.1063 0.4423 -0.3363.0000 0.1000 0.4231 -0.32313.1000 0.0943 0.4038 -0.30963.2000 0.0890 0.3846 -0.29563.3000 0.0841 0.3654 -0.28133.4000 0.0796 0.3462 -0.26653.5000 0.0755 0.3269 -0.25153.6000 0.0716 0.3077 -0.23613.7000 0.0681 0.2885 -0.22043.8000 0.0648 0.2692 -0.20453.9000 0.0617 0.2500 -0.18834.0000 0.0588 0.2308 -0.17194.1000 0.0561 0.2115 -0.15544.2000 0.0536 0.1923 -0.13874.3000 0.0513 0.1731 -0.12184.4000 0.0491 0.1538 -0.10474.5000 0.0471 0.1346 -0.08764.6000 0.0451 0.1154 -0.07034.7000 0.0433 0.0962 -0.05284.8000 0.0416 0.0769 -0.03534.9000 0.0400 0.0577 -0.01775.0000 0.0385 0.0385 0n 2 3 4 5 6 7误差平方136.9209 79.1689 63.334 6.9775 23.7384 0.8329 和n 8 9 10 11 12 13 14 误差平方9.0015 0.5726 3.6152 0.572 1.5676 0.4648 0.7472和n 15 16 17 18 19 20 21误差平方0.3366 0.3945 0.2327 0.2291 0.1593 0.1438 0.1101和。