【笔记】定积分的近似计算 | 您所在的位置:网站首页 › 近似数的步骤和方法 › 【笔记】定积分的近似计算 |
目标:加深理解积分理论中分割、近似、求和、取极限的思想方法,了解定积分近似计算的矩形法、梯形法与抛物线法,会用MATLAB语言编写求定积分近似值的程序,会用MALAB中的命令求定积分。 预备知识 在许多实际问题中,常常需要计算定积分的值。根据微分学基本定理,若被积函数f(x)在区间[a,b]上连续,只需要找到被积函数的一个原函数F(x),就可以用牛顿-莱布尼茨公式求出定积分值。但在工程技术与科学实验中,有些定积分的被积函数的原函数可能求不出来,即使可求出,计算也可能相当复杂。特别地,当被积函数是图形或表格时,更不能用牛顿-莱布尼茨公式计算。因此必须寻求定积分的近似计算方法。大多数实际问题的积分需要用数值积分方法求出近似结果。数值积分原则上可以用于计算各种被积函数的定积分,无论被积函数是解析形式还是数表形式,其基本原理都是用多项式函数近似代替被积函数,用对多项式的积分结果近似代替对被积函数的积分。由于所选多项式的不同,可以有许多种数值积分方法,下面介绍最常用的几种差值型数值积分方法。 1、矩形法 定积分的几何意义是计算曲边梯形的面积,如果将区间[a,b]n等分,每个小区间上都是一个小的曲边梯形,用一个个小矩形代替这些小曲边梯形,然后把所有小矩形的面积加起来就近似等于整个曲边梯形的面积,于是便求出了定积分的近似值,这就是矩形法的基本原理。 2、梯形法 将积分区间[a,b]n等分,用线段依次连接各分点,每段都形成一个小的直角梯形。如果用这些小直角梯形面积之和代替原来的小曲边梯形面积之和,就可求得定积分的近似值。 3、抛物线法 在此不做说明 4、相关的MATLAB命令 命令:sum(a),求数组a的和 例1 调用命令sum求矩阵x的各列元素的累加和。 >> x=[1 2 3;4 5 6;7 8 9]; >> sum(x) ans = 12 15 18 命令:format long,长格式,即屏幕显示15位有效数字。 命令:double(),若输入的是字符则转化为相应的ASCII码;若输入的是整数值则转化为相应的实型数值。 命令:quad(),抛物线法求数值积分。 格式:quad(fun,a,b),注意此处的fun是函数,并且为数值形式,所以使用*、/、^等运算时要加上小数点。 其中,qudal(fun,a,b,…)为用高精度进行计算,在同样的精度下高阶方法quadl要求的节点较少,效率可能比quad更高。 例2 >> Q=quad('1./(x.^3-2*x-5)',0,2) Q = -0.4605 >> F=inline('1./(x.^3-2*x-5)'); >> S1=quad(F,0,2) S1 = -0.4605 命令:trapz(),梯形法求数值积分。 格式:trapz(x,y)。 其中x为带有步长的积分区间,y为数值形式的运算(相当于上面介绍的函数fun)。 例3 >> x=-1:0.1:1; >> y=1./(1+25*x.^2); >> T=trapz(x,y) T = 0.5492 例4 计算 >> x=0:pi/100:pi; >> y=sin(x); >> trapz(x,y) ans = 1.9998 注:梯形法是步长越小精度越高。 命令:dblquad(),抛物线法求二重数值积分。 格式:dblquad(fun,xmin,xmax,ymin,ymax),fun可以用inline定义,也可以通过某个函数文件的句柄传递。 例5 >> Q1=dblquad(inline('y*sin(x)'),pi,2*pi,0,pi) Q1 = -9.8696 >> Q2=dblquad(inline('y*sin(x)'),0,pi,pi,2*pi) Q2 = 29.6088 >> Q3=dblquad(@integrnd,pi,2*pi,0,pi) Q3 = -9.8696 编写M文件: % integrnd.m function z=integrnd(x,y) z=y*sin(x); 命令:fprintf(文件地址,格式,写入的变量),把数据写入指定文件。 例6 >> x=0:.1:1; >> y=[x;exp(x)]; >> % 打开文件 >> fid=fopen('vexp.txt','w'); >> % 写入 >> fprintf(fid,'%6.2f%12.8f\n',y); >> % 关闭文件 >> fclose(fid); 命令:int(f,v,a,b),求f关于v积分,积分区间有[a,b]。 命令:subs(f,’x’,a),将a赋值给符号表达式f中的x,并计算出值。 计算实验:计算定积分近似值 1、矩形法计算定积分近似值 取f(x)=exp(x),求定积分int(exp(x),0,1)的近似值。 >> % 积分区间为[0,1],等距划分为20个子区间。 >> x=linspace(0,1,21); >> % 选取每个子区间的端点并计算端点处函数值 >> % 取区间的左端点处得函数值乘以区间长度全部加起来。 >> y1=y(1:20); >> s1=sum(y1)/20 s1 = 1.6757 若选取右端点: >> y2=y(2:21); >> s2=sum(y2)/20 s2 = 1.7616 2、编程用矩形法计算定积分的近似值 根据 f=input('请输入被积函数f(x)='); qujian=input('区间[a,b]='); n=input('请输入子区间个数n='); s=0; for i=1:n x=qujian(1)+(qujian(2)-qujian(1))/n*i; y=eval(f); s=s+y; end disp('定积分的近似值是:'); s=s*(qujian(2)-qujian(1))/n 存为juxingfa.m。运行结果为: >> juxingfa 请输入被积函数f(x)='x^2' 区间[a,b]=[0,1] 请输入子区间个数n=10 定积分的近似值是: s = 0.3850 >> juxingfa 请输入被积函数f(x)='x^2' 区间[a,b]=[0,1] 请输入子区间个数n=1000 定积分的近似值是: s = 0.3338 注:可见子区间个数较少时精度不够高。 3、编程用梯形法计算定积分的近似值 根据公式 f=input('请输入被积函数f(x)='); qujian=input('区间[a,b]='); n=input('请输入子区间个数n='); s=0; for i=1:n x=qujian(1)+(qujian(2)-qujian(1))/n*i; y=eval(f); s=s+y; end x=qujian(1); y=eval(f); s=s+y/2; x=qujian(2); y=eval(f); s=s+y/2; disp('定积分的近似值是:'); s=s*(qujian(2)-qujian(1))/n 存为tixingfa.m,运行结果如下: 请输入被积函数f(x)='x^2' 区间[a,b]=[0 1] 请输入子区间个数n=10 定积分的近似值是: s = 0.4350 4、导数、单调性与极值 导数是分析充分光滑的函数f(x)的单调性和极值点的有效方法。函数在x0点达到局部极大(或局部极小)的充分条件是f’(x0)=0且f’’(x0)0)。 命令:fminbnd或fminsearch,求极小值 >> [x,f]=fminsearch(inline(fun),0.5) |
CopyRight 2018-2019 实验室设备网 版权所有 |