定积分的近似计算(数学实验报告matlab版)
matlab实验报告--定积分的近似计算
abs((inum2-integrate)/integrate))
fprintf('the relative error between inum3 and real-value is about: %g\n\n',...
abs((inum3-integrate)/integrate)) 【调试结果】
○2 使用函数 quad()
quad('sin(x)./x',0,inf) 【调试结果】 ans =
NaN
○3 程序法
%矩阵法
format long
n=inf;a=0;b=inf;
syms x fx
fx=sin(x)./x;
i=1:n;
xj=a+(i-1)*(b-a)/n; xi=a+i*(b-a)/n;
实验目的:
本实验将主要研究定积分的三种近似计算算法:矩形法、梯形法、抛物线法。对于定 积分的近似数值计算,Matlab 有专门函数可用。
实验原理与数学模型:
1. 矩形法 根据定积分的定义,每一个积分和都可以看作是定积分的一个近似值,即
在几何意义上,这是用一系列小矩形面积近似小曲边梯形的结果,所以把这个近似计 算方法称为矩形法.不过,只有当积分区间被分割得很细时,矩形法才有一定的精确度.
【调试结果】
inum =
0.78539816339745
the relative error between inum and real-value is about: 2.82716e-016
【情况记录】
1、梯形法和抛物线法程序设计较为顺利。但要注意使用 for 循环函数和求和函数时
的不同 matlab 命令,避免混淆出错。使用函数 trapz(),quad()时要注意被积函数是数 值形式,应使用数组计算,应用点除即 ./ ,否则将出错,不能调试出结果。
MATLAB实验三 定积分的近似计算
实验三定积分的近似计算一、问题背景与实验目的利用牛顿—莱布尼兹公式虽然可以精确地计算定积分的值,但它仅适用于被积函数的原函数能用初等函数表达出来的情形.如果这点办不到或者不容易办到,这就有必要考虑近似计算的方法.在定积分的很多应用问题中,被积函数甚至没有解析表达式,可能只是一条实验记录曲线,或者是一组离散的采样值,这时只能应用近似方法去计算相应的定积分.本实验将主要研究定积分的三种近似计算算法:矩形法、梯形法、抛物线法.对于定积分的近似数值计算,Matlab有专门函数可用.二、相关函数(命令)及简介1.sum(a):求数组a的和.2.format long:长格式,即屏幕显示15位有效数字.(注:由于本实验要比较近似解法和精确求解间的误差,需要更高的精度).3.double():若输入的是字符则转化为相应的ASCII码;若输入的是整型数值则转化为相应的实型数值.4.quad():抛物线法求数值积分.格式:quad(fun,a,b) ,注意此处的fun是函数,并且为数值形式的,所以使用*、/、^等运算时要在其前加上小数点,即.*、./、.^等.例:Q = quad('1./(x.^3-2*x-5)',0,2);5.trapz():梯形法求数值积分.格式:trapz(x,y)其中x为带有步长的积分区间;y为数值形式的运算(相当于上面介绍的函数fun)例:计算0sin()dx xπ⎰x=0:pi/100:pi;y=sin(x);trapz(x,y)6.dblquad():抛物线法求二重数值积分.格式:dblquad(fun,xmin,xmax,ymin,ymax),fun可以用inline定义,也可以通过某个函数文件的句柄传递.例1:Q1 = dblquad(inline('y*sin(x)'), pi, 2*pi, 0, pi)顺便计算下面的Q2,通过计算,比较Q1 与Q2结果(或加上手工验算),找出积分变量x、y的上下限的函数代入方法.Q2 = dblquad(inline('y*sin(x)'), 0, pi, pi, 2*pi)例2:Q3 = dblquad(@integrnd, pi, 2*pi, 0, pi)这时必须存在一个函数文件integrnd.m:function z = integrnd(x, y) z = y*sin(x);7.fprintf (文件地址,格式,写入的变量):把数据写入指定文件.例:x = 0:.1:1; y = [x; exp(x)];fid = fopen('exp.txt','w'); %打开文件 fprintf(fid,'%6.2f %12.8f\n',y); %写入 fclose(fid) %关闭文件 8.syms 变量1 变量2 …:定义变量为符号. 9.sym('表达式'):将表达式定义为符号.解释:Matlab 中的符号运算事实上是借用了Maple 的软件包,所以当在Matlab 中要对符号进行运算时,必须先把要用到的变量定义为符号. 10.int(f,v,a,b):求f 关于v 积分,积分区间由a 到b .11.subs(f ,'x',a):将 a 的值赋给符号表达式 f 中的 x ,并计算出值.若简单地使用subs(f),则将f 的所有符号变量用可能的数值代入,并计算出值.三、实验内容1. 矩形法根据定积分的定义,每一个积分和都可以看作是定积分的一个近似值,即1()d ()nbi i ai f x x f x ς==∆∑⎰在几何意义上,这是用一系列小矩形面积近似小曲边梯形的结果,所以把这个近似计算方法称为矩形法.不过,只有当积分区间被分割得很细时,矩形法才有一定的精确度.针对不同i ς的取法,计算结果会有不同,我们以 120d 1xx +⎰为例(取100=n ),(1) 左点法:对等分区间b x i n ab a x x a x n i =<<-+=<<<=ΛΛ10,在区间],[1i i x x -上取左端点,即取1-=i i x ς,12 01d ()1ni i i xf x x ς==∆≈+∑⎰0.78789399673078, 理论值 12 0d 14x x π=+⎰,此时计算的相对误差0.7878939967307840.0031784ππ-=≈(2)右点法:同(1)中划分区间,在区间],[1i i x x -上取右端点,即取i i x =ς,12 01d ()1ni i i xf x x ς==∆≈+∑⎰0.78289399673078, 理论值 12 0d 14x x π=+⎰,此时计算的相对误差 0.7828939967307840.0031884ππ-=≈(3)中点法:同(1)中划分区间,在区间1[,]i i x x -上取中点,即取12i ii x x ς-+=, 12 01d ()1ni i i xf x x ς==∆≈+∑⎰0.78540024673078, 理论值 12 0d 14x x π=+⎰,此时计算的相对误差 60.7854002467307842.653104ππ--=≈⨯如果在分割的每个小区间上采用一次或二次多项式来近似代替被积函数,那么可以期望得到比矩形法效果好得多的近似计算公式.下面介绍的梯形法和抛物线法就是这一指导思想的产物.2. 梯形法等分区间b x i n a b a x x a x n i =<<-+=<<<=ΛΛ10,nab x -=∆ 相应函数值为n y y y ,,,10Λ(n i x f y i i ,,1,0),(Λ==).曲线)(x f y =上相应的点为n P P P ,,,10Λ(n i y x P i i i ,,1,0),,(Λ==)将曲线的每一段弧i i P P 1-用过点1-i P ,i P 的弦i i P P 1-(线性函数)来代替,这使得每个],[1i i x x -上的曲边梯形成为真正的梯形,其面积为x y y ii ∆⨯+-21,n i ,,2,1Λ=. 于是各个小梯形面积之和就是曲边梯形面积的近似值,11 11()d ()22nnbi i i i ai i y y x f x x x y y --==+∆≈⨯∆=+∑∑⎰, 即11 ()d ()22bn n ay y b a f x x y y n --≈++++⎰L , 称此式为梯形公式.仍用 12 0d 1x x +⎰的近似计算为例,取100=n ,10112 0d ()122n n y y x b a y y x n --≈++++=+⎰L 0.78539399673078, 理论值 12 0d 14x x π=+⎰,此时计算的相对误差 60.7853939967307845.305104ππ--=≈⨯很显然,这个误差要比简单的矩形左点法和右点法的计算误差小得多.3. 抛物线法由梯形法求近似值,当)(x f y =为凹曲线时,它就偏小;当)(x f y =为凸曲线时,它就偏大.若每段改用与它凸性相接近的抛物线来近似时,就可减少上述缺点,这就是抛物线法.将积分区间],[b a 作n 2等分,分点依次为b x i n a b a x x a x n i =<<-+=<<<=2102ΛΛ,nab x 2-=∆, 对应函数值为n y y y 210,,,Λ(n i x f y i i 2,,1,0),(Λ==),曲线上相应点为n P P P 210,,,Λ(n i y x P i i i 2,,1,0),,(Λ==).现把区间],[20x x 上的曲线段)(x f y =用通过三点),(000y x P ,),(111y x P ,),(222y x P 的抛物线)(12x p x x y =++=γβα来近似代替,然后求函数)(1x p 从0x 到2x 的定积分:21 ()d x x p x x =⎰22 ()d x x x x x αβγ++=⎰)()(2)(30220223032x x x x x x -+-+-γβα]4)(2)()()[(62022022202002γβαγβαγβα++++++++++-=x x x x x x x x x x 由于2201x x x +=,代入上式整理后得 21 ()d x x p x x ⎰)](4)()[(612122202002γβαγβαγβα++++++++-=x x x x x x x x )4(621002y y y x x ++-=)4(6210y y y nab ++-= 同样也有422 ()d x x p x x ⎰)4(6432y y y n ab ++-=……222 ()d n n x nx p x x -⎰)4(621222n n n y y y nab ++-=-- 将这n 个积分相加即得原来所要计算的定积分的近似值:22222212 11()d ()d (4)6ii nnbx i i i i ax i i b af x x p x x y y y n---==-≈=++∑∑⎰⎰, 即021******* ()d [4()2()]6bn n n ab af x x y y y y y y y y n---≈++++++++⎰L L 这就是抛物线法公式,也称为辛卜生(Simpson )公式.仍用 12 0d 1x x +⎰的近似计算为例,取100=n ,102132124222 0d [4()2()]16n n n x b ay y y y y y y y x n ---≈+++++++++⎰L L=0.78539816339745,理论值 12 0d 14x x π=+⎰,此时计算的相对误差 160.7853981633974542.827104ππ--=≈⨯4. 直接应用Matlab 命令计算结果(1) 数值计算 120d .1xx +⎰ 方法1:int('1/(1+x^2)','x',0,1) (符号求积分)方法2:quad('1./(1+x.^2)',0,1) (抛物线法求数值积分)方法3:x=0:0.001:1; y=1./(1+x.^2);trapz(x,y) (梯形法求数值积分) (2)数值计算 212 01d d x x y y -+⎰⎰方法1:int(int('x+y^2','y',-1,1),'x',0,2) (符号求积分)方法2:dblquad(inline('x+y^2'),0,2,-1,1) (抛物线法二重数值积分)四、自己动手1. 实现实验内容中的例子,即分别采用矩形法、梯形法、抛物线法计算 120d 1xx +⎰,取258=n ,并比较三种方法的精确程度.2. 分别用梯形法与抛物线法,计算 2 1d xx⎰,取120=n .并尝试直接使用函数trapz()、quad()进行计算求解,比较结果的差异.3. 试计算定积分 0sin d xx x+∞⎰.(注意:可以运用trapz()、quad()或附录程序求解吗?为什么?)4. 将 120d 1xx +⎰的近似计算结果与Matlab 中各命令的计算结果相比较,试猜测Matlab 中的数值积分命令最可能采用了哪一种近似计算方法?并找出其他例子支持你的观点.5. 通过整个实验内容及练习,你能否作出一些理论上的小结,即针对什么类型的函数(具有某种单调特性或凹凸特性),用某种近似计算方法所得结果更接近于实际值?6. 学习fulu2sum.m 的程序设计方法,尝试用函数 sum 改写附录1和附录3的程序,避免for 循环.五、附录附录1:矩形法(左点法、右点法、中点法)(fulu1.m ) format long n=100;a=0;b=1;inum1=0;inum2=0;inum3=0; syms x fx fx=1/(1+x^2); for i=1:nxj=a+(i-1)*(b-a)/n; %左点 xi=a+i*(b-a)/n; %右点 fxj=subs(fx,'x',xj); %左点值fxi=subs(fx,'x',xi); %右点值fxij=subs(fx,'x',(xi+xj)/2); %中点值inum1=inum1+fxj*(b-a)/n;inum2=inum2+fxi*(b-a)/n;inum3=inum3+fxij*(b-a)/n;endinum1inum2inum3integrate=int(fx,0,1)integrate=double(integrate)fprintf('The relative error between inum1 and real-value is about: %d\n\n',...abs((inum1-integrate)/integrate))fprintf('The relative error between inum2 and real-value is about: %d\n\n',...abs((inum2-integrate)/integrate))fprintf('The relative error between inum3 and real-value is about: %d\n\n',...abs((inum3-integrate)/integrate))附录2:梯形法(fulu2.m)format longn=100;a=0;b=1;inum=0;syms x fxfx=1/(1+x^2);for i=1:nxj=a+(i-1)*(b-a)/n;xi=a+i*(b-a)/n;fxj=subs(fx,'x',xj);fxi=subs(fx,'x',xi);inum=inum+(fxj+fxi)*(b-a)/(2*n);endinumintegrate=int(fx,0,1)integrate=double(integrate)fprintf('The relative error between inum and real-value is about: %d\n\n',...abs((inum-integrate)/integrate))附录2sum:梯形法(fulu2sum.m),利用求和函数,避免for 循环format longn=100;a=0;b=1;syms x fxfx=1/(1+x^2);i=1:n;xj=a+(i-1)*(b-a)/n; %所有左点的数组xi=a+i*(b-a)/n; %所有右点的数组fxj=subs(fx,'x',xj); %所有左点值fxi=subs(fx,'x',xi); %所有右点值f=(fxi+fxj)/2*(b-a)/n; %梯形面积inum=sum(f) %加和梯形面积求解integrate=int(fx,0,1)integrate=double(integrate)fprintf('The relative error between inum and real-value is about: %d\n\n',...abs((inum-integrate)/integrate))附录3:抛物线法(fulu3.m)format longn=100;a=0;b=1;inum=0;syms x fxfx=1/(1+x^2);for i=1:nxj=a+(i-1)*(b-a)/n; %左点xi=a+i*(b-a)/n; %右点xk=(xi+xj)/2; %中点fxj=subs(fx,'x',xj);fxi=subs(fx,'x',xi);fxk=subs(fx,'x',xk);inum=inum+(fxj+4*fxk+fxi)*(b-a)/(6*n);endinumintegrate=int(fx,0,1)integrate=double(integrate)fprintf('The relative error between inum and real-value is about: %d\n\n',...abs((inum-integrate)/integrate))。
数值积分的方法计算定积分,matlab实验
n 1 n 1 h f ( a ) 2 f ( x ) 4 f ( x 1 ) f (b) i i 6 i 1 i 0 2
编程如下: a=0; b=1; %积分上下限分别为 0,1, 将[0,1]区间 3 等分保证误差小于 e-5 n=3; h=(b-a)/n; %h 为区间长度 m=0;n=0; %定义变量并初始化 for x1=a+h:h:b-h m=m+exp(-x1); %计算以步长为 h 的所有节点函数值之和 end for x2=a+h/2:h:b-h/2 %计算以步长为 h/2 的所有节点函数值之和 n=n+exp(-x2); end S=h*(exp(0)+2*m+4*n+exp(-1))/6*2/sqrt(pi) %复化 simpson 公式求积分并输出 format long;
姓 名
学 号
班 级
成 绩
教师姓名:
时间:
( 教 师 填 写 )
实 验 报 告 要 求
如果步骤较多,请自行加页(A4 幅面)ຫໍສະໝຸດ 2e0
1
x
dx .
(1)复化梯形公式 分析: 由复化梯形公式可知, 余项 RTn
b a 2 '' '' x h f ( ) , 且 a=0, b=1, f ( x) e , 12
h
ba 1 1 x 1 x 1 e e 105 ,所以可 , 所以可得 RTn 2 2 2 n 12 n 12n 12n
n 1 h 得 n 91.29 ,所以 n 取 92, Tn [ f (a) 2 f ( xi ) f (b)] 。 2 i 1
程序: n=92,a=0,b=1;
matlab定积分及应用
实验四 定积分及应用实验的目的1、掌握利用Matlab 进行积分运算;2、掌握积分在计算面积、体积等问题中的应用;3、掌握各种积分指令的区别与特点。
实验的基本理论与方法1、定积分定义:函数)(x f 在区间],[b a 上的定积分定义为:设函数)(x f 在],[b a 上有界,在区间],[b a 上任取1-n 个分点:b x x x x x a n n =<<<<<=-1210 ,把],[b a 分成n 个小区间],[1i i i x x -=∆, n i ,,2,1 =。
这些分点构成对区间],[b a 的一个分割,用T 表示。
小区间i ∆的长度为1--=∆i i i x x x 。
记{}i ni x T ∆=≤≤1ma x ,称为分割T 的模。
在区间],[1i i ix x -=∆上取点i ξ)(1i i i x x ≤≤-ξ,做函数值)(i f ξ与小区间长度i x ∆的乘积),2,1()(n i x f i i =∆ξ,并作和∑=∆=ni i i x f S 1)(ξ。
当0→T 时,和S 总趋于确定的极限,这时这个极限为函数)(x f 在区间],[b a 上的定积分,记作⎰badx x f )(。
即i ni i T bax f dx x f ∆=∑⎰=→1)(lim )(ξ。
2、定积分的应用①计算平面图形的面积:由连续曲线)0)()((≥=x f x f y ,直线)(,b a b x a x <==及x 轴所围成的曲边梯形面积为⎰=badx x f S )(;②计算旋转体的体积:由连续曲线)(x f y =,直线)(,b a b x a x <==及x 轴所围成的曲边梯形绕x 轴旋转一周所成立体的体积为⎰=badx x f V 2)]([π;③计算平面曲线的弧长:设曲线弧由直线坐标方程))((b x a x f y ≤≤=给出,其中)(x f 在],[b a 上具有一阶连续导数,则曲线弧长dx y l ba⎰'+=21;设曲线弧由参数方程⎩⎨⎧≤≤==)(,)()(βαt t y y t x x 给出,其中)(),(t y t x 在],[βα上具有连续导数,则曲线弧长dt t y t x l ⎰'+'=βα22)()(;设曲线弧由极坐标方程))((βθαθ≤≤=r r 给出,其中)(θr 在],[βα上具有连续导数,则曲线弧长θθθβαd r r l ⎰'+=22)()(。
实验二:定积分的近似计算
实验二:定积分的近似计算实验目标:1、 熟悉MATLAB 的程序设计方法;2、 熟悉MATLAB 下命令文件和函数文件的建立和使用;3、 学习定积分的三种近似计算方法:矩形法、梯形法、辛普生法;4、 理解数值计算的误差分析。
问题背景:求定积分的近似值的数值方法就是用被积函数的有限个抽样值的离散或加权平均近似值代替定积分的值。
求某函数的定积分时,在多数情况下,被积函数的原函数很难用初等函数表达出来, 因此能够借助牛顿-莱布尼兹公式计算定积分的机会是不多的。
另外,许多实际问题中的被积函数往往是列表函数或其他形式的非连续函数,对这类函数的定积分,也不能用不定积分方法求解。
由于以上原因,数值积分的理论与方法一直是计算数学研究的基本课题。
定积分的近似解的简单方法包括:矩形公式、梯形公式和辛普生公式。
根据定积分的定义,每个积分和都可以看做是定积分的近似值,即⎰∑=∆=ba ni i i x f dx x f 1)()(ζ在几何意义上,这是用一系列小块区域的面积近似小曲边梯形的面积。
当然,只有当积分区间被分割的很细时,计算结果才具有一定的精确度。
● 矩形法:设积分区间被等分为若干份,第i 份是由][1+i i x x 表示,则该小块区域面积为:)(*1i i i x f x x -+ 或)(*11++-i i i x f x x 或)(*211++-i i i x f x x● 梯形法:设积分区间被等分为若干份,第i 份是由][1+i i x x 表示,取)(i x f 和)(1+i x f 的加权平均值作为平均高度)(i f ζ,则该小块区域面积为:2)()(*11+++-i i i i x f x f x x ● 辛普生法:设积分区间被等分为若干份,第i 份是由][1+i i x x 表示,中点为21+i x ,取函数)(x f 在i x ,1+i x ,21+i x 这是三点的函数值的加权平均值作为平均高度的近似值,则该小块区域面积为:6)()(4)(*1211+++++-i i i i i x f x f x f x x实验内容:1、 试推导定积分的三种近似计算方法的迭代公式(矩形法、梯形法、辛普生法)。
三用MATLAB实现定积分计算
s=s+feval(f,z1(j))+feval(f,z2(j));
0,2*pi,1000)
end
s=
s=s*h/2;
-267.2458
Gauss-lobatto是改进的高斯积分方法,采取自适应求积方法
三 用MATLAB实现定积分计算: 2 sin xdx 0
⑴ 矩形公式与梯形公式 z1 =
形的公求式积代公数式精。度为对于1,f 辛(x)甫=1森, x公, 式x 2的, x代3,数应精该度有为 3。
节成点立我x,ba下i和们依f面系先(次介x数考11)将绍dfA虑f(x的i(,xx节))是d=使点x1取t代数, (x消数xAb,为1对xaa精f22)(2区/bx,度而21x间)尽使3代等可用Ab入2分2能(fa1,(的1高1x1)即2限计的)f可制(算所得a,的谓2b到n积高确给分斯b定定近2公aA后似t式1,)同A值d。2时t有,x确1代,x定数2
这两种用随机模拟的方式求积分近似值的方法 z=sum(y)*pi/2/n
/2
z=
蒙特卡罗方法
sin xdx
1.0010
0
3、蒙特卡罗方法的通用函数与调用格式
均值估计法
随机投点法 (设0≤ f(x) ≤1)
b
a
f
( x)dx
ba n
n i1
f
(a (b a)ui )
直接调用。这里被积函数为内部函数,无需另外定义。
s=gaussinteg(‘sin', 0, pi/2,1000) s=
1.0000
6000
§2 数值积分应用问题举例4000
2000
0
一 求卫星轨道长度
定积分的近似计算
fxj=subs(fx,'x',xj);%符号运算,替换
fxi=subs(fx,'x',xi);%符号运算,替换
sum=sum+(fxi+fxj)*(b-a)/(2*n);%求矩形面积
end
sum%如果没有这句话就算不出最后的积分
结果:sum =0.6932
(2)抛物线法:
fzz=subs(fx,'x',zz);
f=(fxx+4*fyy+fzz)*(b-a)/(6*n);
s=sum(f)
结果:s=0.7854
实验结果报告及实验总结:
1.实验结果报告:
(1)梯形法:结果:sum =0.6932
(2)抛物线法:结果:sum = 0.6931
(3)trapz()法:结果:ans =0.6931
n=120;
a=1;
b=2;
sum=0;
syms x fx;
fx=1/x;
for i=1:n
xx=a+(2*i-2)*(b-a)/(2*n);%第一点的自变量取值后面是除2*n,不是n
yy=a+(2*i-1)*(b-a)/(2*n);%第二点的自变量取值后面是2*i-1
zz=a+(2*i-0)*(b-a)/(2*n);%复制后没改自变量的名字
且发现trapz()的调试结果与梯形法结果相同,故可猜测该Matlab中的数值积分命令函数trapz()采用了梯形法近似计算方法。
2.实验结果报告:
(1)使用int:结果:ans =pi/2
(2)使用函数trapz():
结果:Maximum variable size allowed by the program is exceeded.
定积分的近似计算97786
数学实验报告实验序号:4 日期:2012 年12 月13 日实验名称定积分的近似计算问题背景描述:利用牛顿—莱布尼兹公式虽然可以精确地计算定积分的值,但它仅适用于被积函数的原函数能用初等函数表达出来的情形.如果这点办不到或者不容易办到,这就有必要考虑近似计算的方法.在定积分的很多应用问题中,被积函数甚至没有解析表达式,可能只是一条实验记录曲线,或者是一组离散的采样值,这时只能应用近似方法去计算相应的定积分.实验目的:本实验将主要研究定积分的三种近似计算算法:矩形法、梯形法、抛物线法。
对于定积分的近似数值计算,Matlab有专门函数可用。
实验原理与数学模型:1.矩形法根据定积分的定义,每一个积分和都可以看作是定积分的一个近似值,即在几何意义上,这是用一系列小矩形面积近似小曲边梯形的结果,所以把这个近似计算方法称为矩形法.不过,只有当积分区间被分割得很细时,矩形法才有一定的精确度.针对不同的取法,计算结果会有不同。
(1)左点法:对等分区间,在区间上取左端点,即取。
(2)右点法:同(1)中划分区间,在区间上取右端点,即取。
(3)中点法:同(1)中划分区间,在区间上取中点,即取。
2.梯形法等分区间,相应函数值为().曲线上相应的点为()将曲线的每一段弧用过点,的弦(线性函数)来代替,这使得每个上的曲边梯形成为真正的梯形,其面积为,.于是各个小梯形面积之和就是曲边梯形面积的近似值,,即,称此式为梯形公式。
3.抛物线法将积分区间作等分,分点依次为,,对应函数值为(),曲线上相应点为().现把区间上的曲线段用通过三点,,的抛物线来近似代替,然后求函数从到的定积分:由于,代入上式整理后得同样也有……将这个积分相加即得原来所要计算的定积分的近似值:,即这就是抛物线法公式,也称为辛卜生(Simpson)公式.实验所用软件及版本:Matlab 7.0主要内容(要点):1.分别用梯形法与抛物线法,计算,取.并尝试直接使用函数trapz()、quad()进行计算求解,比较结果的差异.2.试计算定积分.(注意:可以运用trapz()、quad()或附录程序求解吗?为什么?)3.学习fulu2sum.m的程序设计方法,尝试用函数sum 改写附录1和附录3的程序,避免for 循环。
MATLAB复化梯形法及龙贝格法计算定积分
MATLAB复化梯形法及龙贝格法计算定积分复化梯形法是一种数值积分方法,用于计算定积分的近似值。
该方法的基本思想是将积分区间等分成多个子区间,并在每个子区间上使用梯形公式来进行近似计算。
具体步骤如下:1.将积分区间[a,b]等分成n个子区间,每个子区间的长度为h=(b-a)/n。
2.在每个子区间上,使用梯形公式计算近似积分值。
梯形公式可以表示为:T=(f(x0)+f(x1))*h/2,其中x0和x1分别是子区间的左右边界,f(x)是被积函数。
3.对所有子区间的近似积分值进行求和,得到整个积分区间的近似积分值。
复化梯形法的精度可以通过增加子区间的数量来提高,即使n越来越大,积分值的近似精度也会越来越高。
以下是一个用MATLAB实现复化梯形法计算定积分的示例代码:```matlabh=(b-a)/nresult = 0;for i = 0:n-1x0=a+i*h;x1=a+(i+1)*h;result = result + (f(x0) + f(x1)) * h / 2;endend```接下来,我们来介绍龙贝格法,龙贝格法是一种迭代数值积分方法,用于计算定积分的近似值。
该方法的基本思想是在梯形公式的基础上应用Richardson外推技术,通过逐步加密和外推,提高积分值的精度。
具体步骤如下:1.初始化一个矩阵,矩阵的第一列为复化梯形法的近似积分值。
2.逐列递推计算,每一列的元素为由前一列的元素计算得到。
计算公式为:R(j,k+1)=R(j,k)+(R(j,k)-R(j-1,k))/((4^k)-1)其中,R(j,k)是第j次迭代中计算的近似积分值,k表示第k次迭代。
3.判断是否达到预设的精度要求,如果满足要求,则返回最终近似积分值;否则,继续迭代计算。
以下是一个用MATLAB实现龙贝格法计算定积分的示例代码:```matlabfunction result = romberg(f, a, b, epsilon, max_iter)R = zeros(max_iter, max_iter);h=b-a;R(1,1)=h*(f(a)+f(b))/2;for k = 2:max_iterh=h/2;sum = 0;for i = 1:2^(k-2)x=a+(2*i-1)*h;sum = sum + f(x);endR(k, 1) = R(k-1, 1) / 2 + h * sum;for j = 2:kR(k,j)=R(k,j-1)+(R(k,j-1)-R(k-1,j-1))/((4^(j-1))-1); endif abs(R(k, k) - R(k-1, k-1)) < epsilonresult = R(k, k);return;endendresult = R(max_iter, max_iter);end```这个代码定义了一个名为`romberg`的函数,它接受五个参数:被积函数`f`、积分区间的左边界`a`、积分区间的右边界`b`、精度要求`epsilon`和最大迭代次数`max_iter`。
matlab数值计算 实验报告
matlab数值计算实验报告Matlab数值计算实验报告引言:Matlab是一种强大的数值计算软件,广泛应用于科学和工程领域。
本实验旨在通过实际案例,展示Matlab在数值计算中的应用能力。
本报告将从三个方面进行讨论:数值积分、线性方程组求解和最优化问题。
一、数值积分:数值积分是数学中常见的问题,Matlab提供了多种函数和方法来解决这类问题。
我们以求解定积分为例进行讨论。
假设我们要求解函数f(x) = x^2在区间[0, 1]上的定积分。
我们可以使用Matlab中的quad函数来进行计算,代码如下:```matlabf = @(x) x.^2;integral = quad(f, 0, 1);disp(integral);```运行以上代码,我们可以得到定积分的近似值为0.3333。
通过调整积分方法和精度参数,我们可以得到更精确的结果。
二、线性方程组求解:线性方程组求解是数值计算中的重要问题,Matlab提供了多种函数和方法来解决线性方程组。
我们以一个简单的线性方程组为例进行讨论。
假设我们要求解以下线性方程组:```2x + y = 5x - y = 1```我们可以使用Matlab中的linsolve函数来求解,代码如下:```matlabA = [2 1; 1 -1];B = [5; 1];X = linsolve(A, B);disp(X);```运行以上代码,我们可以得到方程组的解为x = 2,y = 3。
通过调整方程组的系数矩阵和右侧向量,我们可以求解更复杂的线性方程组。
三、最优化问题:最优化问题在科学和工程领域中广泛存在,Matlab提供了多种函数和方法来解决这类问题。
我们以求解无约束最优化问题为例进行讨论。
假设我们要求解函数f(x) = x^2的最小值。
我们可以使用Matlab中的fminunc函数来进行计算,代码如下:```matlabf = @(x) x.^2;x0 = 1; % 初始点options = optimoptions('fminunc', 'Display', 'iter');[x, fval] = fminunc(f, x0, options);disp(x);disp(fval);```运行以上代码,我们可以得到最小值的近似解为x = 0,f(x) = 0。
2022年MATLAB数学实验报告定积分的近似计算
MATLAB数学实验报告实验日期:11月20日实验名称定积分旳近似计算姓名:学号:班级:问题背景描述:运用牛顿—莱布尼兹公式虽然可以精确地计算定积分旳值,但它仅合用于被积函数旳原函数能用初等函数体现出来旳情形.如果这点办不到或者不容易办到,这就有必要考虑近似计算旳措施.在定积分旳诸多应用问题中,被积函数甚至没有解析体现式,也许只是一条实验记录曲线,或者是一组离散旳采样值,这时只能应用近似措施去计算相应旳定积分.实验目旳:本实验将重要研究定积分旳三种近似计算算法:矩形法、梯形法、抛物线法。
对于定积分旳近似数值计算,Matlab有专门函数可用。
实验原理与数学模型:1.矩形法根据定积分旳定义,每一种积分和都可以看作是定积分旳一种近似值,即在几何意义上,这是用一系列小矩形面积近似小曲边梯形旳成果,因此把这个近似计算措施称为矩形法.但是,只有当积分区间被分割得很细时,矩形法才有一定旳精确度.针对不同旳取法,计算成果会有不同。
(1)左点法:对等分区间,在区间上取左端点,即取。
(2)右点法:同(1)中划分区间,在区间上取右端点,即取。
(3)中点法:同(1)中划分区间,在区间上取中点,即取。
2.梯形法等分区间,相应函数值为().曲线上相应旳点为()将曲线旳每一段弧用过点,旳弦(线性函数)来替代,这使得每个上旳曲边梯形成为真正旳梯形,其面积为,.于是各个小梯形面积之和就是曲边梯形面积旳近似值,,即,称此式为梯形公式。
3.抛物线法将积分区间作等分,分点依次为,,相应函数值为(),曲线上相应点为().现把区间上旳曲线段用通过三点,,旳抛物线来近似替代,然后求函数从到旳定积分:由于,代入上式整顿后得同样也有……将这个积分相加即得本来所要计算旳定积分旳近似值:,即这就是抛物线法公式,也称为辛卜生(Simpson)公式.实验所用软件及版本:Matlab 7.0重要内容(要点):1.分别用梯形法与抛物线法,计算,取.并尝试直接使用函数trapz()、quad()进行计算求解,比较成果旳差别.2.试计算定积分.(注意:可以运用trapz()、quad()或附录程序求解吗?为什么?)3.学习fulu2sum.m旳程序设计措施,尝试用函数sum 改写附录1和附录3旳程序,避免for 循环。
matlab蒙特卡洛法求定积分
matlab蒙特卡洛法求定积分摘要:一、引言二、蒙特卡洛法简介三、用MATLAB 实现蒙特卡洛法求定积分的方法四、实例:用MATLAB 实现蒙特卡洛法求定积分五、结论正文:一、引言在数学中,定积分是一种常见的计算方式,它可以用来求解曲线下的面积、长度、体积等。
在实际应用中,有些定积分的具体形式难以求解,这时就可以采用蒙特卡洛法来近似求解。
蒙特卡洛法是一种通过随机抽样来估算数学期望值的方法,其在求解定积分方面有着广泛的应用。
本文将介绍如何在MATLAB 中使用蒙特卡洛法求解定积分。
二、蒙特卡洛法简介蒙特卡洛法是一种基于随机抽样的数值计算方法,通过大量模拟实验来近似求解问题。
在求解定积分时,蒙特卡洛法通过随机生成一定数量的样本点,然后计算这些样本点对应的函数值之和,再乘以样本点的权重,最后得到定积分的近似解。
蒙特卡洛法的优点在于它不需要求解积分的具体形式,只需要知道积分函数的表达式即可。
同时,蒙特卡洛法的精度可以通过增加样本点数量来提高。
三、用MATLAB 实现蒙特卡洛法求定积分的方法在MATLAB 中,可以使用内置的random 函数生成随机数,结合定积分的表达式,就可以实现蒙特卡洛法求解定积分。
具体步骤如下:1.创建一个函数文件,输入定积分的表达式;2.使用random 函数生成一定数量的随机数;3.将随机数代入定积分的表达式,计算每个样本点对应的函数值;4.计算所有样本点函数值之和,再乘以样本点的权重,得到定积分的近似解。
四、实例:用MATLAB 实现蒙特卡洛法求定积分例如,求解定积分∫(0~π) sin x dx。
首先,创建一个函数文件,输入以下代码:```matlabfunction integral = montecarlo_integration(n, a, b, f)% n 为样本点数量,a 和b 为积分区间,f 为积分函数x = linspace(a, b, n+1);y = f(x);integral = mean(y);end```然后,在命令窗口中输入以下命令:```matlab= 10000;a = 0;b = pi;f = @(x) sin(x);integral = montecarlo_integration(n, a, b, f);```最后,得到定积分∫(0~π) sin x dx 的近似解为:```matlabintegral =2.5717```五、结论蒙特卡洛法是一种强大的数值计算方法,可以用来求解各种定积分。
定积分计算实验报告(3篇)
第1篇一、实验目的1. 理解定积分的概念,掌握定积分的计算方法。
2. 熟悉数值积分的方法,提高数值计算能力。
3. 通过实验,验证定积分的计算结果,加深对定积分理论的理解。
二、实验原理定积分是数学分析中的一个基本概念,它表示函数在某一区间上的累积效果。
对于给定的函数f(x),在区间[a, b]上的定积分可以表示为:∫[a, b] f(x) dx其中,dx表示无穷小的区间宽度。
在实际计算中,定积分往往采用数值积分的方法进行近似计算。
三、实验仪器与软件1. 仪器:计算机2. 软件:MATLAB四、实验步骤1. 输入函数表达式:在MATLAB中输入待积分函数的表达式,例如f(x) = x^2。
2. 设置积分区间:设定积分的上下限a和b。
3. 选择数值积分方法:MATLAB提供了多种数值积分方法,如梯形法、辛普森法、高斯法等。
根据需要选择合适的方法。
4. 进行数值积分计算:调用MATLAB的数值积分函数,如quad函数,进行积分计算。
5. 结果分析:观察计算结果,与理论值进行对比,分析误差来源。
五、实验数据及结果1. 函数表达式:f(x) = x^22. 积分区间:[0, 1]3. 数值积分方法:辛普森法4. 计算结果:I ≈ 1.1666666667六、误差分析1. 理论值:∫[0, 1] x^2 dx = [x^3/3] |[0, 1] = 1/32. 误差来源:a. 数值积分方法的误差:由于数值积分方法是一种近似计算方法,其计算结果与真实值存在一定的误差。
b. 计算过程中的舍入误差:在计算过程中,由于计算机的浮点数表示,可能导致舍入误差。
3. 误差分析:计算结果与理论值相差较大,说明数值积分方法的误差较大。
在实际应用中,可以根据需要选择合适的数值积分方法,以减小误差。
七、实验结论1. 通过本次实验,掌握了定积分的计算方法,了解了数值积分的方法及其优缺点。
2. 了解了数值积分方法在计算过程中的误差来源,为实际应用提供了参考。
用递推公式计算定积分(matlab版)
用递推公式计算定积分实验目的:1.充分理解不稳定的计算方法会造成误差的积累,在计算过程中会导致误差的迅速增加,从而使结果产生较大的误差.2.在选择数值计算公式来进行近似计算时,应学会选用那些在计算过程中不会导致误差迅速增长的计算公式.3.理解不稳定的计算公式造成误差积累的来源与具体过程;4.掌握简单的matlab语言进行数值计算的方法.实验题目:对n=0,1,2,…,20,计算定积分:实验原理:由于y<n>= = –在计算时有两种迭代方法,如下:方法一:y<n>=– 5*y<n-1>,n=1,2,3, (20)取y<0>= = ln6-ln5 ≈ 0.182322方法二:利用递推公式:y<n-1>=-*y<n>,n=20,19, (1)而且,由 = * ≤≤* =可取:y<20>≈*<>≈0.008730.实验内容:对算法一,程序代码如下:function [y,n]=funa<>syms k n t;t=0.182322;n=0;y=zeros<1,20>;y<1>=t;for k=2:20y<k>=1/k-5*y<k-1>;n=n+1;endy<1:6>y<7:11>对算法二,程序代码如下:%计算定积分;%n--表示迭代次数;%y用来存储结果;function [y,n]=f<>;syms k y_20;y=zeros<21,1>;n=1;y_20=<1/105+1/126>/2;y<21>=y_20;for k=21:-1:2y<k-1>=1/<5*<k-1>>-y<k>/5;n=n+1;end实验结果:由于计算过程中,前11个数字太小,后9个数字比较大,造成前面几个数字只显示0.0000的现象,所以先输出前6个,再输出7—11个,这样就能全部显示出来了.算法一结果:[y,n]=funa%先显示一y<1>—y<6>ans =0.1823-0.41162.3914-11.706958.7346-293.5063%再显示y<7>—y<11>ans =1.0e+005 *0.0147-0.07340.3669-1.83469.1728y =1.0e+012 *Columns 1 through 11 0.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000Columns 12 through 20 -0.00000.0000-0.00010.0006-0.00290.0143-0.07170.3583-1.7916n = 19算法二结果:>> [y,b]=fy =0.18230.08840.05800.04310.03430.02850.02430.02120.01880.01690.01540.01410.01300.01200.01120.01050.00990.00930.00890.00830.0087b =21实验分析:从两题的计算结果可以看出来,算法一是不稳定的,而算法二是稳定的.对算法一:由于y<1>本身具有一定的误差 ,设为a_1,则由于y<n>=1/n-5y<n-1>=1/n-5<1/<n-1>-5y<n-1>>=……=1/n-5/<n-1>-5^2/<n-2>-…-<5^n>*y<0>所以经过多次迭代后会使误差增大很多倍.由此可知:在实际应用过程中应尽量避免使用数值不稳定的公式.。
matlab计算方法实验报告5(数值积分)
算法分析与源程序( 50%),实验结果及分析(30%),实验报告(20%)
指导老师签名:
t1=t2;s1=s2;
elseif k>=3
r1=r2;
c2=s2+(1/15)*(s2-s1);
r2=c2+(1/63)*(c2-c1);
k=k+1;h=h/2;
t1=t2;s1=s2;
c1=c2;
end
end
实验结果与分析
函数xex在区间[1,2]对x进行积分求值,要求误差为0.5*10-7,并与精确值进行比较。(精确值:7.38905609893065)
计算方法实验报告(5)学生姓名杨邦学号指导教师
吴明芬
实验时间
2014.4.16
地点
综合实验大楼203
实验题目
数值积分方法
实验目的
利用复化梯形、辛普森公式和龙贝格数值积分公式计算定积分的近似植。
实验内容
梯形、辛普森、柯特斯法及其Matlab实现;
变步长的梯形、辛普森、柯特斯法及其Matlab实现。
题目由同学从学习材料中任意选两题
梯形:>> jifeng_tixing(1,2,7019,fun)
ans =7.38905612723022
辛普生:>> jifeng_xingpu(1,2,24,fun)
ans =7.38905612621471
龙贝格:>> jifeng_long(fun,1,2,10e-7)
ans =7.38905609893079
算
法
分
析
与
源
程
序
梯形:function y=jifeng_tixing(a,b,n,fun)
数学实验“几种常见的求积分近似解的方法”实验报告(内含matlab程序)
数学实验“几种常见的求积分近似解的方法”实验报告(内含matlab程序)西京学院数学软件实验任务书实验二十一实验报告一、实验名称:Romberg 积分法,Gauss 型积分法,高斯-勒让德积分法,高斯-切比雪夫积分法,高斯-拉盖尔积分法,高斯-埃尔米特积分法。
二、实验目的:进一步熟悉Romberg 积分法,Gauss 型积分法,高斯-勒让德积分法,高斯-切比雪夫积分法,高斯-拉盖尔积分法,高斯-埃尔米特积分法。
三、实验要求:运用Matlab/C/C++/Java/Maple/Mathematica 等其中一种语言完成程序设计。
四、实验原理:1.Romberg 积分法:龙贝格积分法是用里查森外推算法来加快复合梯形求积公式的收敛速度,它的算法如下,其中()i m T 是通过一系列逼近原定积分的龙贝格分值.计算(0)1[()()]2b aT f a f b -=+ 对1,2,3,k n = ,计算下列各步:21()(1)11111(21)()[()]222k k k k k j b a j b a TT f a ---=---=++∑对1,2,,m k = 和,1,2,,1i k k k =-- ,计算111441m i i i m m m m T T T--+-=-随着计算的步骤的增加,()i mT 越来越逼近积分()ba f x dx ?。
2.Gauss 型积分法:高斯积分公式的思想是用n 个不等距的节点123,,,nx x x x 对被积函数进行插值,然后对插值后的函数进行积分,其积分公式为:111()()nk k k f x dx A f x -=≈∑?如果积分区间不是[1,1]-,则需转换到此区间:11()()222bab a b a b af x dx f t dt ---+=+?其中系数k A 、节点k x 与n 的关系如下表所示: 3.高斯-切比雪夫积分法:第一类切比雪夫积分形式为:11()()nk k k f x dx A f x -=≈∑?其中k A n π=,21cos2k k x nπ-= 4.高斯-拉盖尔积分法:高斯-拉盖尔公式有两种形式:1()()nxk k k e f x dx A f x +∞-=≈∑?1()()k nx k k k f x dx A e f x +∞=≈∑?下面编制的程序是针对第一种形式的高斯-拉盖尔公式,即1()()nxk k k e f x dx A f x +∞-=≈∑?因此程序的第一个输入参数——被积函数,是上式中的()f x 。
《数学实验》实验报告——定积分的近似求解
3
2 梯形法程序如下: f=input('请输入被积函数 f(x)='); qujian=input('请输入积分区间[a,b]='); n=input('请输入子区间个数 n='); s=0; for i=1:n-1 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
《数学实验》实验报告
班级 试验 内容 **** 学号 **** 姓名 试验 类别
自选试验
****
成绩 试验 时间 2011 年 5 月 20 日—22 日
定积分的近似求解
试验问题:
用梯形法与抛物法,通过 MATLAB,计算 x 2 dx 的近似值,取 n=10,比较结果的差异,研究
0 1
定积分的两种近似计算方法。
1 1 1 2 ph 3 6rh h(2 ph 2 6r ) h( y 0 4 y1 y 2 ) 3 3 3 。 ba n ,则上面所求的 S 等于区间 [ x0 , x2 ] 上以抛物线为曲边的曲边梯形的面积。同理可
取
以得到区间 [ xi 1 , xi 1 ] 上以抛物线为曲边的曲边梯形的面积:
试验目的:
通过分别用梯形法与抛物线法计算定积分的近似值, 进而熟练掌握运用 MATLAB 来解决 定积分的近似求解,体会 MATLAB 的强大功能。
matlab定积分的近似计算
arctan x
1 0
π 4
左点法相对误差:0.78789399673078 / 4 0.003178 /4
右点法相对误差:0.78289399673078 / 4 0.003188 /4
中点法相对误差:0.78540024673078 / 4 2.653 10-6 /4
yi
xi
yi f ( xi ), i 1, 2,
,n
整个曲边梯形的面积:
b
S a f ( x)dx
n
Si
i 1
n i 1
yi1 2
yi xi
Si
15
梯形法
如果我们 n 等分区间 [a,b],即令:
x1 x2 xn
h ba n
则 S
的抛物线来近似原函数 f (x) 。
用抛物线代替该直线, 计算精度是否会更好?
18
抛物线法
设过以上三点的抛物线方程为:
y = x2 + x + = p1(x)
则在区间 [x0, x2] 上,有
x2 f (x)dx x0
x2 x0
p1( x)dx
x2 ( x2 x )dx
b
n
f(x )dx
a
lim
max(xi )0
i 1
f(i
)xi
其中 a=x0<x1<…<xn=b,
xi=xi-xi-1, i(xi-1,xi), i=1,2,…,n
若在[a,b]上, F’(x)=f(x), 则
b
f(x )dx F(b) F(a)
MATLAB实验三 定积分的近似计算
实验三定积分的近似计算一、问题背景与实验目的利用牛顿—莱布尼兹公式虽然可以精确地计算定积分的值,但它仅适用于被积函数的原函数能用初等函数表达出来的情形.如果这点办不到或者不容易办到,这就有必要考虑近似计算的方法.在定积分的很多应用问题中,被积函数甚至没有解析表达式,可能只是一条实验记录曲线,或者是一组离散的采样值,这时只能应用近似方法去计算相应的定积分.本实验将主要研究定积分的三种近似计算算法:矩形法、梯形法、抛物线法.对于定积分的近似数值计算,Matlab有专门函数可用.二、相关函数(命令)及简介1.sum(a):求数组a的和.2.format long:长格式,即屏幕显示15位有效数字.(注:由于本实验要比较近似解法和精确求解间的误差,需要更高的精度).3.double():若输入的是字符则转化为相应的ASCII码;若输入的是整型数值则转化为相应的实型数值.4.quad():抛物线法求数值积分.格式:quad(fun,a,b) ,注意此处的fun是函数,并且为数值形式的,所以使用*、/、^等运算时要在其前加上小数点,即.*、./、.^等.例:Q = quad('1./(x.^3-2*x-5)',0,2);5.trapz():梯形法求数值积分.格式:trapz(x,y)其中x为带有步长的积分区间;y为数值形式的运算(相当于上面介绍的函数fun)例:计算0sin()dx xπ⎰x=0:pi/100:pi;y=sin(x);trapz(x,y)6.dblquad():抛物线法求二重数值积分.格式:dblquad(fun,xmin,xmax,ymin,ymax),fun可以用inline定义,也可以通过某个函数文件的句柄传递.例1:Q1 = dblquad(inline('y*sin(x)'), pi, 2*pi, 0, pi)顺便计算下面的Q2,通过计算,比较Q1 与Q2结果(或加上手工验算),找出积分变量x、y的上下限的函数代入方法.Q2 = dblquad(inline('y*sin(x)'), 0, pi, pi, 2*pi)例2:Q3 = dblquad(@integrnd, pi, 2*pi, 0, pi)这时必须存在一个函数文件integrnd.m:function z = integrnd(x, y) z = y*sin(x);7.fprintf (文件地址,格式,写入的变量):把数据写入指定文件.例:x = 0:.1:1; y = [x; exp(x)];fid = fopen('exp.txt','w'); %打开文件 fprintf(fid,'%6.2f %12.8f\n',y); %写入 fclose(fid) %关闭文件 8.syms 变量1 变量2 …:定义变量为符号. 9.sym('表达式'):将表达式定义为符号.解释:Matlab 中的符号运算事实上是借用了Maple 的软件包,所以当在Matlab 中要对符号进行运算时,必须先把要用到的变量定义为符号. 10.int(f,v,a,b):求f 关于v 积分,积分区间由a 到b .11.subs(f ,'x',a):将 a 的值赋给符号表达式 f 中的 x ,并计算出值.若简单地使用subs(f),则将f 的所有符号变量用可能的数值代入,并计算出值.三、实验内容1. 矩形法根据定积分的定义,每一个积分和都可以看作是定积分的一个近似值,即1()d ()nbi i ai f x x f x ς==∆∑⎰在几何意义上,这是用一系列小矩形面积近似小曲边梯形的结果,所以把这个近似计算方法称为矩形法.不过,只有当积分区间被分割得很细时,矩形法才有一定的精确度.针对不同i ς的取法,计算结果会有不同,我们以 120d 1xx +⎰为例(取100=n ),(1) 左点法:对等分区间b x i n ab a x x a x n i =<<-+=<<<=ΛΛ10,在区间],[1i i x x -上取左端点,即取1-=i i x ς,12 01d ()1ni i i xf x x ς==∆≈+∑⎰0.78789399673078, 理论值 12 0d 14x x π=+⎰,此时计算的相对误差0.7878939967307840.0031784ππ-=≈(2)右点法:同(1)中划分区间,在区间],[1i i x x -上取右端点,即取i i x =ς,12 01d ()1ni i i xf x x ς==∆≈+∑⎰0.78289399673078, 理论值 12 0d 14x x π=+⎰,此时计算的相对误差 0.7828939967307840.0031884ππ-=≈(3)中点法:同(1)中划分区间,在区间1[,]i i x x -上取中点,即取12i ii x x ς-+=, 12 01d ()1ni i i xf x x ς==∆≈+∑⎰0.78540024673078, 理论值 12 0d 14x x π=+⎰,此时计算的相对误差 60.7854002467307842.653104ππ--=≈⨯如果在分割的每个小区间上采用一次或二次多项式来近似代替被积函数,那么可以期望得到比矩形法效果好得多的近似计算公式.下面介绍的梯形法和抛物线法就是这一指导思想的产物.2. 梯形法等分区间b x i n a b a x x a x n i =<<-+=<<<=ΛΛ10,nab x -=∆ 相应函数值为n y y y ,,,10Λ(n i x f y i i ,,1,0),(Λ==).曲线)(x f y =上相应的点为n P P P ,,,10Λ(n i y x P i i i ,,1,0),,(Λ==)将曲线的每一段弧i i P P 1-用过点1-i P ,i P 的弦i i P P 1-(线性函数)来代替,这使得每个],[1i i x x -上的曲边梯形成为真正的梯形,其面积为x y y ii ∆⨯+-21,n i ,,2,1Λ=. 于是各个小梯形面积之和就是曲边梯形面积的近似值,11 11()d ()22nnbi i i i ai i y y x f x x x y y --==+∆≈⨯∆=+∑∑⎰, 即11 ()d ()22bn n ay y b a f x x y y n --≈++++⎰L , 称此式为梯形公式.仍用 12 0d 1x x +⎰的近似计算为例,取100=n ,10112 0d ()122n n y y x b a y y x n --≈++++=+⎰L 0.78539399673078, 理论值 12 0d 14x x π=+⎰,此时计算的相对误差 60.7853939967307845.305104ππ--=≈⨯很显然,这个误差要比简单的矩形左点法和右点法的计算误差小得多.3. 抛物线法由梯形法求近似值,当)(x f y =为凹曲线时,它就偏小;当)(x f y =为凸曲线时,它就偏大.若每段改用与它凸性相接近的抛物线来近似时,就可减少上述缺点,这就是抛物线法.将积分区间],[b a 作n 2等分,分点依次为b x i n a b a x x a x n i =<<-+=<<<=2102ΛΛ,nab x 2-=∆, 对应函数值为n y y y 210,,,Λ(n i x f y i i 2,,1,0),(Λ==),曲线上相应点为n P P P 210,,,Λ(n i y x P i i i 2,,1,0),,(Λ==).现把区间],[20x x 上的曲线段)(x f y =用通过三点),(000y x P ,),(111y x P ,),(222y x P 的抛物线)(12x p x x y =++=γβα来近似代替,然后求函数)(1x p 从0x 到2x 的定积分:21 ()d x x p x x =⎰22 ()d x x x x x αβγ++=⎰)()(2)(30220223032x x x x x x -+-+-γβα]4)(2)()()[(62022022202002γβαγβαγβα++++++++++-=x x x x x x x x x x 由于2201x x x +=,代入上式整理后得 21 ()d x x p x x ⎰)](4)()[(612122202002γβαγβαγβα++++++++-=x x x x x x x x )4(621002y y y x x ++-=)4(6210y y y nab ++-= 同样也有422 ()d x x p x x ⎰)4(6432y y y n ab ++-=……222 ()d n n x nx p x x -⎰)4(621222n n n y y y nab ++-=-- 将这n 个积分相加即得原来所要计算的定积分的近似值:22222212 11()d ()d (4)6ii nnbx i i i i ax i i b af x x p x x y y y n---==-≈=++∑∑⎰⎰, 即021******* ()d [4()2()]6bn n n ab af x x y y y y y y y y n---≈++++++++⎰L L 这就是抛物线法公式,也称为辛卜生(Simpson )公式.仍用 12 0d 1x x +⎰的近似计算为例,取100=n ,102132124222 0d [4()2()]16n n n x b ay y y y y y y y x n ---≈+++++++++⎰L L=0.78539816339745,理论值 12 0d 14x x π=+⎰,此时计算的相对误差 160.7853981633974542.827104ππ--=≈⨯4. 直接应用Matlab 命令计算结果(1) 数值计算 120d .1xx +⎰ 方法1:int('1/(1+x^2)','x',0,1) (符号求积分)方法2:quad('1./(1+x.^2)',0,1) (抛物线法求数值积分)方法3:x=0:0.001:1; y=1./(1+x.^2);trapz(x,y) (梯形法求数值积分) (2)数值计算 212 01d d x x y y -+⎰⎰方法1:int(int('x+y^2','y',-1,1),'x',0,2) (符号求积分)方法2:dblquad(inline('x+y^2'),0,2,-1,1) (抛物线法二重数值积分)四、自己动手1. 实现实验内容中的例子,即分别采用矩形法、梯形法、抛物线法计算 120d 1xx +⎰,取258=n ,并比较三种方法的精确程度.2. 分别用梯形法与抛物线法,计算 2 1d xx⎰,取120=n .并尝试直接使用函数trapz()、quad()进行计算求解,比较结果的差异.3. 试计算定积分 0sin d xx x+∞⎰.(注意:可以运用trapz()、quad()或附录程序求解吗?为什么?)4. 将 120d 1xx +⎰的近似计算结果与Matlab 中各命令的计算结果相比较,试猜测Matlab 中的数值积分命令最可能采用了哪一种近似计算方法?并找出其他例子支持你的观点.5. 通过整个实验内容及练习,你能否作出一些理论上的小结,即针对什么类型的函数(具有某种单调特性或凹凸特性),用某种近似计算方法所得结果更接近于实际值?6. 学习fulu2sum.m 的程序设计方法,尝试用函数 sum 改写附录1和附录3的程序,避免for 循环.五、附录附录1:矩形法(左点法、右点法、中点法)(fulu1.m ) format long n=100;a=0;b=1;inum1=0;inum2=0;inum3=0; syms x fx fx=1/(1+x^2); for i=1:nxj=a+(i-1)*(b-a)/n; %左点 xi=a+i*(b-a)/n; %右点 fxj=subs(fx,'x',xj); %左点值fxi=subs(fx,'x',xi); %右点值fxij=subs(fx,'x',(xi+xj)/2); %中点值inum1=inum1+fxj*(b-a)/n;inum2=inum2+fxi*(b-a)/n;inum3=inum3+fxij*(b-a)/n;endinum1inum2inum3integrate=int(fx,0,1)integrate=double(integrate)fprintf('The relative error between inum1 and real-value is about: %d\n\n',...abs((inum1-integrate)/integrate))fprintf('The relative error between inum2 and real-value is about: %d\n\n',...abs((inum2-integrate)/integrate))fprintf('The relative error between inum3 and real-value is about: %d\n\n',...abs((inum3-integrate)/integrate))附录2:梯形法(fulu2.m)format longn=100;a=0;b=1;inum=0;syms x fxfx=1/(1+x^2);for i=1:nxj=a+(i-1)*(b-a)/n;xi=a+i*(b-a)/n;fxj=subs(fx,'x',xj);fxi=subs(fx,'x',xi);inum=inum+(fxj+fxi)*(b-a)/(2*n);endinumintegrate=int(fx,0,1)integrate=double(integrate)fprintf('The relative error between inum and real-value is about: %d\n\n',...abs((inum-integrate)/integrate))附录2sum:梯形法(fulu2sum.m),利用求和函数,避免for 循环format longn=100;a=0;b=1;syms x fxfx=1/(1+x^2);i=1:n;xj=a+(i-1)*(b-a)/n; %所有左点的数组xi=a+i*(b-a)/n; %所有右点的数组fxj=subs(fx,'x',xj); %所有左点值fxi=subs(fx,'x',xi); %所有右点值f=(fxi+fxj)/2*(b-a)/n; %梯形面积inum=sum(f) %加和梯形面积求解integrate=int(fx,0,1)integrate=double(integrate)fprintf('The relative error between inum and real-value is about: %d\n\n',...abs((inum-integrate)/integrate))附录3:抛物线法(fulu3.m)format longn=100;a=0;b=1;inum=0;syms x fxfx=1/(1+x^2);for i=1:nxj=a+(i-1)*(b-a)/n; %左点xi=a+i*(b-a)/n; %右点xk=(xi+xj)/2; %中点fxj=subs(fx,'x',xj);fxi=subs(fx,'x',xi);fxk=subs(fx,'x',xk);inum=inum+(fxj+4*fxk+fxi)*(b-a)/(6*n);endinumintegrate=int(fx,0,1)integrate=double(integrate)fprintf('The relative error between inum and real-value is about: %d\n\n',...abs((inum-integrate)/integrate))。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
fxij=(fxi+fxj)/2;
m=fxj*(b-a)/n;
p=fxi*(b-a)/n;k=fxij*(b-a)/n;
inum1=sum(m)
inum2=sum(p)
inum3=sum(k)
改写抛物线法
ex6gaixiepaowuxianfa.m
format long
syms x fx
fx=1./x;
for k=1:n
i=2*k-1; j=2*k;
xi=a+(b-a)*i/(2*n);
xj=a+(b-a)*j/(2*n);
fxi=subs(fx,'x',xi);
fxj=subs(fx,'x',xj);
m=m+fxi;p=p+fxj;
end
fxa=subs(fx,'x',a);
2,int使用较广泛,对被积函数的要求较弱,对广义积分同样适用。
3,Matlab编程中出错时可以根据提示修改,当不会修改时,换一种方法也行。可以尝试
多种方法求解。
4,sum()可以避免for循环
思考与深入:
Matlab的计算功能强大,可以方便的进行一些很复杂的计算,而且运算效率极高。
我对Matlab了解不够透彻,有时不知自己哪里出错,根据提示也不能修改,我会在今后的
学习应用中不断了解。
教师评语:
abs(inum-integrate)/integrate)
改写矩形法
ex6gaixiejuxingfa.m
format long
n=100;a=0;b=1;
syms x fx
fx=1/(1+x^2);
i=1:n;
xj=a+(i-1)*(b-a)/n;
xi=a+i*(b-a)/n;
fxj=subs(fx,'x',xj);
trapz(x,y)
(2)抛物线法:
formatlong
n=120;a=1;b=2;inum=0;
symsxfx
fx=1/x;
fori=1:n
xj=a+(i-1)*(b-a)/n;
xi=a+i*(b-a)/n;
xk=(xi+xj)/2;
fxj=subs(fx,'x',xj);
fxi=subs(fx,'x',xi);
针对不同 的取法,计算结果会有不同。
(1) 左点法:对等分区间
,在区间 上取左端点,即取 。
(2)右点法:同(1)中划分区间,在区间 上取右端点,即取 。
(3)中点法:同(1)中划分区间,在区间 上取中点,即取 。
2.梯形法
等分区间
, 相应函数值为 ( ).曲线 上相应的点为 ( )将曲线的每一段弧 用过点 , 的弦 (线性函数)
ex2zhi2paowuxianfa.m
format long
n=120;
a=1;b=2;m=0;p=0;
syms x fx
fx=1./x;
for k=1:n
i=2*k-1; j=2*k;
xi=a+(b-a)*i/(2*n);
xj=a+(b-a)*j/(2*n);
fxi=subs(fx,'x',xi);
fxb=subs(fx,'x',b);
inum=(b-a)*(fxa+fxb+4*m+2*p)/(6*n)
integrate=int(fx,1,2)
integrate=double(integrate)
fprintf('The relative error between inum and real-value is about:%g\n\n',...
fxk=subs(fx,'x',xk);
inum=inum+(fxj+4*fxk+fxi)*(b-a)/(6*n);
end
inum
quad('1./x',1,2)
2、(1)符号求积分:
int('sin(x)/x','x',0,inf)
(2)quad('sin(x)./x',0,inf)
3ex2zhi2tixingfa.m
来代替,这使得每个 上的曲边梯形成为真正的梯形,其面积为
, .于是各个小梯形面积之和就是曲边梯形面积的近似值,
,
即 ,称此式为梯形公式。
3.抛物线法
将积分区间 作 等分,分点依次为
, ,对应函数值为
( ),曲线上相应点为
( ).
现把区间 上的曲线段 用通过三点 , , 的抛物线
来近似代替,然后求函数 从 到 的定积分:
run('E:\matlab\ex2zhi2trapz.m')
inum =
0.69315152080005
integrate =
log(2)
integrate =
0.69314718055995
The relative error between inum and real-value is about:6.26164e-006
应用近似方法去计算相应的定积分,以便更好地熟练地掌握定积分的近似计算.
实验目的:
1加深理解积分理论中分割、近似、求和、取极限的思想方法;
2了解定积分近似计算的矩形法、梯形法与抛物线法;
3会用MATLAB语言编写求定积分近似值的程序,会用MALAB中的命令求定积分。
实验内容
1分别用梯形法与抛物线法,计算 ,取n=120.要求使用函数trapz()、quad()进行计算求解,并比较结果的差异;
format long
n=120;a=1;b=2;inum=0;
syms x fx
fx=1./x;
for i=1:n
xj=a+(i-1)*(b-a)/n; xi=a+i*(b-a)/n;
fxj=subs(fx,'x',xj); fxi=subs(fx,'x',xi);
inum=inum+(fxj+fxi)*(b-a)/(2*n);
symsxfx
fx=1/x;
fori=1:n
xj=a+(i-1)*(b-a)/n;
xi=a+i*(b-a)/n;
fxj=subs(fx,'x',xj);
fxi=subs(fx,'x',xi);
inum=inum+(fxj+fxi)*(b-a)/(2*n);
end
inum
x=1:1/120:2;
y=1./x;
由于 ,代入上式整理后得
同样也有
……
将这 个积分相加即得原来所要计算的定积分的近似值:
,
即
这就是抛物线法公式,也称为辛卜生(Simpson)公式.
实验所用软件及版本:2012B
主要内容:
1,分别用梯形法与抛物线法,计算 ,将积分区间[1,2]作120等分。并尝试用函数trapz(),quad()进行算求解,比较结果的差异。
2,试计算定积分.
(注意:可以运用trapz()、quad()、或附录程序求解吗?为什么?)
3,学习fuluBsum.m的程序设计方法,尝试用函数sum改写矩形法和抛物线法的程序,避免for循环。
实验过程记录(含基本步骤、主要程序清单及异常情况记录等):
1、(1)梯形法:
、formatlong
n=120;a=1;b=2;inum=0;
n=100;a=0;b=1;
syms x fx
fx=1/(1+x^2);
i=0:(n-1);
xj=a+(2*i)*(b-a)/(2*n);
xi=a+(2*i+1)*(b-a)/(2*n);
xk=a+(2*i+2)*(b-a)/(2*n);
fxj=subs(fx,'x',xj);
fxi=subs(fx,'x',xi);
数学实验报告
实验序号:2日期:2013年12月5日
班级
2011应数一班
姓名
孙婉婉
学号
1101114143
实验名称
定积分的近似计算
问题背景描述:
牛顿—莱布尼兹公式仅使用于被积函数的原函数能用初等函数表达出来的情形。如果这
点办不到或者不容易办到,就有必要考虑近似计算的方法。在定积分的许多应用问题中,被
积函数甚至没有解析表达式,可能只是一条实验记录曲线,或者是一组离散的采样值,这时
fxj=subs(fx,'x',xj);
m=m+fxi;p=p+fxj;
end
fxa=subs(fx,'x',a);
fxb=subs(fx,'x',b);
inum=(b-a)*(fxa+fxb+4*m+2*p)/(6*n)
integrate=int(fx,1,2)
integrate=double(integrate)
2试计算定积分 .(注意:可以运用trapz( )、quad( )或附录程序求解吗?为什么?);
3学习fuluBsum.m的程序设计方法,尝试用函数sum改写附录C的程序,避免for循环。
实验原理与数学模型:
实验原理与数学模型:
1. 矩形法:根据定积分的定义,每一个积分和都可以看作是定积分的一个近似值,即 在几何意义上,这是用一系列小矩形面积近似小曲边梯形的结果,所以把这个近似计算方法称为矩形法.不过,只有当积分区间被分割得很细时,矩形法才有一定的精确度.