1、实验三生日蛋糕问题【实验目的】1介绍了数值积分方法，加深对积分概念的理解。2通过实例学习用数值积分知识解决面积、体积等实际应用问题。3学习使用MATLAB软件中有关积分计算的命令。【实验内容】一个数学家即将要迎来他90岁的生日，有很多的学生要来为他祝寿，所以要订做一个特大的蛋糕。为了纪念他提出的一项重要成果口腔医学的悬链线模型，他的弟子要求蛋糕店的老板将蛋糕边缘圆盘半径作成下列悬链线函数：2/5，01（单位：）由于蛋糕店从来没有做过这么大的蛋糕，蛋糕店的老板必须要计算一下成本。这主要涉及两个问题的计算：一个是蛋糕的质量，由此可以确定需要多少鸡蛋和面料；另一个是蛋糕表面积（底面除外），由此确定需
2、要多少奶油。【实验准备】1数值积分用数值方法近似地求一个函数在区间（，）上定积分的基本思路，可以归结到定积分的定义：（1）其中，（，），1，2，。从几何意义上说，对于区间，上非负函数，积分值是与直线，及轴所围成的曲边梯形的面积。在应用问题中，常常是利用微元法进行分析，而问题最终的解归结为微分的和（即积分）。多元函数的积分称为多重积分。二重积分定义为 （2）当非负时，积分值几何上表示曲顶柱体的体积，二重积分的计算主要是转换为两次单积分来解决，无论是解析方法还是数值方法，如何实现这种转换是解决问题的关键。2定积分计算公式由定积分的几何意义我们知道，定积分表示曲线下的面积，不妨先从图形上看看如何近似
3、计算这块面积。先将区间（，）等分，称为积分步长，其中，在小区间上用矩形面积近似下面的梯形的面积，可取左端点函数值为小矩形的高（图中的虚线），或取右端点的函数值为小矩形的高（图中的虚线），于是在整个区间（，）内构成台阶形，容易看出，两个台阶的面积分别为，（3），（4）在图3.1中，两个台阶形分别小于和大于所求面积，（3）、（4）就是计算定积分的矩形公式。 图3.1定积分的矩形公式和梯形公式示意图若将二者平均，则每个小区间上的小矩形变为小梯形（如图中粗实线），整个区间上的梯形面积为，（5）它相当于用分段线性插值函数作为的近似，（5）式就是常用的梯形求积分公式。为提高精度可以用分段二次插值函数代替，
4、由于每段要用到相邻两个小区间端点的三个函数值，所以小区间的数目必须是偶数，记，0，1，。在第段的两个小区间上用三个节点（，），（，），（，）作二次插值函数，然后积分可得求段之和就得到整个区间上的近似积分，（6）（6）称为辛普森求积公式（抛物线公式）。3求积分的MATLAB命令z = trapz( y )如果函数y是一个向量，则z返回其积分值；若是矩阵，则z返回一个行向量，每个分量表示对矩阵函数每列的积分；z = trapz( x , y )x是表示积分区间离散化向量；y是与x对应的同维数的向量，表示被积函数；z返回积分的近似值；trapz是最基本的数值积分方法，采用梯形积分法，对积分区间n等分
5、，固定步长，精度低，适用于数值函数和光滑性不好的函数。例1，求积分在MATLAB命令窗口输入下述命令： x = 0:pi/100:pi; y = sin(x); a= trapz(x,y)a = 1.9998 b= pi/100*trapz(y)%步长乘以函数y的积分值上面分别使用了trapz命令的两种形式，可以得到正弦函数在区间0，上相同的积分值。对于变步长的数值积分，可用quadl、quad或quad8（MATLAB5.3以上版本中被quadl替代），其使用形式如下：q = quadl( fun , a , b )给出函数fun在区间a，b内数值积分的近似值；q = quadl( fun
6、, a , b , tol )tol为容许误差精度，默认为1e6；q = quadl( fun , a , b , tol , trace )在有轭迹的情况下，显示循环时函数计算次数、a、ba和积分值q； q , fcnt = quadl( fun , a , b , . )返回函数计算次数给fcnt；fun可以是内在字符串对象或M函数文件名，定义函数fun时使用.*，./或.运算符，以便它能对向量或数组进行计算；quadl使用自适应步长Gauss/Lobatto法，quad使用自适应步长Simpson法，建议使用quadl，它不但精度较高，且对假收敛和假奇异积分具有一定的适应性，而quad较
7、差；trapz，quadl，quad都不能用于求广义积分；有关quadl和quad的详细用法可查阅帮助。例2，求积分为解例2，先写M函数文件myfun.mfunction y = myfun(x) y = 1./(x.3-2*x-5);然后在命令窗口输入命令： q = quadl(myfun,0,2)q = -0.4605也可以直接在命令窗口输入下述命令： F=inline(1./(x.3-2*x-5) q=quadl(F,0,2)可得到相同的结果。如求二重积分：则使用dblquad，其使用形式如下：q = dblquad( fun , xmin , xmax , ymin , ymax )f
8、un的定义与quadl中fun定义一样，其余分别是变量x、y的上下限。如求三重积分：q = dblquad( fun , xmin , xmax , ymin , ymax , zmin, zmax )fun的定义与quadl中fun定义一样，其余分别是变量x、y和z的上下限。若求不定积分，可使用下面的命令：R = int( S )符号表达式S的不定积分；R = int( S , v )符号表达式S关于变量v的不定积分；R = int( S , a , b )符号表达式S的定积分，a，b分别为下、上限；R = int( S , v , a , b ) 符号表达式S关于变量v从a到b的定积分；当
9、int求不出不定积分时，会自动转求数值积分；int也可用来求重积分和广义积分，但做数值计算时速度慢，效率不好；有关int的详细用法可查阅MATLAB帮助。 syms x%指明x为符号对象 R=int(-2*x/(1+x2)2)R =1/(1+x2) a=int(x*log(1+x),0,1)%求表达式从0到1的定积分a =1/4 t=int(exp(-xsin(x),0,1)Warning: Explicit integral could not be found. In D:MATLAB6p5toolboxsymbolicsymint.m at line 58t =int(exp(-xsin
10、(x),x = 0 . 1)%该段文字说明无法求得解析解 t=vpa(t,4)%以4位小数的精度给出符号表达式的值t =.4549%得到近似解例3，计算广义积分 z=int(exp(-x)*sin(x),1,inf)%其中inf表示无穷大z =1/2*exp(-1)*cos(1)+1/2*exp(-1)*sin(1)【实验方法与步骤】1引例问题的分析对于一个圆盘形的单层蛋糕，如图2（a）绕水平中心轴旋转而成，若高为（），半径为（），密度为（），则蛋糕的质量（）和表面积（）为：2如果蛋糕是双层圆盘的，如图2（b）绕水平中心轴线旋转而成，每层高/2，下层蛋糕半径为，上层蛋糕半径为，此时蛋糕的质量和
11、表面积为：/2/2（）/22/22/2（） /2 /2 （a） （b） （c） 图3.2不同层数的生日蛋糕依次类推，如果蛋糕是层的，每层高为/，半径分别为，则蛋糕的质量和表面积为2事实上，蛋糕边缘圆盘半径2/5，01那么当，1时22数学家的生日蛋糕问题转化为求上面两个数值积分2MATLAB命令求解蛋糕边缘级圆盘半径是高度的函数2/5，01现要求积分和2在MATLAB命令窗口输入下述命令： syms h%定义h为自变量 r=2-(exp(2*h)+exp(-2*h)/5; quadl(pi*(2-(exp(2*h)+exp(-2*h)/5).2,0,1)%求积分ans =5.4171 r0=su
12、bs(r,h,0)%计算的值r0 =1.6000 quadl(2*pi*(2-(exp(2*h)+exp(-2*h)/5),0,1)+pi*r02ans =16.0512 %求得2的值求得该数学家的生日大蛋糕的质量和表面积为5.4171（），16.0512（）实际上，问题的分析部分使用的是数值积分中的矩形法，求解时直接使用MATLAB积分命令求解，使用不定积分int也可以得到相同的解答结果。【练习与思考】1人造地球卫星可视为平面上的椭圆。我国第一颗人造地球卫星近地点距地球表面439，远地点距地球表面2384，地球半径为6371，求该卫星的轨道长度。2炮弹射击的目标为一正椭圆开区域，当瞄准目标的
13、中心发射时，在众多因素的影响下，炮弹着地点与目标中心有随机偏差。可以合理地假设炮弹着地点围绕中心呈二维正态分布，且偏差在方向和方向相互独立。若椭圆形区域在方向半轴长120，方向半轴长80，设炮弹着地点偏差的均方差在方向和方向均为100，求炮弹落在椭圆形区域内的概率。3某个阿拉伯国家有一座著名的伊斯兰教堂，它以中央大厅的金色巨大拱形圆顶名震遐迩。因年久失修，国王下令将教堂顶部重新贴金箔装饰。据档案记载，大厅的顶部形状为半球面，其半径为30。考虑到可能的损耗和其他技术因素，实际用量将会比教堂顶部面积多1.5%。据此，国王的财政大臣拨出了可制造5800有规定厚度金箔的黄金。建筑商人哈桑略通数学，他计算了一下，觉得黄金会有盈余。于是，他以较低的承包价得到了这项装饰工程。但在施工前的测量中，工程师发现教堂顶部实际上并非是一个精确的半球面，而是一个半椭球碳，其半立轴恰是30，而半长轴和半短轴分别是30.6和29.6。这一来哈桑犯了愁，他担心黄金是否还有盈余？甚至可能短缺，最后的结果究竟如何呢？157
