Maab数值积分与数值微分
(完整版)(数值积分与微分的MATLAB命令)
0.9008 n=
77
例题2求定积分
1 ex2 dx
1
>>quad(inline(‘exp(-x.^2)’),-1,1,0.5e-4) 或者 quad(@(x)exp(-x.^2),-1,1,0.5e-4)
7.triplequad 用于计算长方体区域上的三重积分
例 计算
dy
1 4xzex2 yz2 dzdydx
0
00
>>fun=@(x,y,z)4*x.*z.*exp(-x.^2-y.*z.^2)
>>S=triplequad(fun,0,pi,0,pi,0,1)
8.int 符号函数积分
例6 x=-1:0.1:1; y=exp(-x.^2); trapz(x,y)
4.自适应Lobatto方法 quadl命令适用于光滑函数， quadl调用被积函
数的次数明显少于quad命令，而且精度比quad函 数更高，其调用格式与quad一样 [I,n]=quad('fname',a,b,tol,trace)
8.1.2 数值积分的实现方法 1．变步长辛普生法(自适应simpson积分法) 基于变步长辛普生法，MATLAB给出了quad函数
来求定积分。该函数的调用格式为：
[I,n]=quad('fname',a,b,tol,trace)
其中fname是被积函数名。a和b分别是定积分的下 限和上限。tol用来控制积分精度，缺省时取
2．牛顿－柯特斯法 基于牛顿－柯特斯法，MATLAB给出了
第7章 MATLAB数值积分与数值微分PPT教学课件
MATLAB数值积分 MATLAB数值微分
2020/12/10
1
7.1 数值积分
7.1.1 数值积分基本原理 求解定积分的数值方法多种多样，如简单方法：
梯形法、辛普生(Simpson)法、牛顿－柯特斯 (Newton-Cotes)法等方法。它们的基本思想都是 将整个积分区间[a,b]分成n个子区间[xi,xi+1]， i=1,2,…,n，其中x1=a，xn+1=b。求定积分问题就分 解为求和问题。
其中，fname是被积函数名。a和b分别是定积分的下 限和上限。tol用来控制积分精度，默认时取tol=10-6。 trace控制是否展现积分过程，若取非0则展现积分过 程，取0则不展现，默认时取trace=0。返回参数I即 定积分值，n为被积函数的调用次数。
2020/12/10
3
例7-1 求定积分。
2020/12/10
10
例7-4 用trapz函数计算定积分。 命令如下：
X=1:0.01:2.5; Y=exp(-X); %生成函数关系数据向量 trapz(X,Y) ans =
0.2858
2020/12/10
Байду номын сангаас11
7.1.3 多重定积分的数值求解
MATLAB提供的dblquad函数用于求二重积分的数值 解，triplequad函数用于求三重积分的数值解。函 数的调用格式为
2020/12/10
2
7.1.2 数值积分的实现方法
1．变步长辛普生法
基于变步长辛普生法，MATLAB给出了quad函数和 quadl函数来求定积分。函数的调用格式为：
[I,n]=quad(@fname,a,b,tol,trace)
Matlab中常用的积分和微分算法解析
Matlab中常用的积分和微分算法解析积分和微分是数学中重要的概念和工具，广泛应用于科学、工程和计算领域。
在Matlab中，提供了丰富的积分和微分算法，可以方便地进行数值计算和符号计算。
在本文中，我们将解析Matlab中常用的积分和微分算法，并探讨其应用。
一、数值积分算法数值积分是通过将求和转化为积分的方式，对函数在一定区间内的近似计算。
在Matlab中，有许多数值积分算法可供选择，包括梯形法则、辛普森法则和高斯求积法等。
1. 梯形法则梯形法则是一种基本的数值积分算法。
它将区间分成多个小梯形，并将每个小梯形的面积近似表示为梯形的面积，然后将这些面积相加得到最终的近似积分值。
在Matlab中，可以使用trapz函数来实现梯形法则的计算。
例如，对函数f(x)在区间[a, b]上进行积分，可以使用如下代码：```matlaba = 0;b = 1;x = linspace(a, b, 100);y = f(x);integral_value = trapz(x, y);```其中，linspace函数用于生成均匀分布的点，f(x)是待积分的函数。
trapz函数可以自动计算积分值。
2. 辛普森法则辛普森法则是一种更精确的数值积分算法。
它将区间分成多个小三角形，并将每个小三角形的面积近似表示为一个带有二次多项式的面积，然后将这些面积相加得到最终的近似积分值。
在Matlab中，可以使用quad函数来实现辛普森法则的计算。
例如，对函数f(x)在区间[a, b]上进行积分，可以使用如下代码：```matlaba = 0;b = 1;integral_value = quad(@f, a, b);```其中，@f表示函数句柄，quad函数可以自动计算积分值。
3. 高斯求积法高斯求积法是一种更高精度的数值积分算法。
它利用多个节点和权重，通过插值的方式来近似积分值。
在Matlab中，可以使用gaussquad函数来实现高斯求积法的计算。
MATLAB_简介(5)MATLAB数值积分与微分
二重积分dblquad()与三重积分
0
2
[ y sin(x ) x cos(y )]dxdy
fun=inline('y*sin(x)+x*cos(y)')
Q=dblquad(fun,pi,2*pi,0,pi) Q = -9.8698 [x,y]=meshgrid(pi:.1:2*pi,0:.1:pi); z=fun(x,y); mesh(x,y,z)
>> plot(x,y,'o',x,y)
>> title('y(x) data plot') >> ylabel('y(x)'), xlabel('x') >> dy=diff(y)./diff(x); >> xd=x(1:length(x)-1); >> plot(xd,dy) >> title('Approximate derivative using diff') >> ylabel('dy/dx'), xlabel('x')
注意二者皆以后向差分计算且数据点只剩 4 个 而不是5个。而 dy/dx 的 数值微分则为 dy/dx=diff(y)./diff(x)。 因此要计算下列多项式在 [-4, 5] 区间的微分
>> x=linspace(-4,5); % 产生100个x的离散点
>> p=[1 -3 -11 27 10 -24]; %被积函数各项的系数
%函数计算的次数n
q= -0.4605 n= 53
再来看一个积分式
经典matlab数值积分与微分PPT培训课件
辛普森法
将积分区间划分为若干个等宽 的小区间，每个小区间上取抛 物线面积近似代替函数面积。
牛顿-莱布尼茨法
利用定积分的几何意义，通过 求和式来近似计算定积分的值 。
Matlab中的数值积分函数
quad
trapz
simpson
quadl
使用矩形法进行数值积 分。
使用梯形法进行数值积 分。
使用辛普森法进行数值 积分。
案例二
计算函数$f(x) = sin(x)$在区间 [0, π]的积分。
描述
使用Matlab的数值积分函数计 算函数$f(x) = sin(x)$在区间[0, π]的积分，并分析结果的准确性。
数值微分案例分析
案例一
计算函数$f(x) = x^3$在点x=2的导 数值。
描述
使用Matlab的数值微分函数（如diff 或gradient）计算函数$f(x) = x^3$ 在点x=2的导数值。
书籍推荐
《Matlab从入门到精通》、 《Matlab数值分析》等，适 合有一定基础的读者深入学 习。
Matlab的未来发展展望
云端化 随着云计算技术的发展，Matlab 可能会推出云端版本，让用户无 需安装软件即可使用Matlab的功 能。
可视化增强 Matlab在数据可视化方面具有优 势，未来可能会进一步增强其可 视化功能，提供更多样化的图表 和可视化效果。
截断误差
由于差分代替导数，会产生截断误差，这种误差的大小取决于差 分的阶数和步长。
舍入误差
由于计算机的浮点运算精度限制，会导致舍入误差，这种误差的 大小取决于计算机的浮点精度。
04
经典案例分析
数值积分案例分析
案例一
计算函数$f(x) = x^2$在区间 [0, 2]的积分。
Matlab数值积分和微分(一重、多重都有),简单实用
第8章MA TLAB数值积分与微分8.1 数值积分8.2 数值微分8.1 数值积分8.1.1 数值积分基本原理求解定积分的数值方法多种多样，如简单的梯形法、辛普生(Simpson)•法、牛顿－柯特斯(Newton-Cotes)法等都是经常采用的方法。
它们的基本思想都是将整个积分区间[a,b]分成n 个子区间[xi,xi+1]，i=1,2,…,n，其中x1=a，xn+1=b。
这样求定积分问题就分解为求和问题。
8.1.2 数值积分的实现方法1．变步长辛普生法基于变步长辛普生法，MA TLAB给出了quad函数来求定积分。
该函数的调用格式为：[I,n]=quad('fname',a,b,tol,trace)其中fname是被积函数名。
a和b分别是定积分的下限和上限。
tol用来控制积分精度，缺省时取tol=0.001。
trace控制是否展现积分过程，若取非0则展现积分过程，取0则不展现，缺省时取trace=0。
返回参数I即定积分值，n为被积函数的调用次数。
例8-1 求定积分。
(1) 建立被积函数文件fesin.m。
functionf=fesin(x)f=exp(-0.5*x).*sin(x+pi/6);(2) 调用数值积分函数quad求定积分。
[S,n]=quad('fesin',0,3*pi)S = 0.9008n = 772．牛顿－柯特斯法基于牛顿－柯特斯法，MA TLAB给出了quad8函数来求定积分。
该函数的调用格式为：[I,n]=quad8('fname',a,b,tol,trace)其中参数的含义和quad函数相似，只是tol的缺省值取10-6。
•该函数可以更精确地求出定积分的值，且一般情况下函数调用的步数明显小于quad函数，从而保证能以更高的效率求出所需的定积分值。
例8-2 求定积分。
(1) 被积函数文件fx.m。
function f=fx(x)f=x.*sin(x)./(1+cos(x).*cos(x));(2) 调用函数quad8求定积分。
Matlab中的数值积分和微分方法
Matlab中的数值积分和微分方法在数学和工程领域，数值积分和微分是解决问题的常见方法之一。
而在计算机科学中， Matlab作为一种强大的数值计算软件，提供了许多数值积分和微分的函数，使得这两个问题的解决变得更加简单和高效。
本文将探讨 Matlab 中常用的数值积分和微分方法，包括不定积分、定积分、数值微分和高阶数值微分。
我们将逐一讨论这些方法的原理和使用方法，并展示一些实际的应用案例，以帮助读者更好地理解和应用这些技术。
一、不定积分不定积分是指求一个函数的原函数。
在 Matlab 中，我们可以使用 `int` 函数来实现不定积分的计算。
例如，如果我们想求解函数 f(x) = x^2 的不定积分，可以使用下面的代码：```syms x;F = int(x^2);```这里的 `syms x` 表示将 x 定义为一个符号变量，`int(x^2)` 表示求解函数 x^2 的不定积分。
得到的结果 F 将是一个以 x 为变量的符号表达式。
除了求解简单函数的不定积分外，Matlab 还支持求解复杂函数的不定积分，例如三角函数、指数函数等。
我们只需要将函数表达式作为 `int` 函数的参数即可。
二、定积分定积分是指求函数在一个闭区间上的积分值。
在 Matlab 中，我们可以使用`integral` 函数来计算定积分。
例如，如果我们想计算函数 y = x^2 在区间 [0, 1] 上的积分值，可以使用下面的代码：```y = @(x) x^2;result = integral(y, 0, 1);```这里的 `@(x)` 表示定义一个匿名函数，`integral(y, 0, 1)` 表示求解函数 y = x^2 在区间 [0, 1] 上的积分。
得到的结果 result 将是一个数值。
与不定积分类似，Matlab 还支持对复杂函数求解定积分，只需要将函数表达式作为 `integral` 函数的第一个参数，并指定积分的区间。
第7章数值积分与数值微分
当步长h充分小时，有
f ( x) f ( x) h f ( x ) f ( x) h f ( x) f ( x) h
分别称为函数在x点处以h(h>0)为步长的向前差商、 向后差商和中心差商。当步长h(h>0)充分小时，函数f 在点x的微分接近于函数在该点的任意一种差分，而f 在点x的导数接近于函数在该点的任意一种差商。
7.2.2 数值微分的实现
计算任意函数f(x)在给定点x处数 值导数的方式有两种，一种是用多项 式或样条函数g(x)对f(x)进行逼近（插 值或拟合），然后用逼近函数g(x)在点 x处的导数作为f(x)在点x处的导数；另 一种是用f(x)在点x处的某种差商作为 其导数。
在MATLAB中，没有直接提供求数值导 数的函数，只有计算向前差分的函数diff，其 调用格式为： DX=diff(X)：计算向量X的向前差分 DX(i)=X(i+1)-X(i)，i=1,2,…,n-1。 DX=diff(X,n)：计算X的n阶向前差分 例如，diff(X,2)=diff(diff(X))。 DX=diff(A,n,dim)：计算矩阵A的n阶差分 dim=1时(缺省状态)，按列计算差分；dim=2， 按行计算差分。
[S,n]=quad(‘fesin’,0,3*pi) 或者@fesin S= 0.9008 n= 77
也可不建立关于被积函数的函数文件，而 使用语句函数（内联函数）求解， 命令如下： g=inline('exp(-0.5*x).*sin(x+pi/6)'); %定义语 句函数 [S,n]=quad(g,0,3*pi) %注意函数 名不加‘号或@符号 S= 0.9008 n= 77
2.5
1
第8章 MATLAB数值积分与微分PPT教学课件
2020/12/10
2
8.1.2 数值积分的实现方法 1．变步长辛普生法 基于变步长辛普生法，MATLAB给出了quad函数来
求定积分。该函数的调用格式为：
[I,n]=quad('fname',a,b,tol,trace) 其中fname是被积函数名。a和b分别是定积分的下
限和上限。tol用来控制积分精度，缺省时取 tol=0.001。trace控制是否展现积分过程，若取非0 则展现积分过程，取0则不展现，缺省时取 trace=0。返回参数I即定积分值，n为被积函数的 调用次数。
算向前差分的函数diff，其调用格式为： DX=diff(X)：计算向量X的向前差分，DX(i)=X(i+1)-X(i)，
i=1,2,…,n-1。 DX=diff(X,n)：计算X的n阶向前差分。例如，
diff(X,2)=diff(diff(X))。 DX=diff(A,n,dim)：计算矩阵A的n阶差分，dim=1时(缺省状
第8章 MATLAB数值积分与微分 8.1 数值积分 8.2 数值微分
2020/12/10
1
8.1 数值积分 8.1.1 数值积分基本原理
求解定积分的数值方法多种多样，如简单 的梯形法、辛普生(Simpson)•法、牛顿－柯 特斯(Newton-Cotes)法等都是经常采用的方 法。它们的基本思想都是将整个积分区间 [中 为a,x求b1]=分和a，成问xn题n个+。1=子b区。间这[样xi,求xi+定1]，积i分=1问,2,题…就,n分，解其
调用函数quad求定积分： format long; fx=inline('exp(-x)'); [I,n]=quad(fx,1,2.5,1e-10) I=
MATLAB教程第8章MATLAB数值积分与微分
MATLAB教程第8章MATLAB数值积分与微分1.数值积分数值积分是计算函数的定积分值的近似方法。
在MATLAB中，有几个函数可以帮助我们进行数值积分。
(1) quad函数quad函数是MATLAB中用于计算一维定积分的常用函数。
它的语法如下：I = quad(fun, a, b)其中，fun是被积函数的句柄，a和b分别是积分区间的下界和上界，I是近似的积分值。
例如，我们可以计算函数y=x^2在区间[0,1]内的积分值：a=0;b=1;I = quad(fun, a, b);disp(I);(2) integral函数integral函数是在MATLAB R2024a版本引入的新函数，它提供了比quad函数更稳定和准确的积分计算。
integral函数的语法如下：I = integral(fun, a, b)其中fun、a和b的含义与quad函数相同。
例如，我们可以使用integral函数计算函数y = x^2在区间[0, 1]内的积分值：a=0;b=1;I = integral(fun, a, b);disp(I);2.数值微分数值微分是计算函数导数的近似方法。
在MATLAB中，可以使用diff 函数计算函数的导数。
(1) diff函数diff函数用于计算函数的导数。
它的语法如下：derivative = diff(fun, x)其中，fun是需要计算导数的函数，x是自变量。
例如，我们可以计算函数y=x^2的导数：syms x;fun = x^2;derivative = diff(fun, x);disp(derivative);(2) gradient函数gradient函数可以计算多变量函数的梯度。
它的语法如下：[g1, g2, ..., gn] = gradient(fun, x1, x2, ..., xn)其中fun是需要计算梯度的函数，x1, x2, ..., xn是自变量。
例如，我们可以计算函数f=x^2+y^2的梯度：syms x y;fun = x^2 + y^2;[gx, gy] = gradient(fun, x, y);disp(gx);disp(gy);以上是MATLAB中进行数值积分和微分的基本方法和函数。
第8章 MATLAB数值积分与微分PPT课件
8.1.2 数值积分的实现方法 1．变步长辛普生法 基于变步长辛普生法，MATLAB给出了quad函数
来求定积分。该函数的调用格式为：
[I,n]=quad('fname',a,b,tol,trace) 其中fname是被积函数名。a和b分别是定积分的下
限和上限。tol用来控制积分精度，缺省时取 tol=0.001。trace控制是否展现积分过程，若取非0 则展现积分过程，取0则不展现，缺省时取 trace=0。返回参数I即定积分值，n为被积函数的 调用次数。
调用函数quad求定积分： format long; fx=inline('exp(-x)'); [I,n]=quad(fx,1,2.5,1e-10) I=
0.28579444254766 n=
65
调用函数quad8求定积分： format long; fx=inline('exp(-x)'); [I,n]=quad8(fx,1,2.5,1e-10) I=
例8-2 求定积分。 (1) 被积函数文件fx.m。 function f=fx(x) f=x.*sin(x)./(1+cos(x).*cos(x)); (2) 调用函数quad8求定积分。 I=quad8('fx',0,pi) I=
2.4674
例8-3 分别用quad函数和quad8函数求定积分的近 似值，并在相同的积分精度下，比较函数的调用 次数。
i=1,2,…,n-1。 DX=diff(X,n)：计算X的n阶向前差分。例如，
diff(X,2)=diff(diff(X))。 DX=diff(A,n,dim)：计算矩阵A的n阶差分，dim=1时(缺省状
态)，按列计算差分；dim=2，按行计算差分。
Ma ab 数值积分与数值微分
Matlab数值积分与数值微分Matlab数值积分1.一重数值积分的实现方法变步长辛普森法、高斯-克朗罗德法、梯形积分法1.1变步长辛普森法Matlab提供了quad函数和quadl函数用于实现变步长辛普森法求数值积分.调用格式为:[I,n]=Quad(@fname,a,b,tol,trace)[I,n]=Quadl(@fname,a,b,tol,trace)Fname是函数文件名，a,b分别为积分下限、积分上限；tol 为精度控制，默认为1.0×10-6，trace控制是否展开积分过程，若为0则不展开，非0则展开，默认不展开.返回值I为积分数值；n为调用函数的次数.--------------------------------------------------------------------- 例如：求∫e0.5x sin(x+π)dx3π的值.先建立函数文件fesin.mfunction f=fesin(x)f=exp(-0.5*x).*sin(x+(pi/6));再调用quad函数[I,n]=quad(@fesin,0,3*pi,1e-10)I=0.9008n=365--------------------------------------------------------------------- 例如：分别用quad函数和quadl函数求积分∫e0.5x sin(x+π6)dx3π的近似值，比较函数调用的次数.先建立函数文件fesin.mfunction f=fesin(x)f=exp(-0.5*x).*sin(x+(pi/6));formatlong[I,n]=quadl(@fesin,0,3*pi,1e-10)I=n=198[I,n]=quad(@fesin,0,3*pi,1e-10)I=n=365--------------------------------------------------------------------- 可以发现quadl函数调用原函数的次数比quad少，并且比quad函数求得的数值解更精确.1.2高斯-克朗罗德法Matlab提供了自适应高斯-克朗罗德法的quadgk函数来求震荡函数的定积分，函数的调用格式为：[I,err]=quadgk(@fname,a,b)Err返回近似误差范围，其他参数的意义与quad函数相同，积分上下限可以是-Inf或Inf，也可以是复数，若为复数则在复平面上求积分.--------------------------------------------------------------------- 例如：求积分∫xsinx1+cos2xdxπ的数值.先编写被积函数的m文件fsx.mfunction f=fsx(x)f=x.*sin(x)./(1+cos(x).^2);再调用quadgk函数I=quadgk(@fsx,0,pi)I=2.4674--------------------------------------------------------------------- 例如：求积分∫xsinx1+cos2xdx+∞−∞的值.先编写被积函数的m文件fsx.mfunction f=fsx(x)f=x.*sin(x)./(1+cos(x).^2); 再调用quadgk函数I=quadgk(@fsx,-Inf,Inf)I=-9.0671e+017---------------------------------------------------------------------1.3梯形积分法对于一些不知道函数关系的函数问题，只有实验测得的一组组样本点和样本值，由表格定义的函数关系求定积分问题用梯形积分法，其函数是trapz函数，调用格式为：I=Traps(X,Y)X,Y为等长的两组向量，对应着函数关系Y=f(X)X=(x1,x2,…,x n)(x1<x2<…<x n),Y=(y1,y2,…,y n)，积分区间是[x1,x n]--------------------------------------------------------------------- 例如：已知某次物理实验测得如下表所示的两组样本点.现已知变量x和变量y满足一定的函数关系，但此关系未知，设y=f(x)，求积分13.39∫f(x)dx1.38的数值.X=[1.38,1.56,2.21,3.97,5.51,7.79,9.19,11.12,13.39];Y=[3.35,3.96,5.12,8.98,11.46,17.63,24.41,29.83,32.21];I=trapz(X,Y)I=217.1033---------------------------------------------------------------------例如：用梯形积分法求积分：∫e −x dx 2.51的数值.x=1:0.01:2.5; y=exp(-x); I=trapz(x,y) I= 0.2858---------------------------------------------------------------------2. 多重数值积分的实现重积分的积分函数一般是二元函数f(x,y)或三元函数f(x,y,z)；形如：∫∫f (x,y )dxdy ba dc∫∫∫f(x,y,z)dxdydz b a dc f eMatlab 中有dblquad 函数和triplequad 函数来对上述两个积分实现.调用格式为：I=dblquad(@fun,a,b,c,d,tol)I=triplequad(@fun,a,b,c,d,e,f,tol)Fun 为被积函数，[a,b]为x 的积分区间；[c,d]为y 的积分区间；[e,f]为z 的积分区间.Dblquad 函数和triplequad 函数不允许返回调用的次数，如果需要知道函数调用的次数，则在定义被积函数的m 文件中增加一个计数变量，统计出被积函数被调用的次数.---------------------------------------------------------------------例如:计算二重积分I =∫∫√x 2+y 2dxdy π2−π2π2−π2的值.先编写函数文件fxy.mfunction f=fxy(x,y) global k; k=k+1;f=sqrt(x.^2+y.^2);再调用函数dblquadglobalk; k=0;I=dblquad(@fxy,-pi/2,pi/2,-pi/2,pi/2,1.0e-10) I= 11.8629 k k= 37656---------------------------------------------------------------------例如：求三重积分∫∫∫4xze −z2y−x 2dxdydz ππ1的值.编写函数文件fxyz1.mfunction f=fxyz1(x,y,z) global j; j=j+1;f=4*x.*z.*exp(-z.*z.*y-x.*x); 调用triplequad 函数 edit globalj;j=0;I=triplequad(@fxyz1,0,pi,0,pi,0,1,1.0e-10)I=1.7328jj=1340978---------------------------------------------------------------------Matlab数值微分1.数值微分与差商导数的三种极限定义f′(x)=limn→0f(x+h)−f(x)f′(x)=limn→0f(x)−f(x−h)hf′(x)=limn→0f(x+h2)−f(x−h2)h上述公式中假设h>0，引进记号：?f(x)=f(x+h)−f(x)?f(x)= f(x)−f(x−h)δf(x)= f(x+h2)−f(x−h2)称上述?f(x)、?f(x)、δf(x)为函数在x点处以h(h>0)为步长的向前差分、向后差分、中心差分，当步长h足够小时，有：f′(x)≈?f(x) hf′(x)≈?f(x)f′(x)≈δf(x) h?f(x) h 、?f(x)h、δf(x)h也分别被称为函数在x点处以h(h>0)为步长的向前差商、向后差商、中心差商.当h足够小时，函数f(x)在x点处的导数接近于在该点的任意一种差商，微分接近于在该点的任意一种差分.2.函数导数的求法2.1用多项式或样条函数g(x)对函数f(x)进行逼近(插值或拟合)，然后用逼近函数g(x)在点x处的导数作为f(x)在该点处的导数.2.2用f(x)在点x处的差商作为其导数.3.数值微分的实现方法Matlab中，只有计算向前差分的函数diff，其调用格式为：·DX=diff(X):计算向量X的向前差分，DX(i)=X(i+1)-X(i),i=1,2,…,n-1·DX=diff(X,n):计算向量X的n阶向前差分，例如diff(X,2)=diff(diff(X))·DX=diff(A,n,dim):计算矩阵A的n阶向前差分，dim=1(默认值)按列计算差分，dim=2按行计算差分.--------------------------------------------------------------------- 例如：生成6阶范德蒙德矩阵，然后分别按行、按列计算二阶向前差分A=vander(1:6)A=111111321684212438127931102425664164131256251252551777612962163661D2A1=diff(A,2,1)D2A1=180501220057011018200132019424200255030230200D2A2=diff(A,2,2)D2A2=000084211083612457614436920004008016540090015025--------------------------------------------------------------------- 例如：设f(x)=√x3+2x2−x+12+√(x+5)56+5x+2求函数f(x)的数值导数，并在同一坐标系中作出f’(x)的图像.已知函数f(x)的导函数如下：f′(x)=3x2+4x−12√x3+2x2−x+12+1√()56+5编辑函数文件fun7.m和fun8.mfunctionf=fun7(x)f=sqrt(x.^3+2*x.^2-x+12)+(x+5).^(1/6)+5*x+2;functionf=fun8(x)f=(3*x.^2+4*x-1)/2./sqrt(x.^3+2*x.^2-x+12)+1/6./(x+5).^(5/6)+5 ;x=-3:0.01:3;p=polyfit(x,fun7(x),5);用5次多项式拟合曲线dp=polyder(p);对拟合多项式进行求导dpx=polyval(dp,x);对dp在假设点的求函数值dx=diff(fun7([x,3.01]))/0.01;直接对dx求数值导数gx=fun8(x);求函数f的函数在假设点的导数plot(x,dpx,x,dx,'.',x,gx,'-')可以发现，最后得到的三条曲线基本重合.--------------------------------------------------------------------- 练习：A.用高斯-克朗罗德法求积分∫dx1+x2 +∞−∞的值并讨论计算方法的精确度.（该积分值为π）function f=fun9(x)f=1./(1+x.^2);formatlong[I,err]=quadgk(@fun9,-Inf,Inf)I=err=B.设函数f(x)=sin xx+cos2x用不同的办法求该函数的数值导数，并在同一坐标系中作出f′(x)的图像.已知f′(x)=x cos x+cos x cos2x−sin x+2sin x sin2x(x+cos2x)2function f=fun10(x)f=sin(x)./(x+cos(2*x));function f=fun11(x)f=(x.*cos(x)+cos(x).*cos(2*x)-sin(x)-2*sin(x).*sin(2*x))/(x+cos(2* x)).^2;x=-3:0.01:3;p=polyfit(x,fun10(x),5);dp=polyder(p);dpx=polyval(dp,x);dx=diff(fun10([x,3.01]))/0.01;gx=fun11(x);plot(x,dpx,'r:',x,dx,'.g',x,gx,'-k')。
MATLAB编程与工程应用——第8章 MATLAB数值积分与微分
exp81mmatlab数值积分与微分二数值积分的实现方法只是tol的缺省值取106该函数可以更精确地求出定积分的值且一般情况下函数调用的步数明显小于quad函数从而保证能以更高的效率求出所需的定积分值
MATLAB数值积分与微分 第8章 MATLAB数值积分与微分
8.1 数值积分 8.2 数值微分
第8章
第8章 MATLAB数值积分与微分
二、数值积分的实现方法
2．牛顿－柯特斯法 牛顿－ [I,n]=quadl('fname',a,b,tol,trace) 其中参数的含义和quad函数相似， 只是tol的缺省值取10-6。 该函数可以更精确地求出定积分的值，且一般情况下函 数调用的步数明显小于quad函数，从而保证能以更高的 效率求出所需的定积分值。 例8.2 求funx在[0,1]上的定积分。 exp8_2.m funx在[0,1]上的定积分。 上的定积分
第8章
MATLAB数值积分与微分
8.1 数值积分
二、数值积分的实现方法
1．变步长辛普生法 [I,n]=quad('fname',a,b,tol,trace) fname是被积函数名。 a和b分别是定积分的下限和上限。 tol用来控制积分精度，缺省时取tol=0.001。 trace控制是否展现积分过程，若取非0则展现积分过程， 取0则不展现，缺省时取trace=0。 I为定积分值 n为被积函数的调用次数。 例8.1 求funx在[0,1]上的定积分。 exp8_1.m funx在[0,1]上的定积分。 上的定积分
第8章
exp8_3.m
MATLAB数值积分与微分
8.1 数值积分
三、二元函数重积分的数值计算
Matlab数值积分与数值微分修订稿
M a t l a b数值积分与数值微分集团档案编码：[YTTR-YTPT28-YTNTL98-UYTYNN08]M a t l a b数值积分与数值微分Matlab数值积分1.一重数值积分的实现方法变步长辛普森法、高斯-克朗罗德法、梯形积分法1.1变步长辛普森法Matlab提供了quad函数和quadl函数用于实现变步长辛普森法求数值积分.调用格式为:[I,n]=Quad(@fname,a,b,tol,trace)[I,n]=Quadl(@fname,a,b,tol,trace)Fname是函数文件名，a,b分别为积分下限、积分上限；tol为精度控制，默认为1.0×10-6，trace控制是否展开积分过程，若为0则不展开，非0则展开，默认不展开.返回值I为积分数值；n为调用函数的次数.---------------------------------------------------------------------例如：求∫e e.ee eee(e+ee)e e3π的值.先建立函数文件fesin.mfunction f=fesin(x)f=exp(-0.5*x).*sin(x+(pi/6));再调用quad函数[I,n]=quad(@fesin,0,3*pi,1e-10)I=0.9008n=365---------------------------------------------------------------------例如：分别用quad函数和quadl函数求积分∫e e.ee eee(e+ee)e e3π的近似值，比较函数调用的次数.先建立函数文件fesin.mfunction f=fesin(x)f=exp(-0.5*x).*sin(x+(pi/6));formatlong[I,n]=quadl(@fesin,0,3*pi,1e-10)I=n=198[I,n]=quad(@fesin,0,3*pi,1e-10)I=n=365---------------------------------------------------------------------可以发现quadl函数调用原函数的次数比quad少，并且比quad函数求得的数值解更精确.1.2高斯-克朗罗德法Matlab提供了自适应高斯-克朗罗德法的quadgk函数来求震荡函数的定积分，函数的调用格式为：[I,err]=quadgk(@fname,a,b)Err返回近似误差范围，其他参数的意义与quad函数相同，积分上下限可以是-Inf或Inf，也可以是复数，若为复数则在复平面上求积分.---------------------------------------------------------------------例如：求积分∫e eee ee+eee e eee π的数值.先编写被积函数的m文件fsx.mfunction f=fsx(x)f=x.*sin(x)./(1+cos(x).^2);再调用quadgk函数I=quadgk(@fsx,0,pi)I=2.4674---------------------------------------------------------------------例如：求积分∫e eee ee+eee e eee +∞−∞的值.先编写被积函数的m文件fsx.mfunction f=fsx(x)f=x.*sin(x)./(1+cos(x).^2);再调用quadgk函数I=quadgk(@fsx,-Inf,Inf)I=-9.0671e+017---------------------------------------------------------------------1.3梯形积分法对于一些不知道函数关系的函数问题，只有实验测得的一组组样本点和样本值，由表格定义的函数关系求定积分问题用梯形积分法，其函数是trapz函数，调用格式为：I=Traps(X,Y)X,Y为等长的两组向量，对应着函数关系Y=f(X)X=(x1,x2,…,x n)(x1<x2<…<x n),Y=(y1,y2,…,y n)，积分区间是[x1,x n]---------------------------------------------------------------------例如：已知某次物理实验测得如下表所示的两组样本点.现已知变量x和变量y满足一定的函数关系，但此关系未知，设y=f(x)，求积分13.39∫e(e)ee1.38的数值.X=[1.38,1.56,2.21,3.97,5.51,7.79,9.19,11.12,13.39];Y=[3.35,3.96,5.12,8.98,11.46,17.63,24.41,29.83,32.21];I=trapz(X,Y)I=217.1033---------------------------------------------------------------------例如：用梯形积分法求积分：∫e −e ee 2.51的数值.x=1:0.01:2.5; y=exp(-x); I=trapz(x,y) I=0.2858---------------------------------------------------------------------2. 多重数值积分的实现重积分的积分函数一般是二元函数f(x,y)或三元函数f(x,y,z)；形如：∫∫e (e ,e )eeee ee ee∫∫∫e (e ,e ,e )eeeeee e ee ee eMatlab 中有dblquad 函数和triplequad 函数来对上述两个积分实现.调用格式为： I=dblquad(@fun,a,b,c,d,tol) I=triplequad(@fun,a,b,c,d,e,f,tol)Fun 为被积函数，[a,b]为x 的积分区间；[c,d]为y 的积分区间；[e,f]为z 的积分区间.Dblquad 函数和triplequad 函数不允许返回调用的次数，如果需要知道函数调用的次数，则在定义被积函数的m 文件中增加一个计数变量，统计出被积函数被调用的次数.---------------------------------------------------------------------例如:计算二重积分I =∫∫√e e +e e eeee e e−eee e−ee 的值.先编写函数文件fxy.mfunction f=fxy(x,y) global k; k=k+1;f=sqrt(x.^2+y.^2);再调用函数dblquadglobalk; k=0;I=dblquad(@fxy,-pi/2,pi/2,-pi/2,pi/2,1.0e-10) I= 11.8629 k k= 37656---------------------------------------------------------------------例如：求三重积分∫∫∫eeee −ee e −e eeeeeee eeeeee的值.编写函数文件fxyz1.m function f=fxyz1(x,y,z) global j; j=j+1;f=4*x.*z.*exp(-z.*z.*y-x.*x); 调用triplequad 函数 edit globalj; j=0;I=triplequad(@fxyz1,0,pi,0,pi,0,1,1.0e-10) I= 1.7328 jj=1340978---------------------------------------------------------------------Matlab数值微分1.数值微分与差商导数的三种极限定义e′(e)=eeee→e e(e+e)−e(e)ee′(e)=eeee→e e(e)−e(e−e)ee′(e)=eeee→e e(e+ee)−e(e−e e)e上述公式中假设h>0，引进记号：e(e)=e(e+e)−e(e)e(e)=e(e)−e(e−e)ee(e)=e(e+ee)−e(e−ee)称上述e(e)、e(e)、ee(e)为函数在x点处以h(h>0)为步长的向前差分、向后差分、中心差分，当步长h足够小时，有：e′(e)≈e(e) ee′(e)≈e(e) ee′(e)≈ee(e)ee(e) e 、e(e)e、ee(e)e也分别被称为函数在x点处以h(h>0)为步长的向前差商、向后差商、中心差商.当h足够小时，函数f(x)在x点处的导数接近于在该点的任意一种差商，微分接近于在该点的任意一种差分.2.函数导数的求法2.1用多项式或样条函数g(x)对函数f(x)进行逼近(插值或拟合)，然后用逼近函数g(x)在点x处的导数作为f(x)在该点处的导数.2.2用f(x)在点x处的差商作为其导数.3.数值微分的实现方法Matlab中，只有计算向前差分的函数diff，其调用格式为：·DX=diff(X):计算向量X的向前差分，DX(i)=X(i+1)-X(i),i=1,2,…,n-1·DX=diff(X,n):计算向量X的n阶向前差分，例如diff(X,2)=diff(diff(X))·DX=diff(A,n,dim):计算矩阵A的n阶向前差分，dim=1(默认值)按列计算差分，dim=2按行计算差分.---------------------------------------------------------------------例如：生成6阶范德蒙德矩阵，然后分别按行、按列计算二阶向前差分A=vander(1:6)A=111111321684212438127931102425664164131256251252551777612962163661D2A1=diff(A,2,1)D2A1=180501220057011018200132019424200255030230200D2A2=diff(A,2,2)D2A2=000084211083612457614436920004008016540090015025---------------------------------------------------------------------例如：设e(e)=√e e+ee e−e+ee+√(e+e)e+ee+e求函数f(x)的数值导数，并在同一坐标系中作出f’(x)的图像.已知函数f(x)的导函数如下：e′(e)=ee e+ee−ee e e+ee e−e+ee+ee+e编辑函数文件fun7.m和fun8.mfunctionf=fun7(x)f=sqrt(x.^3+2*x.^2-x+12)+(x+5).^(1/6)+5*x+2;functionf=fun8(x)f=(3*x.^2+4*x-1)/2./sqrt(x.^3+2*x.^2-x+12)+1/6./(x+5).^(5/6)+5; x=-3:0.01:3;p=polyfit(x,fun7(x),5);用5次多项式拟合曲线dp=polyder(p);对拟合多项式进行求导dpx=polyval(dp,x);对dp在假设点的求函数值dx=diff(fun7([x,3.01]))/0.01;直接对dx 求数值导数gx=fun8(x);求函数f 的函数在假设点的导数plot(x,dpx,x,dx,'.',x,gx,'-')可以发现，最后得到的三条曲线基本重合.--------------------------------------------------------------------- 练习：A. 用高斯-克朗罗德法求积分∫ee e +e e +∞−∞ 的值并讨论计算方法的精确度.（该积分值为π） function f=fun9(x)f=1./(1+x.^2);formatlong[I,err]=quadgk(@fun9,-Inf,Inf)I=err=B.设函数 e (e )=eee e + 用不同的办法求该函数的数值导数，并在同一坐标系中作出e ′(e )的图像. 已知e ′(e )=e eee e +eee e eee ee −eee e +e eee e eee ee (e eee ee ) function f=fun10(x)f=sin(x)./(x+cos(2*x));function f=fun11(x)f=(x.*cos(x)+cos(x).*cos(2*x)-sin(x)-2*sin(x).*sin(2*x))/(x+cos(2*x)).^2; x=-3:0.01:3;p=polyfit(x,fun10(x),5);dp=polyder(p);dpx=polyval(dp,x);dx=diff(fun10([x,3.01]))/0.01;gx=fun11(x);plot(x,dpx,'r:',x,dx,'.g',x,gx,'-k')。
Maab数值积分与数值微分
Matlab数值积分与数值微分Matlab数值积分1.一重数值积分的实现方法变步长辛普森法、高斯-克朗罗德法、梯形积分法1.1变步长辛普森法Matlab提供了quad函数和quadl函数用于实现变步长辛普森法求数值积分.调用格式为:[I,n]=Quad(@fname,a,b,tol,trace)[I,n]=Quadl(@fname,a,b,tol,trace)Fname是函数文件名，a,b分别为积分下限、积分上限；tol 为精度控制，默认为1.0×10-6，trace控制是否展开积分过程，若为0则不展开，非0则展开，默认不展开.返回值I为积分数值；n为调用函数的次数.--------------------------------------------------------------------- 例如：求∫e0.5x sin(x+π)dx3π的值.先建立函数文件fesin.mfunction f=fesin(x)f=exp(-0.5*x).*sin(x+(pi/6));再调用quad函数[I,n]=quad(@fesin,0,3*pi,1e-10)I=0.9008n=365--------------------------------------------------------------------- 例如：分别用quad函数和quadl函数求积分∫e0.5x sin(x+π6)dx3π的近似值，比较函数调用的次数.先建立函数文件fesin.mfunction f=fesin(x)f=exp(-0.5*x).*sin(x+(pi/6));formatlong[I,n]=quadl(@fesin,0,3*pi,1e-10)I=n=198[I,n]=quad(@fesin,0,3*pi,1e-10)I=n=365--------------------------------------------------------------------- 可以发现quadl函数调用原函数的次数比quad少，并且比quad函数求得的数值解更精确.1.2高斯-克朗罗德法Matlab提供了自适应高斯-克朗罗德法的quadgk函数来求震荡函数的定积分，函数的调用格式为：[I,err]=quadgk(@fname,a,b)Err返回近似误差范围，其他参数的意义与quad函数相同，积分上下限可以是-Inf或Inf，也可以是复数，若为复数则在复平面上求积分.--------------------------------------------------------------------- 例如：求积分∫xsinx1+cos2xdxπ的数值.先编写被积函数的m文件fsx.mfunction f=fsx(x)f=x.*sin(x)./(1+cos(x).^2);再调用quadgk函数I=quadgk(@fsx,0,pi)I=2.4674--------------------------------------------------------------------- 例如：求积分∫xsinx1+cos2xdx+∞−∞的值.先编写被积函数的m文件fsx.mfunction f=fsx(x)f=x.*sin(x)./(1+cos(x).^2); 再调用quadgk函数I=quadgk(@fsx,-Inf,Inf)I=-9.0671e+017---------------------------------------------------------------------1.3梯形积分法对于一些不知道函数关系的函数问题，只有实验测得的一组组样本点和样本值，由表格定义的函数关系求定积分问题用梯形积分法，其函数是trapz函数，调用格式为：I=Traps(X,Y)X,Y为等长的两组向量，对应着函数关系Y=f(X)X=(x1,x2,…,x n)(x1<x2<…<x n),Y=(y1,y2,…,y n)，积分区间是[x1,x n]--------------------------------------------------------------------- 例如：已知某次物理实验测得如下表所示的两组样本点.现已知变量x和变量y满足一定的函数关系，但此关系未知，设y=f(x)，求积分13.39∫f(x)dx1.38的数值.X=[1.38,1.56,2.21,3.97,5.51,7.79,9.19,11.12,13.39];Y=[3.35,3.96,5.12,8.98,11.46,17.63,24.41,29.83,32.21];I=trapz(X,Y)I=217.1033---------------------------------------------------------------------例如：用梯形积分法求积分：∫e −x dx 2.51的数值.x=1:0.01:2.5; y=exp(-x); I=trapz(x,y) I= 0.2858---------------------------------------------------------------------2. 多重数值积分的实现重积分的积分函数一般是二元函数f(x,y)或三元函数f(x,y,z)；形如：∫∫f (x,y )dxdy ba dc∫∫∫f(x,y,z)dxdydz b a dc f eMatlab 中有dblquad 函数和triplequad 函数来对上述两个积分实现.调用格式为：I=dblquad(@fun,a,b,c,d,tol)I=triplequad(@fun,a,b,c,d,e,f,tol)Fun 为被积函数，[a,b]为x 的积分区间；[c,d]为y 的积分区间；[e,f]为z 的积分区间.Dblquad 函数和triplequad 函数不允许返回调用的次数，如果需要知道函数调用的次数，则在定义被积函数的m 文件中增加一个计数变量，统计出被积函数被调用的次数.---------------------------------------------------------------------例如:计算二重积分I =∫∫√x 2+y 2dxdy π2−π2π2−π2的值.先编写函数文件fxy.mfunction f=fxy(x,y) global k; k=k+1;f=sqrt(x.^2+y.^2);再调用函数dblquadglobalk; k=0;I=dblquad(@fxy,-pi/2,pi/2,-pi/2,pi/2,1.0e-10) I= 11.8629 k k= 37656---------------------------------------------------------------------例如：求三重积分∫∫∫4xze −z2y−x 2dxdydz ππ1的值.编写函数文件fxyz1.mfunction f=fxyz1(x,y,z) global j; j=j+1;f=4*x.*z.*exp(-z.*z.*y-x.*x); 调用triplequad 函数 edit globalj;j=0;I=triplequad(@fxyz1,0,pi,0,pi,0,1,1.0e-10)I=1.7328jj=1340978---------------------------------------------------------------------Matlab数值微分1.数值微分与差商导数的三种极限定义f′(x)=limn→0f(x+h)−f(x)f′(x)=limn→0f(x)−f(x−h)hf′(x)=limn→0f(x+h2)−f(x−h2)h上述公式中假设h>0，引进记号：?f(x)=f(x+h)−f(x)?f(x)= f(x)−f(x−h)δf(x)= f(x+h2)−f(x−h2)称上述?f(x)、?f(x)、δf(x)为函数在x点处以h(h>0)为步长的向前差分、向后差分、中心差分，当步长h足够小时，有：f′(x)≈?f(x) hf′(x)≈?f(x)f′(x)≈δf(x) h?f(x) h 、?f(x)h、δf(x)h也分别被称为函数在x点处以h(h>0)为步长的向前差商、向后差商、中心差商.当h足够小时，函数f(x)在x点处的导数接近于在该点的任意一种差商，微分接近于在该点的任意一种差分.2.函数导数的求法2.1用多项式或样条函数g(x)对函数f(x)进行逼近(插值或拟合)，然后用逼近函数g(x)在点x处的导数作为f(x)在该点处的导数.2.2用f(x)在点x处的差商作为其导数.3.数值微分的实现方法Matlab中，只有计算向前差分的函数diff，其调用格式为：·DX=diff(X):计算向量X的向前差分，DX(i)=X(i+1)-X(i),i=1,2,…,n-1·DX=diff(X,n):计算向量X的n阶向前差分，例如diff(X,2)=diff(diff(X))·DX=diff(A,n,dim):计算矩阵A的n阶向前差分，dim=1(默认值)按列计算差分，dim=2按行计算差分.--------------------------------------------------------------------- 例如：生成6阶范德蒙德矩阵，然后分别按行、按列计算二阶向前差分A=vander(1:6)A=111111321684212438127931102425664164131256251252551777612962163661D2A1=diff(A,2,1)D2A1=180501220057011018200132019424200255030230200D2A2=diff(A,2,2)D2A2=000084211083612457614436920004008016540090015025--------------------------------------------------------------------- 例如：设f(x)=√x3+2x2−x+12+√(x+5)56+5x+2求函数f(x)的数值导数，并在同一坐标系中作出f’(x)的图像.已知函数f(x)的导函数如下：f′(x)=3x2+4x−12√x3+2x2−x+12+1√()56+5编辑函数文件fun7.m和fun8.mfunctionf=fun7(x)f=sqrt(x.^3+2*x.^2-x+12)+(x+5).^(1/6)+5*x+2;functionf=fun8(x)f=(3*x.^2+4*x-1)/2./sqrt(x.^3+2*x.^2-x+12)+1/6./(x+5).^(5/6)+5 ;x=-3:0.01:3;p=polyfit(x,fun7(x),5);用5次多项式拟合曲线dp=polyder(p);对拟合多项式进行求导dpx=polyval(dp,x);对dp在假设点的求函数值dx=diff(fun7([x,3.01]))/0.01;直接对dx求数值导数gx=fun8(x);求函数f的函数在假设点的导数plot(x,dpx,x,dx,'.',x,gx,'-')可以发现，最后得到的三条曲线基本重合.--------------------------------------------------------------------- 练习：A.用高斯-克朗罗德法求积分∫dx1+x2 +∞−∞的值并讨论计算方法的精确度.（该积分值为π）function f=fun9(x)f=1./(1+x.^2);formatlong[I,err]=quadgk(@fun9,-Inf,Inf)I=err=B.设函数f(x)=sin xx+cos2x用不同的办法求该函数的数值导数，并在同一坐标系中作出f′(x)的图像.已知f′(x)=x cos x+cos x cos2x−sin x+2sin x sin2x(x+cos2x)2function f=fun10(x)f=sin(x)./(x+cos(2*x));function f=fun11(x)f=(x.*cos(x)+cos(x).*cos(2*x)-sin(x)-2*sin(x).*sin(2*x))/(x+cos(2* x)).^2;x=-3:0.01:3;p=polyfit(x,fun10(x),5);dp=polyder(p);dpx=polyval(dp,x);dx=diff(fun10([x,3.01]))/0.01;gx=fun11(x);plot(x,dpx,'r:',x,dx,'.g',x,gx,'-k')。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
Matlab数值积分与数值微分Matlab数值积分1.一重数值积分的实现方法变步长辛普森法、高斯-克朗罗德法、梯形积分法1.1变步长辛普森法Matlab提供了quad函数和quadl函数用于实现变步长辛普森法求数值积分.调用格式为:[I,n]=Quad(@fname,a,b,tol,trace)[I,n]=Quadl(@fname,a,b,tol,trace)Fname是函数文件名，a,b分别为积分下限、积分上限；tol 为精度控制，默认为1.0×10-6，trace控制是否展开积分过程，若为0则不展开，非0则展开，默认不展开.返回值I为积分数值；n为调用函数的次数.--------------------------------------------------------------------- 例如：求∫e0.5x sin(x+π)dx3π的值.先建立函数文件fesin.mfunction f=fesin(x)f=exp(-0.5*x).*sin(x+(pi/6));再调用quad函数[I,n]=quad(@fesin,0,3*pi,1e-10)I=0.9008n=365--------------------------------------------------------------------- 例如：分别用quad函数和quadl函数求积分∫e0.5x sin(x+π6)dx3π的近似值，比较函数调用的次数.先建立函数文件fesin.mfunction f=fesin(x)f=exp(-0.5*x).*sin(x+(pi/6));formatlong[I,n]=quadl(@fesin,0,3*pi,1e-10)I=n=198[I,n]=quad(@fesin,0,3*pi,1e-10)I=n=365--------------------------------------------------------------------- 可以发现quadl函数调用原函数的次数比quad少，并且比quad函数求得的数值解更精确.1.2高斯-克朗罗德法Matlab提供了自适应高斯-克朗罗德法的quadgk函数来求震荡函数的定积分，函数的调用格式为：[I,err]=quadgk(@fname,a,b)Err返回近似误差范围，其他参数的意义与quad函数相同，积分上下限可以是-Inf或Inf，也可以是复数，若为复数则在复平面上求积分.--------------------------------------------------------------------- 例如：求积分∫xsinx1+cos2xdxπ的数值.先编写被积函数的m文件fsx.mfunction f=fsx(x)f=x.*sin(x)./(1+cos(x).^2);再调用quadgk函数I=quadgk(@fsx,0,pi)I=2.4674--------------------------------------------------------------------- 例如：求积分∫xsinx1+cos2xdx+∞−∞的值.先编写被积函数的m文件fsx.mfunction f=fsx(x)f=x.*sin(x)./(1+cos(x).^2); 再调用quadgk函数I=quadgk(@fsx,-Inf,Inf)I=-9.0671e+017---------------------------------------------------------------------1.3梯形积分法对于一些不知道函数关系的函数问题，只有实验测得的一组组样本点和样本值，由表格定义的函数关系求定积分问题用梯形积分法，其函数是trapz函数，调用格式为：I=Traps(X,Y)X,Y为等长的两组向量，对应着函数关系Y=f(X)X=(x1,x2,…,x n)(x1<x2<…<x n),Y=(y1,y2,…,y n)，积分区间是[x1,x n]--------------------------------------------------------------------- 例如：已知某次物理实验测得如下表所示的两组样本点.现已知变量x和变量y满足一定的函数关系，但此关系未知，设y=f(x)，求积分13.39∫f(x)dx1.38的数值.X=[1.38,1.56,2.21,3.97,5.51,7.79,9.19,11.12,13.39];Y=[3.35,3.96,5.12,8.98,11.46,17.63,24.41,29.83,32.21];I=trapz(X,Y)I=217.1033---------------------------------------------------------------------例如：用梯形积分法求积分：∫e −x dx 2.51的数值.x=1:0.01:2.5; y=exp(-x); I=trapz(x,y) I= 0.2858---------------------------------------------------------------------2. 多重数值积分的实现重积分的积分函数一般是二元函数f(x,y)或三元函数f(x,y,z)；形如：∫∫f (x,y )dxdy ba dc∫∫∫f(x,y,z)dxdydz b a dc f eMatlab 中有dblquad 函数和triplequad 函数来对上述两个积分实现.调用格式为：I=dblquad(@fun,a,b,c,d,tol)I=triplequad(@fun,a,b,c,d,e,f,tol)Fun 为被积函数，[a,b]为x 的积分区间；[c,d]为y 的积分区间；[e,f]为z 的积分区间.Dblquad 函数和triplequad 函数不允许返回调用的次数，如果需要知道函数调用的次数，则在定义被积函数的m 文件中增加一个计数变量，统计出被积函数被调用的次数.---------------------------------------------------------------------例如:计算二重积分I =∫∫√x 2+y 2dxdy π2−π2π2−π2的值.先编写函数文件fxy.mfunction f=fxy(x,y) global k; k=k+1;f=sqrt(x.^2+y.^2);再调用函数dblquadglobalk; k=0;I=dblquad(@fxy,-pi/2,pi/2,-pi/2,pi/2,1.0e-10) I= 11.8629 k k= 37656---------------------------------------------------------------------例如：求三重积分∫∫∫4xze −z2y−x 2dxdydz ππ1的值.编写函数文件fxyz1.mfunction f=fxyz1(x,y,z) global j; j=j+1;f=4*x.*z.*exp(-z.*z.*y-x.*x); 调用triplequad 函数 edit globalj;j=0;I=triplequad(@fxyz1,0,pi,0,pi,0,1,1.0e-10)I=1.7328jj=1340978---------------------------------------------------------------------Matlab数值微分1.数值微分与差商导数的三种极限定义f′(x)=limn→0f(x+h)−f(x)f′(x)=limn→0f(x)−f(x−h)hf′(x)=limn→0f(x+h2)−f(x−h2)h上述公式中假设h>0，引进记号：?f(x)=f(x+h)−f(x)?f(x)= f(x)−f(x−h)δf(x)= f(x+h2)−f(x−h2)称上述?f(x)、?f(x)、δf(x)为函数在x点处以h(h>0)为步长的向前差分、向后差分、中心差分，当步长h足够小时，有：f′(x)≈?f(x) hf′(x)≈?f(x)f′(x)≈δf(x) h?f(x) h 、?f(x)h、δf(x)h也分别被称为函数在x点处以h(h>0)为步长的向前差商、向后差商、中心差商.当h足够小时，函数f(x)在x点处的导数接近于在该点的任意一种差商，微分接近于在该点的任意一种差分.2.函数导数的求法2.1用多项式或样条函数g(x)对函数f(x)进行逼近(插值或拟合)，然后用逼近函数g(x)在点x处的导数作为f(x)在该点处的导数.2.2用f(x)在点x处的差商作为其导数.3.数值微分的实现方法Matlab中，只有计算向前差分的函数diff，其调用格式为：·DX=diff(X):计算向量X的向前差分，DX(i)=X(i+1)-X(i),i=1,2,…,n-1·DX=diff(X,n):计算向量X的n阶向前差分，例如diff(X,2)=diff(diff(X))·DX=diff(A,n,dim):计算矩阵A的n阶向前差分，dim=1(默认值)按列计算差分，dim=2按行计算差分.--------------------------------------------------------------------- 例如：生成6阶范德蒙德矩阵，然后分别按行、按列计算二阶向前差分A=vander(1:6)A=111111321684212438127931102425664164131256251252551777612962163661D2A1=diff(A,2,1)D2A1=180501220057011018200132019424200255030230200D2A2=diff(A,2,2)D2A2=000084211083612457614436920004008016540090015025--------------------------------------------------------------------- 例如：设f(x)=√x3+2x2−x+12+√(x+5)56+5x+2求函数f(x)的数值导数，并在同一坐标系中作出f’(x)的图像.已知函数f(x)的导函数如下：f′(x)=3x2+4x−12√x3+2x2−x+12+1√()56+5编辑函数文件fun7.m和fun8.mfunctionf=fun7(x)f=sqrt(x.^3+2*x.^2-x+12)+(x+5).^(1/6)+5*x+2;functionf=fun8(x)f=(3*x.^2+4*x-1)/2./sqrt(x.^3+2*x.^2-x+12)+1/6./(x+5).^(5/6)+5 ;x=-3:0.01:3;p=polyfit(x,fun7(x),5);用5次多项式拟合曲线dp=polyder(p);对拟合多项式进行求导dpx=polyval(dp,x);对dp在假设点的求函数值dx=diff(fun7([x,3.01]))/0.01;直接对dx求数值导数gx=fun8(x);求函数f的函数在假设点的导数plot(x,dpx,x,dx,'.',x,gx,'-')可以发现，最后得到的三条曲线基本重合.--------------------------------------------------------------------- 练习：A.用高斯-克朗罗德法求积分∫dx1+x2 +∞−∞的值并讨论计算方法的精确度.（该积分值为π）function f=fun9(x)f=1./(1+x.^2);formatlong[I,err]=quadgk(@fun9,-Inf,Inf)I=err=B.设函数f(x)=sin xx+cos2x用不同的办法求该函数的数值导数，并在同一坐标系中作出f′(x)的图像.已知f′(x)=x cos x+cos x cos2x−sin x+2sin x sin2x(x+cos2x)2function f=fun10(x)f=sin(x)./(x+cos(2*x));function f=fun11(x)f=(x.*cos(x)+cos(x).*cos(2*x)-sin(x)-2*sin(x).*sin(2*x))/(x+cos(2* x)).^2;x=-3:0.01:3;p=polyfit(x,fun10(x),5);dp=polyder(p);dpx=polyval(dp,x);dx=diff(fun10([x,3.01]))/0.01;gx=fun11(x);plot(x,dpx,'r:',x,dx,'.g',x,gx,'-k')。