利用matlab进行微积分的计算.pptx

1、微积分的计算,微积分的符号运算,Matlab定义符号运算和数值运算。如果运算过程中有符号变量，则运算过程为符号运算。 符号变量可以利用两种方法定义。,(1) 利用函数sym定义符号表达式 语句格式 变量=sym(字符串) sym函数可以定义符号变量、符号常数和一般的符号表达式。, x=sym(1/3) x = 1/3 y=x+1/2 %符号运算 y = 5/6,x=1/3 x= 0.3333 y=1/2+1/3 %数值计算 y = 0.8333,例1：比较符号运算和数值计算,例2：分别利用符号运算和数值运算计算如下的和 并比较其计算速度。,注：Matlab语言中利用tic 和toc 计算时间，

2、语言格式为 tic 程序块 toc 显示器显示tic和toc之间的程序块运行时间,符号运算程序： s=sym(0); tic for k=1:1000 s=s+1/k; end toc s,运行结果： elapsed time is 17.471170 seconds. s = 5336291328229478504559104562404298040965247228038426009710134924845626888949710175750609790198503569140908873155046809837844217211788500946430234432656602250210

3、0278425632852081405544941210442510142672770294774712708917963967779610453224692426866468888281582071984897105110796873249319155529397017508931564519976085734473014183284011724412280649074307703736683170055800293659235088589360235285852808160759574737836655413175508131522517/7128865274665093053166384

4、155714272920668358861885893040452001991154324087581111499476444151913871586911717817019575256512980264067621009251465871004305131072686268143200196609974862745937188343705015434452523739745298963145674982128236956232823794011068809262317708861979540791247754558049326475737829923352751796735248042463

5、638051137034331214781746850878453485678021888075373249921995672056932029099390891687487672697950931603520000,数值运算程序 s=0; tic for k=1:1000 s=s+1/k; end toc s 运行结果： Elapsed time is 0.000015 seconds. s = 7.4855,可以看出，数值计算的计算速度远远高于符号运算。 另外，符号运算只能计算非常简单的问题，复杂问题只能利用数值方法。但由于符号运算得到的是解析式，在许多情况下有利于问题的进一步分析。,利用

