matlab实现数值分析插值及积分

合集下载

MATLAB数值实验一(数据的插值运算及其应用完整版)

MATLAB数值实验一(数据的插值运算及其应用完整版)

佛山科学技术学院实 验 报 告课程名称 数值分析 实验项目 插值法与数据拟合 专业班级 机械工程 姓 名 余红杰 学 号 10 指导教师 陈剑 成 绩 日 期 月 日一、实验目的1、学会Lagrange 插值、牛顿插值和三次样条插值等基本插值方法;2、讨论插值的Runge 现象3、学会Matlab 提供的插值函数的使用方法,会用这些函数解决实际问题。

二、实验原理1、拉格朗日插值多项式2、牛顿插值多项式3、三次样条插值 三、实验步骤1、用MATLAB 编写独立的拉格朗日插值多项式函数2、用MATLAB 编写独立的牛顿插值多项式函数3、用MATLAB 编写独立的三次样条函数(边界条件为第一、二种情形)4、已知函数在下列各点的值为:根据步骤1,2,3编好的程序,试分别用4次拉格朗日多项式4()L x 、牛顿插值多项式4()P x 以及三次样条函数()S x (自然边界条件)对数据进行插值,并用图给出 {(,),0.20.08,0,1,2,,10i i i x y x i i =+=},4()L x 、4()P x 和()S x 。

5、在区间[-1,1]上分别取10,20n =用两组等距节点对龙格函数21(),(11)125f x x x=-≤≤+作多项式插值,对不同n 值,分别画出插值函数及()f x 的图形。

6、下列数据点的插值可以得到平方根函数的近似,在区间[0,64]上作图。

(1)用这9个点作8次多项式插值8()L x 。

(2)用三次样条(第一边界条件)程序求()S x 。

7、对于给函数21()125f x x =+在区间[-1,1]上取10.2(0,1,,10)i x i i =-+=,试求3次曲线拟合,试画出拟合曲线并打印出方程,与第5题的结果比较。

四、实验过程与结果:1、Lagrange 插值多项式源代码:function ya=lag(x,y,xa) %x 所有已知插值点 %y 插值点对应函数值 %xa 所求点,自变量 %ya 所求点插值估计量 ya=0; mu=1; %初始化%循环方式求L 系数,并求和: for i = 1:length(y) for j = 1:length(x) if i ~= jmu = mu * (xa - x(j) ) / ( x(i) - x(j) ); else continue end endya = ya + y(i) * mu ; mu = 1; end2、Newton 源代码:function ya = newton(x,y,xa) %x 所有已知插值点 %y 插值点对应函数值 %xa 所求点,自变量 %ya 所求点插值估计量 %建立系数零矩阵D 及初始化:D = zeros(length(x)-1);ya = y(1);xi = 1;%求出矩阵D,该矩阵第一行为牛顿插值多项式系数:for i=1:(length(x)-1)D(i,1) = (y(i+1) -y(i))/(x(i+1) -x(i));endfor j=2:(length(x)-1)for i=1:(length(x)-j)D(i,j) = (D(i+1,j-1) - D(i,j-1)) / (x(i+j) - x(i)); endend%xi为单个多项式(x-x(1))(x-x(2))...的值for i=1:(length(x)-1)for j=1:ixi = xi*(xa - x(j));endya = ya + D(1,i)*xi;xi = 1;end3、三次样条插值多项式(1)(第一边界条件)源代码:function y=yt1(x0,y0,f_0,f_n,x) _____________(1)%第一类边界条件下三次样条插值;%xi 所求点;%yi 所求点函数值;%x 已知插值点;%y 已知插值点函数值;%f_0左端点一次导数值;%f_n右端点一次导数值;n = length(x0);z = length(y0);h = zeros(n-1,1);k=zeros(n-2,1);l=zeros(n-2,1);S=2*eye(n);for i=1:n-1h(i)= x0(i+1)-x0(i);endfor i=1:n-2k(i)= h(i+1)/(h(i+1)+h(i));l(i)= 1-k(i);end%对于第一种边界条件:k = [1;k]; _______________________(2)l = [l;1]; _______________________(3)%构建系数矩阵S:for i = 1:n-1S(i,i+1) = k(i);S(i+1,i) = l(i);end%建立均差表:F=zeros(n-1,2);for i = 1:n-1F(i,1) = (y0(i+1)-y0(i))/(x0(i+1)-x0(i));endD = zeros(n-2,1);for i = 1:n-2F(i,2) = (F(i+1,1)-F(i,1))/(x0(i+2)-x0(i));D(i,1) = 6 * F(i,2);end%构建函数D:d0 = 6*(F(1,2)-f_0)/h(1); ___________(4)dn = 6*(f_n-F(n-1,2))/h(n-1); ___________(5)D = [d0;D;dn]; ______________(6)m= S\D;%寻找x所在位置,并求出对应插值:for i = 1:length(x)for j = 1:n-1if (x(i)<=x0(j+1))&(x(i)>=x0(j))y(i) =( m(j)*(x0(j+1)-x(i))^3)/(6*h(j))+...(m(j+1)*(x(i)-x0(j))^3)/(6*h(j))+...(y0(j)-(m(j)*h(j)^2)/6)*(x0(j+1)-x(i))/h(j)+... (y0(j+1)-(m(j+1)*h(j)^2)/6)*(x(i)-x0(j))/h(j) ; break;else continue;endendend(2)(自然边界条件)源代码:仅仅需要对上面部分标注的位置做如下修改:__(1):function y=yt2(x0,y0,x)__(2):k=[0;k]__(3):l=[l;0]__(4)+(5):删除—(6):D=[0:D:0]4、——————————————PS:另建了一个f方程文件,后面有一题也有用到。

