Matlab教学课件：第三讲 3数值计算3.ppt

Matlab基础第四讲数值积分q 问题背景数值积分数值积分u 定积分计算的基本公式是牛顿－莱布尼兹公式但当被积函数的原函数不知道时，如何计算？这时就需要利用近似计算特别是在许多实际应用中，被积函数甚至没有解析表达式，而是一条实验记录曲线，或一组离散的采样值，此时只能用近似方法计算定积分u 本节主要研究定积分的三种近似计算算法：矩形法、梯形法和抛物线法同时介绍 Matlab 计算定积分的相关函数q 矩形法u 定积分的定义：实验二实验二、、定积分的近似计算定积分的近似计算矩形矩形法法n 充分大，△x 充分小u 定积分的近似：l 通常我们取左点法右点法中点法l 点 可以任意选取，常见的取法有： 左端点 ，右端点 和中点 步长节点u 右点法：u 中点法：u 左点法：左点左点法、右点法和中点法法、右点法和中点法解：矩形法举例矩形法举例==>==>==>==> h =1/100=0.01, xi = i*h, a=0, b=1, n=100 u 例：用不同的矩形法计算下面的定积分 ( 取 n=100 )， 并比较这三种方法的相对误差。

l 左点法：l 右点法：l 中点法：(i = 0,1,2,...,100)l 理论值：l 左点法相对误差：u 误差分析矩形法举例矩形法举例l 右点法相对误差：l 中点法相对误差：不同的方法有不同的计算精度有没有更好的近似计算定积分的方法有没有更好的近似计算定积分的方法 ?定积分几何意义定积分几何意义u 曲边小梯形的面积可以由直边小梯形的面积来近似u 整个曲边梯形的面积：梯形法梯形法u 如果我们 n 等分区间 [a,b]，，即令：则则==>==>==>==>梯形公式梯形公式梯形法梯形法梯形公式与中点公式有什么区别梯形公式与中点公式有什么区别 ?解：==>==>==>==>u 例：用梯形法计算下面定积分 ( 取 n=100 )， 并计算相对误差梯形法举例梯形法举例a=0, b=1, n=100, f (x) = 1/( 1+x2 ) ==>==>==>==> h =1/100=0.01, xi = i*h, yi = f (xi) l 相对误差：u 2n 等分区间 [a,b] ，得该直线用抛物线代替，计算精度是否会更好？u 计算每个节点上的函数值：抛物线法抛物线法u 在区间 [x0, x2] 上，用过以下三点的抛物线来近似原函数 f (x) 。

u 设过以上三点的抛物线方程为：则在则在区间 [x0, x2] 上，有上，有怎麽来的？怎麽来的？y =   x2 +   x +   = p1(x) 抛物线法抛物线法u 同理可得：u 相加即得：抛物线法抛物线法u 整理后可得：或辛普森 (Simpson) 公式抛物线法公式抛物线法抛物线法==>==>==>==>u 例：用抛物线法计算下面定积分 ( 取 n=100 )， 并计算相对误差解：a=0, b=1, n=100, yi = f (xi) = 1/( 1+xi2 ) l 相对误差：抛物线法抛物线法u 梯形法：trapztrapz(x,y)x 为分割点（节点）组成的向量，y 为被积函数在节点上的函数值组成的向量q Matlab 近似计算定积分的相关函数Matlab 计算定积分函数介绍计算定积分函数介绍前面的做法u 例：用梯形法计算下面定积分 ( 取 n=100) 解：a=0, b=1, n=100, yi = f (xi) = 1/( 1+xi2 ) >>>> x=0:1/100:1;>>>> y=1./(1+x.^2);>>>> trapz(x, y)trapz函数trapz(x,1./(1+x.^2))trapz 举例举例quad(f,a,b,tol)f = f(x) 为被积函数，[a,b] 为积分区间，tol 为计算精度将自变量看成是将自变量看成是向量向量u 抛物线法：quadl 不用自己分割积分区间l 可以指定计算精度，若不指定，缺省精度是 10-6l 精度越高，函数运行的时间越长l 此处的函数 f 是数值形式，应该使用数组运算，即 点运算：.*，./ ，.\ ，.^ 注：抛物线法抛物线法解：解：>>>> quad('1./(1+x.^2)',0,1)>>>> quad('1./(1+x.^2)',0,1,10e-10)>>>> quad('1./(1+x.^2)',0,1,10e-16)函数表达式一定要用函数表达式一定要用 单引号单引号 括起来！括起来！涉及的运算一定要用涉及的运算一定要用 数组运算数组运算！！u 例：用 quad 计算定积分：quad 举例举例q 抛物线法计算二重积分： dblquad abcd的意思的意思怎麽了？怎麽了？dblquad(f,a,b,c,d,tol)u tol 为计算精度，若不指定，则缺省精度为 10-6u f(x,y) 可以由 inline 定义，或通过一个函数句柄传递u [a,b] 是第一积分变量的积分区间，[c,d] 是第二积分变量 的积分区间按字母顺序，大写字母排在小写字母的前面二重积分的计算二重积分的计算>>>> f=inline('4*x*y+3*y^2');>>>> I=dblquad(f, -1,1,0,2)u f(x,y) 中关于第一自变量的运算是数组运算， 即把 x 看成是向量，y 看成是标量。