6、语句 syms 定义符号变量 语句形式 syms x y z syms语句 一次可定义多个符号变量，但不能定义符号常数和表达式,微积分计算问题的符号运算,极限运算 语句格式 G= limit(F); %表达式F中变量趋向于0的极限 G= limit(F,a); %表达式F中变量趋向于a的极限 G= limit(F,v,a); %表达式F中变量v趋向于a的极限 G= limit(F,v,a,right); %表达式F中变量趋向于a的右极限。,例3：计算 syms x f=(1-x)/(1+x)(2/x) f = (-(x - 1)/(x + 1)(2/x) limit(f,0) ans = 1/

7、exp(4),例4：计算 syms x y z f=sin(x*y)/sin(x*z) f = sin(x*y)/sin(x*z) limit(f,x,0) ans = y/z,导数 求导数的语句为 df=diff(F) %求表达式F的一阶导数 df=diff(F,n) %求表达式F的n阶导数 df=diff(F,u,n) %求表达式F的关于变量u的n阶导数,例5：求函数 f=sin(x)/x 的二阶导数。 f=sin(x)/x f = sin(x)/x diff(f,2) ans = (2*sin(x)/x3 - sin(x)/x - (2*cos(x)/x2,例6:求函数 关于x和y的二阶

8、偏导数。 f=(x-y)/(x+y) f = (x - y)/(x + y) dfxy=diff(diff(f,y),x) dfxy = (2*(x - y)/(x + y)3,注：如果表达式F有多个符号变量，且没有指明对哪一个变量求导，则 Matlab规定对离x最近的变量求导。其他运算也有类似规定。,例7：求 f=x2*sin(x) exp(x)*cos(x) log(1+x2) atan(x); diff(f) ans = x2*cos(x) + 2*x*sin(x), exp(x)*cos(x) - exp(x)*sin(x) (2*x)/(x2 + 1), a/(x2 + 1),Mat

9、lab的微积分符号运算都可以对数组进行。,函数的积分 积分符号运算的基本语句 int(F)； %求函数表达式F的不定积分 int(F,v)； %求函数表达式F关于变量v的不定积分 int(F,a,b)； %求函数表达式F在区间a,b上的定积分 int(F,v,a,b); %求函数表达式F在区间a,b上的关于变量v的定积分,例8：求不定积分 syms a b x I=int(exp(a*x)*cos(b*x),x) I = (exp(a*x)*(a*cos(b*x) + b*sin(b*x)/(a2 + b2),例9：求定积分 syms x I=int(sin(x)15,x,0,pi/2) I

10、= 2048/6435,对含参数的运算，符号计算的优势在于可以提供一个用于进一步分析的解析表达式。,级数的符号计算 级数符号计算的基本语句 S=symsum(g)； % 求通项g的无穷和 S=symsum(g,n)； % 求通项g的对整变量n的无穷和 S=symsum(g,a,b)； % 求通项g的下标在a和b之间的项的和 S=symsum(g,n,a,b)； % 求通项g的下标变量n在a和b之间的项的和,例10:求级数的和 syms n symsum(1/n2,1,inf) ans = pi2/6,例11: 求幂级数的和 syms x n f=xn/n f = xn/n S=symsum(f

11、,n,1,inf) S = piecewise(1 1, -log(1 - x),Taylor展式 Taylor展式的基本语句 F=taylor(f); %求函数f 的马克劳林级数的前6项 F=taylor(f,n); %求函数f 的马克劳林级数的前n项 F=taylor(f,n,x0); %求函数f 的在x0展开的Taylor级数的前n项 F=taylor(f,v,n,v0); %求函数f对变量v在v0展开的Taylor级数的前n项,例12: 求函数 的在0点的taylor展式前8项 syms x S=taylor(1-x)/(1+x),8) S = 2*x6 - 2*x7 - 2*x5 +

12、 2*x4 - 2*x3 + 2*x2 - 2*x + 1,微积分的数值计算,符号运算有很大的局限性。第一是计算是纯推理过程，效率非常低。第二是大量的问题无法利用解析的形式给出，例如我们经常使用的概率积分 就无法利用牛顿-莱布尼兹公式推导得到。,相比之下，数值计算在实际计算中有符号运算无法比拟的优点。一是效率高，对公式稍微复杂一些或规模稍微大一些的问题，只能利用数值计算。二是我们经常遇到的计算问题，都可以利用数值计算得到满意的解，而符号运算只能处理其中一小部分问题。,数值计算需要对连续变量离散求解！ 数值计算需要设计数值算法,数值导数,数值计算是在离散节点上的近似计算，通常对连续的区间要首先进

13、行离散。对计算区间a.b上的导数。,首先将区间离散成等距的节点 其中步长 与下标k无关。,一阶和二阶导数的计算由公式 给出, 误差为 O(h2)。,Matlab的数值导数计算可以利用语句diff实现。diff语句对非符号计算表示向前差分，即对 函数diff(x)输出 diff(x,2)输出 diff(diff(x),即 一般, diff(x,k) 的计算可依次类推。,注： diff(x) 的元素个数比x少 1！,一阶导数的计算 一阶导数可以利用 计算。由公式 df的第i个元素与 f(x) 在 xi+1/2点的值最接近，误差是O(h2),例13：利用diff语句计算函数 y=sin(x) 的导数

14、，并计算误差。 参考程序： x=0:pi/100:2*pi; y=sin(x); x1=pi/200:pi/100:2*pi; y1=diff(y)./diff(x); y11=sin(x1); subplot(1,2,1) plot(x1,y1,x1,y11) subplot(1,2,2) y12=cos(x1); plot(x1,y1-y12) legend(df,f) legend(error) 结果图：,二阶导数的计算 语句 输出二阶向前差分，即 可以利用 来计算y(x)的二阶导数。,计算的是n-2个内点上的二阶导数，不包含两个端点上的导数！,例14：计算 y=sin(x)的二阶导数，

15、并与 sin(x)比较。 参考程序： x=0:pi/100:2*pi; y=sin(x); dx=diff(x); d2y=diff(y,2)./dx(2:end).2; subplot(1,2,1) plot(x,y,x(2:end-1),d2y) legend(func,numerd2f) subplot(1,2,2) plot(x(2:end-1),d2y+sin(x(2:end-1) legend(error),数值积分 一元数值积分 Matlab中的数值积分函数quad采用自适应复化Simpson公式，语句格式为 I=quad(fun,a,b) %计算函数fun在a,b上的定积分。

16、I=quad(fun,a,b,tol) %计算函数fun在a,b上精度tol的定积分。,fun的 取法： 过去的取法是一个字符串作为所调用的函数名字。如 quad(sin, 0, pi/2) 在matlab7以上的版本中引入了函数句柄，可以在函数quad中利用指示一个函数句柄。例如 sin %调用以sin为函数名的函数子程序 (x,y)(x.3+y.3) %建立并调用一个函数句柄 (x)myfun(x,a) %调用含有参数的子程序myfun,例15：计算积分 k=1,2,3,4,5 并与精确值比较。,I=zeros(1,5); for k=1:5 g=strcat(sin(x).,num2st

17、r(k); I(k)=quad(g,0,pi/2); end I0=1,pi/4,2/3,3*pi/16,8/15; subplot(1,2,1) plot(I-I0,p) legend(errors of quad) for k=1:5 g=strcat(sin(x).,num2str(k); I(k)=quad(g,0,pi/2,1e-12); end subplot(1,2,2) plot(I-I0,p) legend(errors of quad),积分函数quadv 积分函数quadv类似于函数quad，但它可以对含有向量参数的函数积分，而且可以对矩阵函数积分。,例16：计算积分 k

18、=1,2,3,4,5,这里k 可以看做参数向量!,参考语句： format rat I=quadv(sin(x).1:5,0,pi/2) I = 1 355/452 2/3 1065/1808 8/15,例17：计算矩阵函数 的积分,参考语句： format rat I=quadv(x)exp(-x).*cos(x),(1-x)./sqrt(1+x.2); x.*log(1+x),sin(x)./x,1,2) I = 249/7261 -548/2111 1075/769 689/1045,可用于广义积分的积分函数quadgk 函数quadgk采用Gauss-Kronrod数值积分方法，可以解决无穷区间、无穷间断点等问题。其使用方法与quad类似，但积分限可以取inf,例18：计算 参考语句： I=quadgk