数值积分和微分方程数值解.pptVIP

5.3节数值积分和微分方程

数值解一．数值定积分求面积【例5-3-1】用数值积分法求由 ，y=0, x=0与x=10围成的图形面积，并讨论步长和积分方法对精度的影响。解:◆原理用矩形法和梯形法分别求数值积分并作比较，步长的变化用循环语句实现。MATLAB中的定积分有专门的函数QUAD,QUADL等实现。为了弄清原理，我们先用直接编程的方法来计算，然后再介绍定积分函数及其调用方法。设x向量的长度取为n，即将积分区间分为n-1段， 各段长度为 。算出各点的 ，则矩形法数值积分公式为：梯形法的公式为：比较两个公式，它们之间的差别只是 。在MATLAB中，把向量中各元素叠加的命令是sum。把向量中各元素按梯形法叠加的命令是trapz。梯形法的几何意义是把被积分的函数的各计算点以直线相联，形成许多窄长梯形条，然后叠加，我们把两种算法都编入同一个程序进行比较。123矩形和梯形定积分公式求面积的数值积分程序exn531fordx=[2,1,0.5,0.1] %设不同步长 x=0:.1:10;y=-x.*x+115;%取较密的函数样本plot(x,y),holdon %画出被积曲线并保持x1=0:dx:10;y1=-x1.*x1+115; %求取样点上的y1%用矩形（欧拉）法求积分，注意末尾去掉一个点n=length(x1);s=sum(y1(1:n-1))*dx; q=trapz(y1)*dx; %用梯形法求积分stairs(x1,y1), %画出欧拉法的积分区域plot(x1,y1) %画出梯形法的积分区域[dx,s,q],pause(1),holdoff，end程序运行的结果如下：程序exn531运行结果 程序运行的结果如下： 步长dx矩形法解s梯形法解q 2 910 810 1 865 815 .5 841.25 816.25 .1 821.65 816.65用解析法求出的精确解为2450/3=816.6666...。dx=2时矩形法和梯形法的积分面积见图5-4-1.。在曲线的切线斜率为负的情况下，矩形法的积分结果一定偏大，梯形法是由各采样点联线包围的面积，在曲线曲率为负（上凸）时，其积分结果一定偏小，因此精确解在这两者之间。由这结果也能看出，在步长相同时，梯形法的精度比矩形法高。矩形法数字积分的演示程序rsumsMATLAB中有一个矩形法数字积分的演示程序rsums，可以作一个对比。键入rsums(115-x.^2,0,10)就得到右图。图中表示了被积函数的曲线和被步长分割的小区间，并按各区间中点的函数值构成了各个窄矩形面积。用鼠标拖动图下方的滑尺可以改变步长的值，图的上方显示的是这些矩形面积叠加的结果。MATLAB内的数值定积分函数在实际工作中，用MATLAB中的定积分求面积的函数quad和quadl可以得到比自编程序更高的精度，因为quad函数用的是辛普生法，即把被积函数用二次曲线逼近的算法，而quadl函数采用了更高阶的逼近方法。它们的调用格式如下：Q=QUADL(FUN,A,B,TOL)其中，FUN是表示被积函数的字符串，A是积分下限，B是积分上限。TOL是规定计算的容差，其默认值为1e-6例如，键入 S=quad(-x.*x+115,0,10) 得到 二．求两条曲线所围图形的面积【例5-3-2】。设计算区间[0,4]上两曲线所围面积。解：◆原理：先画出图形，dx=input(dx=);x=0:dx:4;f=exp(-(x-2).^2.*cos(pi*x));g=4*cos(x-2);plot(x,f,x,g,:r)得到右图。从图上看到，其中既有f(x)g(x)的区域，也有f(x)g(x)的区域，求两条曲线所围图形的面积(1)若要求两曲线所围总面积（不管正负），则可加一条语句 s=trapz(abs(f-g))*dx,在dx=0.001时，得到s=6.47743996919702若要求两曲线所围的f(x)g(x)的正面积，则需要一定的技巧.◆方法一。先求出交点x1，再规定积分上下限。x1=fzero(exp(-(x-2).^2.*cos(pi*x))-4*cos(x-2),1)%把积分限设定为0~x1，求出积分结果再乘以2：x=0:dx:x1;f=exp(-(x-2).^2.*cos(pi*x));g=4*cos(x-2);s1