实验2函数的极限、导数、积分与M文件的应用 实验目的 1、理解极限的概念 2、掌握用 matlab软件求函数极限的方法; 3、加深理解积分理论中分割、近似、求和、取极限的思想方法 4、学习并掌握用 matlab求导数,不定积分、定积分、二重积分、曲线积分的方法。 实验内容: 1、学习 matlab命令 应用 matlab求极限的命令: 对应数学运算 matlab命令 lim f(x) limit(f imf(x)|m2减或者lm lim f(x) x→日 limit(f,x,a,'let”) lim f(x) limit(, x, a, right) x→a+ 应用 matlab求解方程的函数 solve调用格式: solve(f(x)):表示求方程f(x)=0的根 2、理解极限概念:参照课本 3、M文件的创建与运行: MATLA B是一个功能强大的程序设计语言,有着良好的交互计算的工作环境。包 含 MATLAB语言代码的文件称为M一文件,你可以使用文本编辑器来建立M一文件, 然后可以象使用其它的 MA TLAB函数或命令一样使用它们。 (1)命令M文件及其运行 命令文件没有输入参数,也不返回输出参数,只是一些命令行的组合。命令M文件中 的命令可以访问 MATLAB工作区中的所有变量,而且其中的所有变量也成为工作区的 一部分,命令文件运行结束,命令文件产生的变量保留在工作区,直到关闭 MATLAB 或用命令删除 (2)函数M文件及其调用: 在 MATLAB编辑窗口还可建立函数M文件,我们可以根据需要建立自己的函数文件, 它们能够像库函数一样方便地调用,从而可扩展 MATLAB的功能。如果对于一类特殊 的问题,建立起许多函数M文件,就能形成工具箱。格式如下 function[输出参数列表]=函数名(形式参数列表) 函数体语句 end 4、x=sym(x’)建立符号变量x syms x y Z建立多个符号变量x,y,z diff(f(x)) 求f(x)的一阶导数 diff(f(x), n) 求f(x)的n阶导数 diff(f(x, y), x) 求二元函数f(x,y)的关于x的偏导数
实验 2 函数的极限、导数、积分与 M 文件的应用 实验目的： 1、理解极限的概念； 2、掌握用 matlab 软件求函数极限的方法； 3、加深理解积分理论中分割、近似、求和、取极限的思想方法； 4、学习并掌握用 matlab 求导数，不定积分、定积分、二重积分、曲线积分的方法。 实验内容： 1、学习 matlab 命令： 应用 matlab 求极限的命令： 对应数学运算 matlab 命令 lim ( ) 0 f x x→ limit(f) lim f (x) x→a limit(f,x,a)或者 limit(f,a) lim f (x) x→a− limit(f,x,a,’left’) lim f (x) x→a+ limit(f,x,a,’right’) 应用 matlab 求解方程的函数 solve 调用格式： solve(f(x)) ：表示求方程 f(x)=0 的根。 2、 理解极限概念： 参照课本 3、 M 文件的创建与运行： MATLAB 是一个功能强大的程序设计语言，有着良好的交互计算的工作环境。包 含 MATLAB 语言代码的文件称为 M－文件，你可以使用文本编辑器来建立 M－文件， 然后可以象使用其它的 MA TLAB 函数或命令一样使用它们。 （1） 命令 M 文件及其运行： 命令文件没有输入参数，也不返回输出参数，只是一些命令行的组合。命令 M 文件中 的命令可以访问 MATLAB 工作区中的所有变量，而且其中的所有变量也成为工作区的 一部分，命令文件运行结束，命令文件产生的变量保留在工作区，直到关闭 MATLAB 或用命令删除。 （2）函数 M 文件及其调用： 在 MATLAB 编辑窗口还可建立函数 M 文件，我们可以根据需要建立自己的函数文件， 它们能够像库函数一样方便地调用，从而可扩展 MATLAB 的功能。如果对于一类特殊 的问题，建立起许多函数 M 文件，就能形成工具箱。格式如下： function [输出参数列表]=函数名（形式参数列表） 函数体语句； end 4、x=sym('x') 建立符号变量 x syms x y z 建立多个符号变量 x,y,z diff(f(x)) 求 f(x)的一阶导数 diff(f(x),n) 求 f(x)的 n 阶导数 diff(f(x,y),x) 求二元函数 f(x,y)的关于 x 的偏导数
diff(f(x, y), x, n) 求二元函数f(x,y)的关于x的n阶偏导数 x为向量时给出x的各个元素之和,若x为矩阵时给出x的列和 symsum(s, n) 表示求∑S 表示求∑S t(f(x)) 计算不定积分f(x)x int(r(,y,x)计算不定积分f(x,y)dtx int(f(x), a, b) 计算定积分」f(x)dhx int(f(,y),x,a,b)计算定积分」f(x,y)dt 例1观察数列 当n→>∞时的变化趋势。 解:输入 matlab命令: n=1:100 xn=n./(n+1); for主=1:100; plot(n(i), xn(i),'b') 其中 for end是循环语句,循环体内的语句被执行100次,n() 表示n的第i个分量。由图可以看出,随着n的增大,点列与直线y= 无限接近,于是有结论:lm、n n→∞n+1 例2分析函数f(x)=xSn一当x→0时的变化趋势。 解:输入 matlab命令 x=-1:0.01:1 plot(x,y)
diff(f(x,y),x,n) 求二元函数 f(x,y)的关于 x 的 n 阶偏导数 sum(x) x 为向量时给出 x 的各个元素之和，若 x 为矩阵时给出 x 的列和 symsum(s,n) 表示求 n s 1 symsum(s,k,m,n) 表示求 = n k m s int(f(x)) 计算不定积分 f (x)dx int(f(x,y),x) 计算不定积分 f (x, y)dx int(f(x),a,b) 计算定积分 b a f (x)dx int(f(x,y),x,a,b) 计算定积分 b a f (x, y)dx 例1 观察数列 n+! n 当 n → 时的变化趋势。 解：输入 matlab 命令： n=1:100; xn=n./(n+1); for i=1:100; plot(n(i),xn(i),'b') hold on end 其中 for end是循环语句，循环体内的语句被执行100次，n(i) 表示 n 的第 i 个分量。由图可以看出，随着 n 的增大，点列与直线 y = 1 无限接近，于是有结论： 1 1 lim = → n + n n 。 例2 分析函数 x f x x 1 ( ) = sin 当 x →0 时的变化趋势。 解：输入 matlab 命令 x=-1:0.01:1; y=x.*sin(1./x); plot(x,y)
hold on plot(x, x,x, -2 从图上看,fx)随着的减小,振幅越来越接近零,频率越来越高,做无限振荡 例3分析函数f(x)=Sn一当x→>O时的变化趋势 解:输入 matlab命令 x=-1:0.01:1; y=sin(1/x)i plot(x,y) 从图上看,当x→0时,Sn一在-1和1之间无限次振荡,极限不存在。仔细观察图 形,回答为什么图像的有些峰值不是-1和1? 例4分析函数∫(x)=(1+-)当x→∞时的变化趋 解:输入 matlab命令 =1:100; y=(1+1./x).^x plot(x,y) 从图上可以看出当x→时,函数值与某常数无限接近,我们知道这个常数就是e。 例5求极限tanx-snx 解:输入 matlab命令 limit((tan(x)-sin(x))/(sin(x)3) 例6解方程ax2+bx+c=0 解:输入 matlab命令 syms a b cx f=a*x 2+b*x+c; solve(f) 例7解方程x3-5x-1=0 解:输入 malab命令
hold on plot(x,x,x,-x) 从图上看，f(x)随着|x|的减小，振幅越来越接近零，频率越来越高，做无限振荡。 例3 分析函数 x f x 1 ( ) = sin 当 x →0 时的变化趋势。 解：输入 matlab 命令 x=-1:0.01:1; y=sin(1./x); plot(x,y) 从图上看，当 x →0 时， x 1 sin 在-1和1之间无限次振荡，极限不存在。仔细观察图 形，回答为什么图像的有些峰值不是-1和1 ？ 例4 分析函数 x x f x ) 1 ( ) = (1+ 当 x → 时的变化趋势。 解：输入matlab命令 x=1:100; y=(1+1./x).^x; plot(x,y) 从图上可以看出当 x → 时，函数值与某常数无限接近，我们知道这个常数就是 e 。 例5 求极限 x x x x 3 0 sin tan sin lim − → 。 解：输入matlab命令 syms x limit((tan(x)-sin(x))/(sin(x)^3)) 例6 解方程 0 2 ax + bx + c = 。 解：输入matlab命令 syms a b c x f=a*x^2+b*x+c; solve(f) 例7 解方程 5 1 0 5 x − x − = 。 解：输入malab命令 syms x
f=x^5-5*x-1; solve(f) 例8分别建立命令式M文件和函数式M文件求向量的平均值和标准差 解:命令式函数文件 在文件编辑窗口输入如下命令: x=rand (20)i stdev=sqrt(sum(x. 2)/n-mean. 2) 随便以一个文件名命名即可,例如保存为 programme. m,在命令窗口输入 programme 即可得到结果。 函数式函数文件 在文件编辑窗口输入如下命令: function [mean, stdev]=s tat(x) mean=sum(x)/n stdevesgrt(sum(x. 2)/n-mean. 2) 保存为文件 stat m,在命令窗口输入 x=rand(20)i [me an, stdev]=stat(x) 即可得到结 例9画出f(x)=e2在x=0(P(0,1))处的切线及若干条割线,观察割线的变化趋势。 解:在曲线y=e2上另取一点M(h,e"),则PM的方程是 x+1,取h=3,2,10.1,0.01,分别作出几条割线。 h 然后作出y=e在x=0处的切线y=x+1。 h=[3,2,1,0.1,0.01]; a=(exp(h)-1)./h 1:0.1:3; plot(x, exp(x),'r) hold on for i=l: 5 plot(h(i), exp(h(i)),'r') plot(x, a(i)*x+l axis square plot(x, x+l,'r') 从图上看出割线越来越接近切线
f=x^5-5*x-1; solve(f) 例8 分别建立命令式M文件和函数式M文件求向量的平均值和标准差。 解：命令式函数文件： 在文件编辑窗口输入如下命令： x=rand(20); n=length(x); mean=sum(x)/n stdev=sqrt(sum(x.^2)/n-mean.^2) 随便以一个文件名命名即可，例如保存为programme.m，在命令窗口输入programme 即可得到结果。 函数式函数文件： 在文件编辑窗口输入如下命令： function [mean,stdev]=stat(x) n=length(x); mean=sum(x)/n stdev=sqrt(sum(x.^2)/n-mean.^2) 保存为文件stat.m，在命令窗口输入： x=rand(20); [mean,stdev]=stat(x) 即可得到结果。 例 9 画出 x f (x) = e 在 x = 0 （ P(0,1) ）处的切线及若干条割线，观察割线的变化趋势。 解：在曲线 x y = e 上另取一点 ( , ) h M h e ，则 PM 的方程是： 1 1 + − = x h e y h ，取 h = 3,2,1,0.1,0.01 ，分别作出几条割线。 然后作出 x y = e 在 x = 0 处的切线 y = x +1。 h=[3,2,1,0.1,0.01]; a=(exp(h)-1)./h; x=-1:0.1:3; plot(x,exp(x),'r'); hold on for i=1:5 plot(h(i),exp(h(i)),'r.') plot(x,a(i)*x+1) end axis square plot(x,x+1,'r') 从图上看出割线越来越接近切线
例10求隐函数所确定函数的导数hnx+ex=e 解:F(x,y)=nx+ex-e,先求F'(x),再求F(y), dy F'(x) dxF(y)’于是输入 matlab命令 diffx=diff(log(x)+exp(-y/x)-exp(1),x); diffy=diff(log(x)+exp(-y/x)-exp(1),y) ￡fx/ diffy 例11设隐函数n(xy)+cOS(y=)+tan(x=)=0,求 解:利用 matlab命令中的 jacobian函数简化步骤, sys x y Z a=jacobian(sin(x*y)+cos (y*z)+tan (z*x), Lx, y, z]: dz dx=-a(1)/a 3 dz_dy=a(2)/a(3) 例2利用定积分定义计家∫edt 解:等距划分为20个子区间 左矩形法求和近似 x=linspace 0, 1, 21) y=exp(x) y1=y(1:20) s1=sum(y1)/20 求得s1=1.6757,图形演示如下 plot(x, y) for i=1: 20 fi1(x(i),x(i+1),x(i+1),x(i),x(i)],[0,0,y(i),y(i),0],"b d 右矩形法求和近似 x=linspace 0, 1, 21) y2=y(2:21); 求得s2=1.7616,图形演示如下 plot(x, y) hold for i=1: 20
例 10 求隐函数所确定函数的导数 x e e x y + = − ln 。 解： F x y x e e x y = + − − ( , ) ln ，先求 F'(x) ，再求 F'( y) ， 则 '( ) '( ) F y F x dx dy = − ，于是输入 matlab 命令 syms x y diffx=diff(log(x)+exp(-y/x)-exp(1),x); diffy=diff(log(x)+exp(-y/x)-exp(1),y); -diffx/diffy 例 11 设隐函数 sin( x y) + cos( yz) + tan( x z) = 0 ，求 y z x z , 。 解：利用 matlab 命令中的 jacobian 函数简化步骤， syms x y z a=jacobian(sin(x*y)+cos(y*z)+tan(z*x),[x,y,z]); dz_dx=-a(1)/a(3) dz_dy=-a(2)/a(3) 例 12 利用定积分定义计算 1 0 e dx x 。 解：等距划分为 20 个子区间 左矩形法求和近似： x=linspace(0,1,21); y=exp(x); y1=y(1:20); s1=sum(y1)/20 求得 s1 =1.6757，图形演示如下 plot(x,y) hold on for i=1:20 fill([x(i),x(i+1),x(i+1),x(i),x(i)],[0,0,y(i),y(i),0],'b') end 右矩形法求和近似： x=linspace(0,1,21); y=exp(x); y2=y(2:21); s2=sum(y2)/20 求得 s2 =1.7616，图形演示如下 plot(x,y) hold on for i=1:20
fi11([x(i),x(i+1),x(i+1),x(i),x(i)],[0,0,y(i+1),y(i+1),0],b) s1和s2都可以作为该定积分的近似值 直接用积分命令int也可以求,命令如下 sys x int(exp(x),0, 1) 例13判断广义积分」d dx的敛散性,收敛时计 算积分值 解:对第一个积分: syms int(1/xp, x, 1, inf) 由结果可以看出,当p1时, aS=1/(-p+1) 对第二个积分: int(1/(2*i)^(1/2)*exp(-x^2/2),-inf,inf) 对第三个积分: int(1/(1-x)2,0,2) 结果为aHS=nf,所以广义积分不收敛。 练习 1、生物学中,有一个刻画生物群体的个体总量增长情况的著名方程一— Logistic方程: Pn+=kP(1-Pa),其中P为某一生物群体的第n代的个体总量与该群体所能达到的最大个体 总量之比,0≤Pn≤1(m=0.1,2,…)k为比例系数 给定初值Po和比例系数k的值后,有方程就能生成一个数列,Po,P1,P2 生物学家为了预测群体总量的变化情况,就要研究这个数列,他们感兴趣的是:这个数列存 在极限吗?数列中的项会出现周期性的变化情况吗?数列会不会出现无法预测的紊乱情 况?(试给定不同的初值Po和比例系数k,研究生物群体变化情况) 2、完全数是指它的所有因子之和等于该完全数。6是一个完全数,因为6的因子为1,2,3 且1+2+3=6.试分别建立命令式M文件和函数式M文件,找出10000以内的所有完全数, 并对它做素因子分解。你能据此猜测完全数的通式吗? 3、绝对素数问题:一个素数,如果其中的任意两位数字交换后得到的数仍为素数,则将此 素数成为绝对素数。例如:13和31都是素数,则13称为绝对素数。试利用 matlab编程求 出三位数中的所有绝对素数。要求分别使用命令式M文件和函数试M文件编写。 4、卫星轨道长度问题和覆盖面积问题
fill([x(i),x(i+1),x(i+1),x(i),x(i)],[0,0,y(i+1),y(i+1),0],'b') end s1 和 s2 都可以作为该定积分的近似值。 直接用积分命令 int 也可以求，命令如下 syms x int(exp(x),0,1) 例 13 判断广义积分 + 1 1 dx x p ， e dx x + − − 2 2 2 1 ， − 2 0 2 (1 ) 1 dx x 的敛散性，收敛时计 算积分值。 解：对第一个积分： syms p int(1/x^p,x,1,inf) 由 结 果 可 以 看 出 ， 当 p 1 时 ， x^ (− p +1) 为 无 穷 ， 当 p 1 时 ， ans = 1/(− p +1) 。 对第二个积分： int(1/(2*pi)^(1/2)*exp(-x^2/2),-inf,inf) 对第三个积分： int(1/(1-x)^2,0,2) 结果为 ans = inf ，所以广义积分不收敛。 练习： 1、生物学中，有一个刻画生物群体的个体总量增长情况的著名方程——Logistic 方程： Pn+1=kPn(1-Pn)，其中 Pn 为某一生物群体的第 n 代的个体总量与该群体所能达到的最大个体 总量之比， 0 P 1(n = 0,1,2, ) n ,k 为比例系数。 给定初值 P0 和比例系数 k 的值后，有方程就能生成一个数列，P0，P1，P2，……，Pn，…… 生物学家为了预测群体总量的变化情况，就要研究这个数列，他们感兴趣的是：这个数列存 在极限吗？数列中的项会出现周期性的变化情况吗？数列会不会出现无法预测的紊乱情 况？（试给定不同的初值 P0 和比例系数 k，研究生物群体变化情况） 2、完全数是指它的所有因子之和等于该完全数。6 是一个完全数，因为 6 的因子为 1，2，3 且 1+2+3=6. 试分别建立命令式 M 文件和函数式 M 文件，找出 10000 以内的所有完全数， 并对它做素因子分解。你能据此猜测完全数的通式吗？ 3、绝对素数问题：一个素数，如果其中的任意两位数字交换后得到的数仍为素数，则将此 素数成为绝对素数。例如：13 和 31 都是素数，则 13 称为绝对素数。试利用 matlab 编程求 出三位数中的所有绝对素数。要求分别使用命令式 M 文件和函数试 M 文件编写。 4、 卫星轨道长度问题和覆盖面积问题：
(1)人造地球卫星轨道可视为平面上的椭圆,我国第一颗人造地球卫星近地点距地球表 面439km,远地点距地球表面2384km,地球半径为6371km,问该卫星的轨道长度 (2)一颗地球同步轨道通讯卫星的轨道位于地球的赤道平面内,且可近似认为是圆轨道, 通讯卫星运行的角速度与地球自转的角速度相同,所以人们看到它在天空不动。问 卫星距离地面的高度h应该为多少?并计算通讯卫星的覆盖面积。 5、路灯照明问题 在一条20m宽的道路两侧,分别按装了一只2kW和3kW的路灯,它们离地面的高度 分别为5m和6m (1)在漆黑的夜晚,当两只路灯开启时,两只路灯连线的路面上最暗的点和最亮的点在 那里? (2)如果3kW的路灯的高度可以在3m到9m之间变化,如何使路面上最暗点的亮度最 大?作数值计算并绘图加以说明。 Psin a (注:1=r2,其中1是照度,r是光源到照射点的距离,a是光源到照射点的光线 与水平面的夹角,P是路灯的功率)
（1） 人造地球卫星轨道可视为平面上的椭圆，我国第一颗人造地球卫星近地点距地球表 面 439km，远地点距地球表面 2384km，地球半径为 6371km，问该卫星的轨道长度？ （2） 一颗地球同步轨道通讯卫星的轨道位于地球的赤道平面内，且可近似认为是圆轨道， 通讯卫星运行的角速度与地球自转的角速度相同，所以人们看到它在天空不动。问 卫星距离地面的高度 h 应该为多少？并计算通讯卫星的覆盖面积。 5、路灯照明问题 在一条 20m 宽的道路两侧，分别按装了一只 2kW 和 3kW 的路灯，它们离地面的高度 分别为 5m 和 6m。 （1） 在漆黑的夜晚，当两只路灯开启时，两只路灯连线的路面上最暗的点和最亮的点在 那里？ （2） 如果 3kW 的路灯的高度可以在 3m 到 9m 之间变化，如何使路面上最暗点的亮度最 大？作数值计算并绘图加以说明。 （注： 2 sin r P I = ，其中 I 是照度，r 是光源到照射点的距离， 是光源到照射点的光线 与水平面的夹角，P 是路灯的功率）