插值MATLAB实现(牛顿差商插值误差龙格现象切比雪夫插值)

插值MATLAB实现(牛顿差商插值误差龙格现象切比雪夫插值)

插值MATLAB实现(牛顿差商插值误差龙格现象切比雪夫插值)插值是数值分析中的一种方法,通过已知数据点的函数值来估计函数在其他点的值。

MATLAB提供了多种方法来实现插值,包括牛顿差商插值、插值误差分析、龙格现象和切比雪夫插值。

下面将详细介绍这些方法的实现原理和MATLAB代码示例。

1.牛顿差商插值:牛顿差商插值是一种基于多项式插值的方法,其中差商是一个连续性的差分商。

该方法的优势在于可以快速计算多项式的系数。

以下是MATLAB代码示例:```matlabfunction [coeff] = newton_interpolation(x, y)n = length(x);F = zeros(n, n);F(:,1)=y';for j = 2:nfor i = j:nF(i,j)=(F(i,j-1)-F(i-1,j-1))/(x(i)-x(i-j+1));endendcoeff = F(n, :);end```该代码中,输入参数x和y分别表示已知数据点的x坐标和y坐标,返回值coeff表示插值多项式的系数。

2.插值误差分析:插值误差是指插值函数与原始函数之间的差异。

一般来说,通过增加插值节点的数量或使用更高次的插值多项式可以减小插值误差。

以下是MATLAB代码示例:```matlabfunction [error] = interpolation_error(x, y, x_eval)n = length(x);p = polyfit(x, y, n-1);y_eval = polyval(p, x_eval);f_eval = sin(pi*x_eval);error = abs(f_eval - y_eval);end```该代码中,输入参数x和y分别表示已知数据点的x坐标和y坐标,x_eval表示插值节点的x坐标,error表示插值误差。

3.龙格现象:龙格现象是插值多项式在等距插值节点上错误增长的现象。

基于Matlab的数值积分公式问题

基于Matlab的数值积分公式问题

数值分析学号:130080402015学生所在学院:测试与光电工程学院学生姓名:张翀任课教师:郑华盛教师所在学院:数信学院基于Matlab的数值积分公式问题张翀,测试与光电工程学院测试计量技术及仪器,130080402015摘要:在求一些函数的定积分时,由于原函数十分复杂难以求出或用初等函数表达,导致积分很难精确求出,只能设法求其近似值,因此能够直接借助牛顿----莱布尼兹公式计算定积分的情形是不多的。

数值积分就是解决此类问题的一种行之有效的方法。

积分的数值计算是数值分析的一个重要分支;因此,探讨近似计算的数值积分方法是有着明显的实际意义的。

本文介绍了数值积分法的几种计算公式,如矩形求积公式、梯形求积公式和辛普森求积公式及相应的MATLAB命令,并给出了用 MATLAB编程求数值积分的实例。

关键词: MA TLAB;数值积分;矩形求积公式;梯形求积公式;辛普森求积公式目录1引言 (1)2数值积分算法介绍 (1)2.1数值求积公式的构造 (1)2.2求积公式的推导 (2)2.3常见的牛顿-科特斯求积公式 (5)2.4复合求积公式 (7)3关于河流横断面积的数值积分问题 (8)4问题的求解过程 (9)5基于MATLAB编程的各种求积公式对问题的求解 (9)6总结 (13)参考文献 (14)附录 (15)1 引言实际问题当中常常需要计算积分。

有些数值方法,如微分方程和积分方程的求解,也都和积分计算相联系。