也可以全部采用数组运算例2：计算二重积分>>>> dblquad(inline('4*x*y+3*x^2'),-1,1,0,2)>>>> dblquad(inline('4*x*y+3*x.^2'),-1,1,0,2)X例1：计算二重积分dblquad 举例举例例：计算二重积分>>>> dblquad(@(x,y)4*x*y+3*x.^2 , -1,1, 0, 2)指定 x、、y 分别是第一第一和第二第二积分变量>>>> dblquad(inline('4*x*y+3*x.^2'), -1,1, 0, 2)q 被积函数 f (x,y) 的另一种定义方法：匿名函数>>>> dblquad(@(y,x)4*x*y+3*x.^2 , -1,1, 0, 2)下面的命令运行结果和上面的一样吗？dblquad 举例举例int(f,a,b) 计算 f 关于默认自变量 的定积分，积分区间为[a,b]int(f) 计算 f 关于默认自变量 的不定积分int(f,v,a,b) 计算函数 f 关于自变量 v 的定积分，积分区间为 [a, b]int(f,v) 计算函数 f 关于自变量 v 的不定积分findsym(f,1)q 符号积分： intint 符号积分符号积分>>>> syms x y;>>>> f=y*sin(x);>>>> int(f,x)>>>> int(f,y)>>>> int(f)>>>> int('a+b')ans=-y*cos(x)ans=1/2*y^2*sin(x)ans=-y*cos(x)ans=a*b+1/2*b^2u 例：指出下面各条语句的输出结果int 举例举例u 例：用 int 函数计算定积分：解：>>>> syms x;>>>> f=1/(1+x^2); >>>> int(f,x,0,1)>>>> f=sym('1/(1+x^2)'); >>>> int(f,x,0,1)>>>> int('1/(1+x^2)',x,0,1)或>>>> int('1/(1+x^2)',0,1)或或int 举例举例double(a) 将将 a 转化为双精度型，若转化为双精度型，若 a 是字符，则取对应的是字符，则取对应的 ASCII 码码>>>> a=3;>>>> double(a)>>>> double('a')例：ans = 3ans = 97其它相关函数其它相关函数>>>> x=1:0.001:2;>>>> y=exp(x.^(-2));>>>> trapz(x,y)l 梯形法：l 抛物线法：>>>> quad('exp(x.^(-2))',1,2,10e-10)l 符号积分法：>>>> syms x>>>> int('exp(x^(-2))',x,1,2)例 1：用 Matlab 函数近似计算积分数值实验数值实验l 抛物线法：>>>> dblquad(inline('x+y^2')),0,2,-1,1)l 符号积分法：>>>> f=int('x+y^2','y',-1,1);>>>> int(f,0,2)数值实验数值实验例 2：用 Matlab 函数近似计算二重积分q 作业（要求写实验报告）上机作业上机作业u 用不同的方法计算下面的定积分 , 并对所采用的方法进行比较。