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

《数值积分和微分方程数值解》由会员上传分享，免费在线阅读，更多相关内容在教育资源-天天文库。

1、5.3节数值积分和微分方程数值解一．数值定积分求面积【例5-3-1】用数值积分法求由，y=0,x=0与x=10围成的图形面积，并讨论步长和积分方法对精度的影响。解:◆原理用矩形法和梯形法分别求数值积分并作比较，步长的变化用循环语句实现。MATLAB中的定积分有专门的函数QUAD,QUADL等实现。为了弄清原理，我们先用直接编程的方法来计算，然后再介绍定积分函数及其调用方法。设x向量的长度取为n，即将积分区间分为n-1段，各段长度为。算出各点的，则矩形法数值积分公式为：矩形和梯形定积分公式梯形法的公式为：比较两个公式，它们之间的差别只

2、是。在MATLAB中，把向量中各元素叠加的命令是sum。把向量中各元素按梯形法叠加的命令是trapz。梯形法的几何意义是把被积分的函数的各计算点以直线相联，形成许多窄长梯形条，然后叠加，我们把两种算法都编入同一个程序进行比较。求面积的数值积分程序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=

3、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梯形法解q29108101865815.5841.25816.25.1821.65816.65用解析法求出的精确解为2450/3=816.6666...。dx=2时矩形法和梯形法的积分面积见图5-4-1.。在曲线的切线

4、斜率为负的情况下，矩形法的积分结果一定偏大，梯形法是由各采样点联线包围的面积，在曲线曲率为负（上凸）时，其积分结果一定偏小，因此精确解在这两者之间。由这结果也能看出，在步长相同时，梯形法的精度比矩形法高。矩形法数字积分的演示程序rsumsMATLAB中有一个矩形法数字积分的演示程序rsums，可以作一个对比。键入rsums('115-x.^2',0,10)就得到右图。图中表示了被积函数的曲线和被步长分割的小区间，并按各区间中点的函数值构成了各个窄矩形面积。用鼠标拖动图下方的滑尺可以改变步长的值，图的上方显示的是这些矩形面积叠加的结果。

5、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)得到S=8.166666666666666e+002二．求两条曲线所围图形的面积

6、【例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若

7、要求两曲线所围的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=2*trapz(abs(f-g))*dx在设定dx=0.001时，得到s1=2.30330486000857求两条曲线所围图形的面积(2)方法二。调用MATL

8、AB中求面积函数quad。这里的关键是建立一个函数文件，把e1=f(x)-g(x)>0的部分取出来。利用逻辑算式(e1>0)，它在e1>0处取值为1，在e1<0处则为零。让逻辑函数(e1>0)与e1作元素群乘法，正的e1