在一元微积分学中,对于积分 ⎰=badx x f I)(,只要找到被积函数f (x )原函数为F( x) ,求f(x)在该区间上的定积分便可用牛顿 - 莱布尼兹公式求解,即⎰-=baa Fb F dx x f )()()(。

用牛顿 - 莱布尼兹公式计算定积分的方法在理论上和解决实际问题中起到了很大的作用 ,但它并不能解决定积分计算的所有问题。

在工程技术领域常遇到十分复杂的情况而无法用牛顿 - 莱布尼兹公式求解.其可能出现的情况有:(1) 某些被积函数f(x),其原函数无法用初等函数表示 ,如⎰dx e x 2, dx xx⎰sin 等。

第十章MATLAB的数值分析

第十章MATLAB的数值分析

• 第一个问题可归结为“已知函数在x0,x1,
– …,xn处的值,求函数在区间[x0,xn]内其它点处的值”,这 种问题适宜用插值方法解决。 – 插值问题可描述为:已知函数在x0,x1,…,xn处的值 y0,y1,…,yn,求函数p(x),使p(xi) = yi。
• 但对第二个问题不宜用插值方法,因为600米已超出所 给数据范围,用插值函数外推插值区间外的数据会 产生较大的误差。
– Q1=prctile(w,25); – Q3=prctile(w,75); – prctile( )函数实现计算样本的百分位数功能
分布形态的测定
• 只用集中趋势和离中趋势来表示所有数据,难免不 够准确。分析总体次数的分布形态有助于识别整个 总体的数量特征。总体的分布形态可以从两个角度 考虑,一是分布的对称程度,另一个是分布的高低。 前者的测定参数称为偏度或偏斜度,后者的测定参 数称为峰度。 • 峰度是掌握分布形态的另一指标,它能描述分布的 平缓或陡峭程度。如果峰度数值等于零,说明分布 为正态;若峰度数值大于零,说明分布呈陡峭状态; 若峰度数值小于零,说明分布形态趋于平缓。
– 解决第二个问题的常用方法是,根据地面到井下 500 处的 数据求出瓦斯浓度与地面到井下距离x之间的近似函数关 系f(x), 由f(x)求井下600米处的瓦斯浓度。
• 插值函数过已知点,拟合函数不一定过已知点。通 常, 插值主要用于求函数值,而拟合的主要目的是求 函数关系。当然,某些问题既可以用插值也可以用 拟合。
插值方法-概述
• 为什么需要插值?
(1) 函数关系y=f(x)没有明确的表达式
(2) y=f(x)表达式复杂,不便于研究和使用
-20 -15
沉陷量/mm 下沉方向为"+"

数值分析的matlab实现

数值分析的matlab实现

第2章牛顿插值法实现参考文献:[1]岑宝俊. 牛顿插值法在凸轮曲线修正设计中的应用[J]. 机械工程师,2009,10:54-55.求牛顿插值多项式和差商的MA TLAB 主程序:function[A,C,L,wcgs,Cw]=newpoly(X,Y)n=length(X);A=zeros(n,n);A(:,1) =Y';s=0.0;p=1.0;q=1.0;c1=1.0;for j=2:nfor i=j:nA(i,j)=(A(i,j-1)-A(i-1,j-1))/(X(i)-X(i-j+1));endb=poly(X(j-1));q1=conv(q,b);c1=c1*j;q=q1;endC=A(n,n);b=poly(X(n));q1=conv(q1,b);for k=(n-1):-1:1C=conv(C,poly(X(k)));d=length(C);C(d)=C(d)+A(k,k);endL(k,:)=poly2sym(C);Q=poly2sym(q1);syms Mwcgs=M*Q/c1;Cw=q1/c1;(1)保存名为newpoly.m 的M 文件(2)输入MA TLAB 程序>> X=[242,243,249,250];>> Y=[13.681,13.526,13.098,13.095];>> [A,C,L,wcgs,Cw]=newpoly(X,Y)输出3阶牛顿插值多项式L 及其系数向量C 差商的矩阵A ,插值余项wcgs 及其)()()1(ξ+n n f x R 的系数向量Cw 。

A =13.6810 0 0 013.5260 -0.1550 0 013.0980 -0.0713 0.0120 013.0950 -0.0030 0.0098 -0.0003C =1.0e+003 *-0.0000 0.0002 -0.0551 4.7634L =- (23*x^3)/84000 + (2981*x^2)/14000 - (7757472138947345*x)/140737488355328 + 5237382665812919/1099511627776wcgs =(M*(x^4 - 984*x^3 + 363071*x^2 - 59535444*x + 3660673500))/24Cw =1.0e+008 *0.0000 -0.0000 0.0002 -0.0248 1.5253输入MATLAB程序>> x=244;>> y=- (23*x^3)/84000 + (2981*x^2)/14000 - (7757472138947345*x)/140737488355328 + 5237382665812919/1099511627776y =13.3976输入MATLAB程序>> x=[244,245,246,247,248];>> y=- (23.*x.^3)./84000 + (2981.*x.^2)./14000 - (7757472138947345.*x)./140737488355328 + 5237382665812919./1099511627776y =13.3976 13.2943 13.2143 13.1560 13.1178第4章 高斯-勒让德积分公式实现用高斯-勒让德积分公式计算dx e I x ⎰--=112221π,取代数精度为3和5,再根据截断误差公式写出误差公式,并将计算结果与精确值进行比较。

在Matlab中如何进行数据插值与拟合

在Matlab中如何进行数据插值与拟合

在Matlab中如何进行数据插值与拟合引言:数据处理是科学研究与工程开发中不可或缺的环节之一。

而数据插值和拟合则是数据处理中常用的技术手段。

在Matlab这一强大的数值分析工具中,提供了丰富的函数与工具箱,使得数据插值与拟合变得更加便捷高效。

本文将详细阐述在Matlab中如何进行数据插值与拟合,并介绍几个常用的插值与拟合方法。

一、数据插值数据插值是通过已知的有限个数据点,推导出数据点之间未知位置上的数值。

在Matlab中,可以利用interp1函数进行数据插值。

假设我们有一组离散的数据点,存储为两个向量x和y。

那么,可以通过以下步骤进行数据插值:1. 调用interp1函数,并传入x和y作为输入参数。

```matlabxi = linspace(min(x), max(x), n);yi = interp1(x, y, xi, '方法');```其中,xi是插值点的位置,min和max分别是x向量的最小值和最大值,n是插值点的数量。

'方法'是要使用的插值方法,可以选择线性插值(method='linear')、样条插值(method='spline')等。

2. 绘制插值结果曲线。

```matlabplot(x, y, 'o', xi, yi)legend('原始数据','插值结果')```使用plot函数可以绘制原始数据点和插值结果的曲线。

通过设置不同的插值方法和插值点的数量,可以探索不同的插值效果。

二、数据拟合数据拟合是通过已知的一组数据点,找到一个符合数据趋势的函数模型。

在Matlab中,可以利用polyfit函数进行多项式拟合。

假设我们有一组离散的数据点,存储为两个向量x和y。

那么,可以通过以下步骤进行数据拟合:1. 调用polyfit函数,并传入x和y作为输入参数。

```matlabp = polyfit(x, y, n);```其中,n是多项式的次数,p是拟合多项式的系数。

数值分析实验报告matlab

数值分析实验报告matlab数值分析实验报告引言:数值分析是一门研究利用计算机数值方法解决数学问题的学科,它在科学计算、工程设计、金融分析等领域具有重要的应用价值。

本实验报告旨在通过使用MATLAB软件,探索数值分析的基本原理和方法,并通过实际案例加深对数值分析的理解。

一、误差分析在数值计算中,误差是无法避免的。

误差分析是数值分析中的重要一环,它帮助我们了解数值计算的准确性和稳定性。

在实验中,我们通过计算机模拟了一个简单的数学问题,并分别计算了绝对误差和相对误差。

通过比较不同算法的误差大小,我们可以选择最适合的算法来解决实际问题。

二、插值与拟合插值和拟合是数值分析中常用的方法,它们可以通过已知的数据点来推导出未知数据点的近似值。

在本实验中,我们通过MATLAB的插值函数和拟合函数,分别进行了插值和拟合的实验。

通过比较不同插值和拟合方法的结果,我们可以选择最适合的方法来处理实际问题。

三、数值积分数值积分是数值分析中的重要内容,它可以用来计算曲线下的面积或函数的积分值。

在实验中,我们通过MATLAB的数值积分函数,对一些简单的函数进行了积分计算。

通过比较数值积分和解析积分的结果,我们可以评估数值积分的准确性和稳定性,并选择最适合的积分方法来解决实际问题。

四、常微分方程的数值解法常微分方程是数值分析中的重要内容,它可以用来描述许多自然现象和工程问题。

在实验中,我们通过MATLAB的常微分方程求解函数,对一些简单的微分方程进行了数值解法的计算。

通过比较数值解和解析解的结果,我们可以评估数值解法的准确性和稳定性,并选择最适合的数值解法来解决实际问题。

五、线性方程组的数值解法线性方程组是数值分析中的经典问题,它在科学计算和工程设计中广泛应用。

在实验中,我们通过MATLAB的线性方程组求解函数,对一些简单的线性方程组进行了数值解法的计算。

通过比较数值解和解析解的结果,我们可以评估数值解法的准确性和稳定性,并选择最适合的数值解法来解决实际问题。

MATLAB数值运算.pdf

第3章 MATLAB 数值运算教学提示:每当难以对一个函数进行积分或者微分以确定一些特殊的值时,可以借助计算机在数值上近似所需的结果,从而生成其他方法无法求解的问题的近似解。

这在计算机科学和数学领域,称为数值分析。

本章涉及的数值分析的主要内容有插值与多项式拟合、数值微积分、线性方程组的数值求解、微分方程的求解等,掌握这些主要内容及相应的基本算法有助于分析、理解、改进甚至构造新的数值算法。

教学要求:本章主要是让学生掌握数值分析中多项式插值和拟合、牛顿-科茨系列数值求积公式、3种迭代方法求解线性方程组、解常微分方程的欧拉法和龙格-库塔法等具体的数值算法,并要求这些数值算法能在MATLAB 中实现。

3.1 多 项 式在工程及科学分析上,多项式常被用来模拟一个物理现象的解析函数。

之所以采用多项式,是因为它很容易计算,多项式运算是数学中最基本的运算之一。

在高等数学中,多项式一般可表示为以下形式:120121()n n n n n f x a x a x a x a x a −−−=+++++…。

当x 是矩阵形式时,代表矩阵多项式,矩阵多项式是矩阵分析的一个重要组成部分,也是控制论和系统工程的一个重要工具。

3.1.1 多项式的表达和创建在MATLAB 中,多项式表示成向量的形式,它的系数是按降序排列的。

只需将按降幂次序的多项式的每个系数填入向量中,就可以在MATLAB 中建立一个多项式。

例如,多项式43231529s s s s +−−+在MATLAB 中,按下面方式组成一个向量x = [1 3 -15 -2 9]MATLAB 会将长度为n +1的向量解释成一个n 阶多项式。

因此,若多项式某些项系数为零,则必须在向量中相应位置补零。

例如多项式41s +在MATLAB 环境下表示为y = [1 0 0 0 1]3.1.2 多项式的四则运算多项式的四则运算包括多项式的加、减、乘、除运算。

下面以对两个同阶次多项式MATLAB 基础及其应用教程·66··66·32()234a x x x x =+++,32()4916b x x x x =+++做加减乘除运算为例,说明多项式的四则运算过程。

如何在Matlab中进行数值积分和数值解

如何在Matlab中进行数值积分和数值解在数学和工程领域,数值积分和数值解是常见的技术手段,可以帮助我们求解复杂的数学问题和实际工程中的模型。

本文将介绍如何使用Matlab进行数值积分和数值解,以及一些注意事项和常用的方法。

一、数值积分数值积分是计算定积分的近似值的方法,可以通过数值逼近或数值插值来实现。

在Matlab中,有几种常用的函数可以用于数值积分,比如trapz、quad等。

1. trapz函数trapz函数是用梯形法则计算积分的函数。

它的使用方法是将要积分的函数作为输入的第一个参数,x轴上的点作为输入的第二个参数。

例如,要计算函数f(x)在区间[a, b]上的积分,可以使用以下代码:result = trapz(x, f(x));2. quad函数quad函数是使用自适应数值积分算法计算积分的函数。

它的使用方法是将要积分的函数作为输入的第一个参数,积分区间的下限和上限作为输入的第二个和第三个参数。

例如,要计算函数f(x)在区间[a, b]上的积分,可以使用以下代码:result = quad(@(x) f(x), a, b);二、数值解数值解是使用数值方法求解复杂的数学问题或实际工程中的模型的近似解。

在Matlab中,有几种常用的函数可以用于数值解,比如fsolve、ode45等。

1. fsolve函数fsolve函数是用于求解非线性方程组的函数。

它的使用方法是将非线性方程组表示为一个函数,然后将该函数作为输入的第一个参数。

例如,要求解方程组f(x) = 0,可以使用以下代码:x = fsolve(@(x) f(x), x0);其中x0是方程的初始猜测值。

2. ode45函数ode45函数是求解常微分方程初值问题的函数。

它的使用方法是将微分方程表示为一个函数,然后将该函数作为输入的第一个参数。

例如,要求解常微分方程dy/dx = f(x, y),可以使用以下代码:[t, y] = ode45(@(t, y) f(t, y), tspan, y0);其中tspan是时间区间,y0是初始条件。

matlab插值法

样条插值
样条插值是一种分段插值方法,它在每个小区间上使用低 次多项式进行插值,同时保证整个插值函数的连续性和光 滑性。
MATLAB中实现插值法
MATLAB提供的插值函数
MATLAB提供了多种内置函数来实现不同类型的插值,如`interp1`、`interp2`、`interp3` 等,分别用于一维、二维和三维数据的插值。
03
二维数据插值方法
网格数据插值
线性插值
基于已知网格点上的数据,通过 线性插值方法计算未知点的值。 这种方法简单快速,但可能不够 精确。
双三次插值
使用周围的16个网格点上的数据 ,构建一个双三次多项式来逼近 未知点的值。这种方法比线性插 值更精确,但计算量较大。
散点数据插值
最近邻插值
将未知点的值设置为距离其最近的已知点的值。这种方法简 单快速,但可能导致不连续的结果。
信号调制与解调中应用
信号调制
在通信系统中,插值法可用于实现信号的调制处理,将基带信号 转换为适合在信道中传输的已调信号。
信号解调
接收端在接收到已调信号后,可以使用插值法对信号进行解调处理 ,还原出原始的基带信号。
符号同步与定时恢复
在数字通信系统中,插值法可用于实现符号同步和定时恢复,确保 接收端能够准确地提取出传输的符号信息。
07
总结与展望
回顾本次课程重点内容
插值法基本概念
插值法是一种通过已知数据点估算未知数据点的方法,广泛应用于 数据分析和科学计算领域。
MATLAB插值法实现
通过MATLAB提供的插值函数,如`interp1`、`interp2`、`interp3` 等,实现一维、二维和三维数据的插值计算。
插值法应用场景
图像修复与增强中应用
  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

Matlab实现数值分析插值及积分摘要:数值分析(numerical analysis)是研究分析用计算机求解数学计算问题的数值计算方法及其理论的学科,是数学的一个分支,它以数字计算机求解数学问题的理论和方法为研究对象。

在实际生产实践中,常常将实际问题转化为数学模型来解决,这个过程就是数学建模。

学习数值分析这门课程可以让我们学到很多的数学建模方法。

分别运用matlab数学软件编程来解决插值问题和数值积分问题。

题目中的要计算差值和积分,对于问题一,可以分别利用朗格朗日插值公式,牛顿插值公式,埃特金逐次线性插值公式来进行编程求解,具体matlab代码见正文。

编程求解出来的结果为:=+。

其中Aitken插值计算的结果图如下:对于问题二,可以分别利用复化梯形公式,复化的辛卜生公式,复化的柯特斯公式编写程序来进行求解,具体matlab代码见正文。

编程求解出来的结果为: 0.6932其中复化梯形公式计算的结果图如下:问题重述问题一:已知列表函数表格 1分别用拉格朗日,牛顿,埃特金插值方法计算。

问题二:用复化的梯形公式,复化的辛卜生公式,复化的柯特斯公式计算积分,使精度小于5。

问题解决问题一:插值方法对于问题一,用三种差值方法:拉格朗日,牛顿,埃特金差值方法来解决。

一、拉格朗日插值法:拉格朗日插值多项式如下:首先构造1+n 个插值节点n x x x ,,,10 上的n 插值基函数,对任一点i x 所对应的插值基函数)(x l i ,由于在所有),,1,1,,1,0(n i i j x j +-=取零值,因此)(x l i 有因子)())(()(110n i i x x x x x x x x ----+- 。

又因)(x l i 是一个次数不超过n 的多项式,所以只可能相差一个常数因子,固)(x l i 可表示成:)())(()()(110n i i i x x x x x x x x A x l ----=+-利用1)(=i i x l 得:)())(()(1110n i i i i i i x x x x x x x x A ----=+-于是),,2,1,0()())(()()())(()()(110110n i x x x x x x x x x x x x x x x x x l n i i i i i i n i i i =--------=+-+-因此满足i i n y x L =)( ),2,1,0(n i =的插值多项式可表示为:∑==nj j j n x l y x L 0)()(从而n 次拉格朗日插值多项式为:),,2,1,0()()(0n i x l y x L nj i j j i n ==∑=matlab 编程:编程思想:主要从上述朗格朗日公式入手:依靠循环,运用poly ()函数和conv ()函数表示拉格朗日公式,其中的poly (i )函数表示以i 作为根的多项式的系数,例如poly (1)表示x-1的系数,输出为1 -1,而poly (poly (1))表示(x-1)*(x-1)=x^2-2*x+1的系数,输出为1 -2 1;而conv ()表示多项式系数乘积的结果,例如conv (poly (1),poly (1))输出为1 -2 1;所以程序最后结果为x^n+x^n-1+……+x^2+x+1(n 的值据结果的长度为准)的对应系数。

在命令窗口输入edit lagran 来建立lagran.m 文件,文件中的程序如下: function [c,l]=lagran(x,y) w=length(x); n=w-1;l=zeros(w,w); for k=1:n+1 v=1;for j=1:n+1 if k~=jv=conv(v,poly(x(j)))/(x(k)-x(j)); end endl(k,:)=v; endc=y*l;输入:>> x=[0 1 2 3 4]; >> y=[1 2 17 82 257]; >> lagran(x,y)运行结果为 ans =1.0000 -0.0000 -0.0000 0 1.0000 结果为:=+。

如图表1:图表 1二.牛顿插值法newton 插值多项式的表达式如下:010011()()()()()n n n N x c c x x c x x x x x x -=+-+⋅⋅⋅+--⋅⋅⋅-其中每一项的系数ci 的表达式如下:12011010[,,,][,,,][,,,]i i i i i f x x x f x x x c f x x x x x -⋅⋅⋅-⋅⋅⋅=⋅⋅⋅=-即为 f (x)在点01,,,i x x x ⋅⋅⋅处的i 阶差商,([]()i i f x f x =,1,2,,i n =),由差商01[,,,]i f x x x ⋅⋅⋅的性质可知:()01001[,,,]()iii j j k j k k jf x x x f x x x ==≠⋅⋅⋅=-∑∏matlab 编程:编程思想:主要从上述牛顿插值公式入手:依靠循环,运用poly ()函数和conv ()函数表示拉格朗日公式,其中的poly (i )函数表示以i 作为根的多项式的系数,例如poly (1)表示x-1的系数,输出为1 -1,而poly(poly(1))表示(x-1)*(x-1)=x^2-2*x+1的系数,输出为1 -2 1;而conv()表示多项式系数乘积的结果,例如conv(poly(1),poly(1))输出为1 -2 1;所以程序最后结果为x^n+x^n-1+……+x^2+x+1(n的值据结果的长度为准)的对应系数。

在命令窗口输入edit nowpoly来建立newpoly.m文件,文件中的程序如下:function [c,d]=newpoly(x,y)n=length(x);d=zeros(n,n);d(:,1)=y';for j=2:nfor k=j:nd(k,j)=(d(k,j-1)-d(k-1,j-1))/(x(k)-x(k-j+1));endendc=d(n,n);for k=(n-1):-1:1c=conv(c,poly(x(k)));m=length(c);c(m)=c(m)+d(k,k);end输入:>> x=[0 1 2 3 4];>> y=[1 2 17 82 257];>> newpoly(x,y)运行结果为ans =1 0 0 0 1所以=+。

如图表2:图表 2三.埃特金插值法:Aitken插值公式如下:递推表达式为:= +当n=1时,= +当n=2时,= +其中的带入递推表达式求得。

由此递推下去,最终得到的结果。

matlab编程:编程思想:埃特金插值多项式又称作Aitken逐次线性插值多项式,根据公式的特点,可以利用2次嵌套循环将公式表示出来。

在命令窗口输入edit Aitken 来建立Aitken.m文件,文件中的程序如下:function f = Aitken(x,y)syms z;n = length(x);y1(1:n) = z;for i=1:n-1for j=i+1:ny1(j) = y(j)*(z-x(i))/(x(j)-x(i))+y(i)*(z-x(j))/(x(i)-x(j));endy = y1;simplify(y1);endsimplify(y1(n));f = collect(y1(n));输入:>> x=[0 1 2 3 4];>> y=[1 2 17 82 257];>> Aitken(x,y)运行结果为ans =z^4 + 1所以=+。

如图表3:图表 3问题二:复化积分对于问题二来说,用复化的梯形公式,复化的辛卜生公式,复化的柯特斯公式结局问题(计算积分,使精度小于5)。

一复化的梯形公式:复化梯形的迭代公式为:;matlab编程:程序1(求f(x)的n阶导数:)在命令窗口输入edit qiudao 来建立qiudao.m文件,文件中的程序如下:function d=qiudao(x,n)syms x;f=1/x;n=input('输入导数阶数: ');d=diff(f,x,n);输入:qiudao(x,n)输入所求导数阶数:2显示:n = 2ans =2/x^3结果为:f2 =2/x^3如图表4:图表 4程序2:在命令窗口输入edit tixing 来建立tixing.m文件,文件中的程序如下:function y=tixing()syms x ; %定义自变量xf=inline('1/x','x');%定义函数f(x)= 1/xf2=inline('2/x^3','x') ;%定义f(x)的二阶导数,输入程序1里求出的f2即可。

f3='-2/x^3';%因fminbnd()函数求的是表达式的最小值,且要求表达式带引号,故取负号,以便求最大值e=5*10^(-5); %精度要求值a=1;%积分下限b=2;%积分上限x1=fminbnd(f3,1,2); %求负的二阶导数的最小值点,也就是求二阶导数的最大值点对应的x值for n=2:1000000;%求等分数nRn=-(b-a)/12*((b-a)/n)^2*f2(x1); %计算余项if abs(Rn)<e %用余项进行判断break % 符合要求时结束endendh=(b-a)/n; %求hTn1=0;for k=1:n-1 %求连加和xk=a+k*h;Tn1=Tn1+f(xk);endTn=h/2*((f(a)+2*Tn1+f(b)));fprintf('用复化梯形算法计算的结果 Tn=')disp(Tn)fprintf('等分数 n=')disp(n) %输出等分数输入:tixing()运行结果为:用复化梯形计算的结果 Tn= 0.6932等分数 n= 58结果如图表:5图表 5二复化的辛卜生公式:复化simpson迭代公式为:;matlab编程:程序1(求f(x)的n阶导数:)在命令窗口输入edit qiudao 来建立qiudao.m文件,文件中的程序如下:function d=qiudao(x,n)syms x;f=1/x;n=input('输入导数阶数: ');d=diff(f,x,n);输入:qiudao(x,n)输入所求导数阶数:4显示:n = 4ans =24/x^5结果为:f2 = 24/x^5如图表6:图表 6程序2:在命令窗口输入edit xinpusheng 来建立xinpusheng.m文件,文件中的程序如下:function y=xinpusheng()syms x ;f=inline('1/x','x');f2=inline('24/x^5','x');f3='-24/x^5';e=5*10^(-5);a=1;b=2;x1=fminbnd(f3,1,2);for n=2:1000000Rn=-(b-a)/180*((b-a)/(2*n))^4*f2(x1);if abs(Rn)<ebreakendendh=(b-a)/n;Sn1=0;Sn2=0;for k=0:n-1xk=a+k*h;xk1=xk+h/2;Sn1=Sn1+f(xk1);Sn2=Sn2+f(xk);endSn=h/6*(f(a)+4*Sn1+2*(Sn2-f(a))+f(b));fprintf('用Simpson公式计算的结果 Sn=')disp(Sn)fprintf('等分数 n=')disp(n)调用xinpusheng.m中xinpusheng函数:输入:>> clear>> xinpusheng()运行结果为用Simpson 公式计算的结果 Sn= 0.6932等分数 n= 4结果如图表7图表 7三 复化的柯特斯公式:牛顿-柯特斯公式如下:matlab 编程:(1)在命令窗口输入edit NewtonCotes 来建立NewtonCotes.m 文件,文件中的程序如下: function [y,Ck,Ak]=NewtonCotes(fun,a,b,n)if nargin==1[mm,nn]=size(fun);if mm>=8error('为了保证NewtonCotes 积分的稳定性,最多只能有9个等距节点!') elseif nn~=2error('fun 构成应为:第一列为x ,第二列为y ,并且个数为小于10的等距节点!')endxk=fun(1,:);fk=fun(2,:);a=min(xk);0()()n b k k a k f x dx A f x =≈∑⎰011011()()()()()()()()()b b k k n k k a a k k k k k k n x x x x x x x x A l x dx dx x x x x x x x x -+-+----==----⎰⎰b=max(xk);n=mm-1;elseif nargin==4xk=linspace(a,b,n+1);if isa(fun,'function_handle')fx=fun(xk);elseerror('fun积分函数的句柄,且必须能够接受矢量输入!')endelseerror('输入参数错误,请参考函数帮助!')endCk=cotescoeff(n);Ak=(b-a)*Ck;y=Ak*fx';(2)在命令窗口输入edit cotescoeff来建立cotescoeff.m文件,文件中的程序如下:function Ck=cotescoeff(n)for i=1:n+1k=i-1;Ck(i)=(-1)^(n-k)/factorial(k)/factorial(n-k)/n*quadl((t)intfun(t,n,k),0,n); end(3)在命令窗口输入edit intfun来建立intfun.m文件,文件中的程序如下:function f=intfun(t,n,k)f=1;for i=[0:k-1,k+1:n]f=f.*(t-i);end% fun,积分表达式,这里有两种选择%(1)积分函数句柄,必须能够接受矢量输入,比如fun=(x)1./x% (2)x,y坐标的离散点,第一列为x,第二列为y,必须等距,且节点的个数小于9,比如: fun=[1:8;1./(1:8)]% 如果fun的表采用第二种方式,那么只需要输入第一个参数即可,否则还要输入a,b,n 三个参数% a,积分下限% b,积分上限% n,牛顿-科特斯数公式的阶数,必须满足1<n<7,因为n>=8时不能保证公式的稳定性% (1)n=1,即梯形公式% (2)n=2,即辛普森公式% (3)n=4,即科特斯公式在显示窗口输入:fun=(x)1./x;a=-1;b=1;n=4;NewtonCotes(fun,a,b,n)运行结果为:ans =0.6932结果如图表8图表 8。

相关文档
最新文档