《实验积分计算幻灯片.ppt》由会员分享，可在线阅读，更多相关《实验积分计算幻灯片.ppt（55页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、实验积分计算第1页，共55页，编辑于2022年，星期五实验目的：1.通过实验加深理解积分理论中分割、近似、求和、取极限的思想方法；2.学习并掌握matlab求不定积分、定积分、二重积分、曲线积分的方法；3.学习matlab命令sum，symsum，int。第2页，共55页，编辑于2022年，星期五l学习 Matlab 命令l计算不定积分l定积分的概念l不定积分和广义积分l二重积分的计算l曲线积分实验内容：第3页，共55页，编辑于2022年，星期五1.学习Matlab命令1).求和命令 sum 调用格式：sum(x)给出向量x的各元素的累加和。若x为矩阵，则是一个元素为每列列和的行向量。例1x=
2、1,2,3,4,5,6,7,8,9,10;sum(x)ans=55第4页，共55页，编辑于2022年，星期五例2x=1,2,3;4,5,6;7,8,9 x=1 2 3 4 5 6 7 8 9sum(x)ans=12 15 18第5页，共55页，编辑于2022年，星期五symsum(s,n)symsum(s,k,m,n)当x的元素很有规律，例如为s(k)时，可用symsum求得x的各项和，即 2).求和命令 symsum 调用格式：symsum(s(k),1,n)=s(1)+s(2)+.+s(n)symsum(s(k),k,m,n)=s(m)+s(m+1)+.+s(n)第6页，共55页，编辑于2
3、022年，星期五例3.syms k nsymsum(k,1,10)ans=55symsum(k2,k,1,n)ans=1/3*(n+1)3-1/2*(n+1)2+1/6*n+1/6第7页，共55页，编辑于2022年，星期五int(f(x)int(f(x)int(f(x,y),x)int(f(x,y),x)3).积分命令 int 调用格式：int(f(x),a,b)int(f(x),a,b)int(f(x,y),x,a,b)int(f(x,y),x,a,b)第8页，共55页，编辑于2022年，星期五2.计算不定积分例4 解：输入命令：syms x;y=x2*log(x);y1=diff(y);%
4、对函数求导对函数求导y0=int(y1);y2=int(y);y2=1/3*x3*log(x)-1/9*x3第9页，共55页，编辑于2022年，星期五例5 解：输入命令：syms x a;y=sqrt(a2-x2),(x+1)/(3*x-1)(1/3),x2*asin(x);int(y,x);第10页，共55页，编辑于2022年，星期五ans=1/2*x*(a2-x2)(1/2)+1/2*a2*asin(1/a2)(1/2)*x),1/15*(3*x-1)(5/3)-1/3*(3*x-1)(2/3),1/3*x3*asin(x)+1/9*x2*(1-x2)(1/2)+2/9*(1-x2)(1/
5、2)第11页，共55页，编辑于2022年，星期五syms a x;f=simple(int(x3*(cos(a*x)2);f=1/16*(4*a3*x3*sin(2*a*x)+2*x4*a4+6*x2*a2*cos(2*a*x)-6*a*x*sin(2*a*x)+3-3*cos(2*a*x)/a4 第12页，共55页，编辑于2022年，星期五f1=x4/8+(x3/(4*a)-3*x/(8*a3)*sin(2*a*x)+(3*x2/(8*a2)-3/(16*a4)*cos(2*a*x);r=simple(f-f1)r=3/16/a4第13页，共55页，编辑于2022年，星期五3 定积分的概念x
6、=linspace(0,1,21);定积分为一个和式极限，取f(x)=exp(x)，积分区间为0,1，等距划分为20个子区间y=exp(x);选取每个子区间的的端点，计算端点处的函数值y1=y(1:20);s1=sum(y1)/20;取区间的左端点乘以区间长度全部加起来s1=1.6757第14页，共55页，编辑于2022年，星期五y2=y(2:21);s2=sum(y2)/20;取区间的右端点乘以区间长度全部加起来s2=1.7616plot(x,y);hold onfor i=1:20fill(x(i),x(i+1),x(i+1),x(i),x(i),0,0,y(i),y(i),0,b)end
7、第15页，共55页，编辑于2022年，星期五第16页，共55页，编辑于2022年，星期五for i=1:20fill(x(i),x(i+1),x(i+1),x(i),x(i),0,0,y(i+1),y(i+1),0,r)end若取右端点，则第17页，共55页，编辑于2022年，星期五第18页，共55页，编辑于2022年，星期五从图上可以看出：当点取得越来越多时，s2-s1的值会越来越小，可试取50个点计算，看结果如何。下面按等分区间计算syms k n;s=symsum(exp(k/n)/n,k,1,n);limit(s,n,inf)ans=exp(1)-1 第19页，共55页，编辑于2022
8、年，星期五结果与上面一样。例7解：输入指令syms x;I=int(exp(x),0,1)得结果:I=exp(1)-14 计算定积分和广义积分第20页，共55页，编辑于2022年，星期五例8解：输入指令syms x;int(abs(x-1),0,2)ans=1.第21页，共55页，编辑于2022年，星期五int()int()还可以求函数边界的定积分问题。还可以求函数边界的定积分问题。syms x t;f=(-2*(x2)+1)/(2*(x2)-3*x+1)2;I=simple(int(f,x,cos(t),exp(-2*t)第22页，共55页，编辑于2022年，星期五I=-(exp(2*t)*
9、cos(t)-1)*(-2*cos(t)+exp(2*t)/(-2+exp(2*t)/(-1+exp(2*t)/(2*cos(t)-1)/(cos(t)-1)第23页，共55页，编辑于2022年，星期五例9解：对第一个积分输入指令：syms x;syms p real;int(1/xp,x,1,inf)ans=limit(-(x-exp(p*log(x)/(p-1)/exp(p*log(x),x=inf)第24页，共55页，编辑于2022年，星期五有结果看出，当p1时，ans=1/(p-1).syms x;int(1/(x-1)2,0,2)ans=inf对第二个积分，输入命令：syms x;i
10、nt(1/(2*pi)(1/2)*exp(-x2/2),-inf,inf)对第三个积分，输入命令：ans=7186705221432913/18014398509481984*2(1/2)*pi(1/2)第25页，共55页，编辑于2022年，星期五例10解：输入指令：syms x t;int(sin(x)/x,0,t)ans=sinint(t)help sinint SININT Sine integral function.SININT(x)=int(sin(t)/t,t,0,x).See also COSINT.Overloaded methods help sym/sinint.m这类积
11、分无法用初等函数或其值来表示。第26页，共55页，编辑于2022年，星期五5 数值积分（1）梯形积分法 -先将积分区间划分为几个小区间，用每个小区间上梯形面积之和作为积分的近似值。z=trapz(x,y)x表示积分区间的离散化向量，y是与x同维数的被积函数向量，z返回积分的近似值.例11.求积分解：我们先用具有较大步长的积分区间的离散化向量求解（高精度近似值为1.49364826）第27页，共55页，编辑于2022年，星期五 clear;x=-1:0.5:1;y=exp(-x.2);z=trapz(x,y)（高精度近似值为1.49364826）z=1.46274050 xx=-1:0.05:1
12、;yy=exp(-xx.2);z=trapz(xx,yy)z=1.49334167 plot(x,y,r-o,xx,yy)第28页，共55页，编辑于2022年，星期五第29页，共55页，编辑于2022年，星期五练习：用梯形法计算积分第30页，共55页，编辑于2022年，星期五（2）Simpson公式 -f 在a，b区间上积分的Simpson公式。用通过三点(a,f(a),(a+b)/2,f(a+b)/2),(b,f(b)的抛物线围成的曲边梯形的面积来代替由 f 围成的曲边梯形的面积，由此获得积分的近似值.-将区间a，b 分为n等分，步长为h=(ba)/n，分点为xk=a+kh(k=0,.,n)
13、。在每个小区间xk,xk+1上采用辛普森公式，得到小区间上 f 积分的近似值。对每个小区间上积分的近似值求和，得到f 在a，b区间积分的近似值，这就是复合辛普森公式：第31页，共55页，编辑于2022年，星期五第32页，共55页，编辑于2022年，星期五(3)adaptive Simpson公式 -通常被积函数f 在整个区间a，b 上变化不是很均匀的，如在某点附近函数变化非常急剧，而在其余地方上的变化比较平缓。这种情况用等距剖分小区间的复合求积公式不很合适。为了使积分达到预定的精度又要节省工作量，则可以在函数变化急剧的部分增多节点，即子区间分得细，而在函数变化平缓的地方减少节点，即子区间分得大
14、，这个方法称为自适应积分法(adaptive quadrature methods)第33页，共55页，编辑于2022年，星期五 -对于自适应辛普森方法，采用逐次将区间二等分的方法。为写法统一，将区间a，b记为a，a+h，其中h=ba为区间的长度，称为0级区间。在区间a,a+h上采用辛普森公式，把结果记作 将区间分成两个相等的子区间a,a+h/2和a+h/2,a+h,这两个子区间称为1级子区间，其长度为h/2.在每个子区间上 采用辛普森公式计算积分，然后相加并令第34页，共55页，编辑于2022年，星期五再将1级子区间中的一个或所有的两个二等分，所得的子区间称为2级子区间，其长度为h/22,.
15、。如此继续下去，最后将区间 a,a+h分成n个子区间ai,ai+1(i=0,1,.,n 1)，这样有子区间的长度一般是不同的，如果子区间 ai,ai+1是r级，则其长度为实际上，区间的划分是根据被积函数的变化而定的，函数变化平缓的地方，子区间大，函数变化急剧的地方，子区间小。第35页，共55页，编辑于2022年，星期五设Sa,a+h表示积分I(f)的近似值，那么有 第36页，共55页，编辑于2022年，星期五 -自适应辛普森方法的matlab的命令为 z=quad(f,a,b,tol,trace)f(x)为被积函数a 为积分下限 b 为积分上限tol 为计算精度，缺省为0.001trace非0
16、时，以动态点图的形式实现积分的整个过程。-注意：调用quad函数时，先要建立一个描述被积函数f(x)的函数文件或语句函数。第37页，共55页，编辑于2022年，星期五例12.用自适应辛普森方法求积分解：使用语句函数（内联函数）（高精度近似值为1.49364826）clear;g=inline(exp(-x.2);z=quad(g,-1,1)z=quad(g,-1,1)z=1.49364827第38页，共55页，编辑于2022年，星期五例12.用自适应辛普森方法求积分（高精度近似值为1.49364826）clear;g=inline(exp(-x.2);z=quad(g,-1,1,1e-6)z=
17、1.49364828第39页，共55页，编辑于2022年，星期五例12.用自适应辛普森方法求积分解：建立函数文件（高精度近似值为1.49364826）function g=li12(x)g=exp(-x.2);z=quad(li12,-1,1)z=1.49364828第40页，共55页，编辑于2022年，星期五练习：用自适应辛普森方法法计算积分第41页，共55页，编辑于2022年，星期五6 二重积分解：输入指令：syms x y;int(int(x*y,y,2*x,x2+1),x,0,1)ans=1/12 第42页，共55页，编辑于2022年，星期五syms x y;int(int(sin(p
18、i*(x2+y2),y,-sqrt(1-x2),sqrt(1-x2),x,-1,1)解：积分区域可用不等式表成：二重积分可化为二次积分第43页，共55页，编辑于2022年，星期五ans=sum(1/2*2(-3-2*_k1)*(-1)_k1*pi(2+2*_k1)/(1+_k1)*hypergeom(1/4,3/4,1/2,2+_k1,3/2+_k1,-1/4*pi2)*2(1/2)/pi(1/2)*2(3/2+2*_k1)/(3/4+_k1)*gamma(5/2+2*_k1)/gamma(2+2*_k1)2+1/2*2(-2-2*_k1)*(-1)_k1*pi(2+2*_k1)*hyperg
19、eom(3/4,5/4,3/2,2+_k1,3/2+_k1,-1/4*pi2)*2(1/2)/pi(1/2)*2(2*_k1-1/2)/(1/4+_k1)/_k1*gamma(3/2+2*_k1)/gamma(2*_k1)/gamma(3+2*_k1),_k1=0.inf)第44页，共55页，编辑于2022年，星期五syms r a;int(int(r*sin(pi*r2),r,0,1),a,0,2*pi)输入命令：结果中仍带有 int，表明 matlab 求不出这一积分，所以采用极坐标化为二次积分：ans=2第45页，共55页，编辑于2022年，星期五7.数值计算矩形区域上二重积分 z=db
20、lquad(f,a,b,c,d)求得二元函数f(x,y)的重积分，其中a,b为变量 x 的积分下、上限；c,d 为变量 y 的积分下、上限;z=tripquad(fun,a,b,c,d,e,f)求得three元函数fun(x,y,z)的重积分，其中a,b为变量 x 的积分下、上限；c,d 为变量 y 的积分下、上限;e,f 为变量 z 的积分下、上限第46页，共55页，编辑于2022年，星期五例15.计算重积分解：fun=inline(x.*exp(x.2+y.2),x,y);dblquad(fun,0,2,-2,2)ans=8.818304115675463e+002第47页，共55页，编辑
21、于2022年，星期五例16.计算重积分解：fun=inline(y.*sin(x)+z.*cos(x),x,y,z)triplequad(fun,0,pi,0,1,-1,1)ans=1.99999999436264第48页，共55页，编辑于2022年，星期五8 曲线积分例17解：曲线的参数方程为曲线积分可化为：第49页，共55页，编辑于2022年，星期五输入命令ans=1/2syms t;int(cos(t)*sin(t),0,pi/2)第50页，共55页，编辑于2022年，星期五当曲线由参数方程给出当曲线由参数方程给出则曲线长为则曲线长为做函数文件做函数文件pmqxhc.mpmqxhc.mfunction y=pmqxhc(x,y,t,a,b)y=int(sqrt(diff(x,t)2+diff(y,t)2),a,b);第51页，共55页，编辑于2022年，星期五例18输入命令：输入命令：syms t a;x=a*(t-sin(t);y=a*(1-cos(t);s=pmqxhc(x,y,t,0,2*pi)pretty(simple(s)8 a第52页，共55页，编辑于2022年，星期五作业4：第53页，共55页，编辑于2022年，星期五第54页，共55页，编辑于2022年，星期五第55页，共55页，编辑于2022年，星期五