利用matlab进行微积分计算

1、MatlabMatlab定义符号运算和数值运算。如果运算过程中有符号变量，则定义符号运算和数值运算。如果运算过程中有符号变量，则运算过程为符号运算运算过程为符号运算。符号变量可以利用两种方法定义。符号变量可以利用两种方法定义。(1) (1) 利用利用函数函数symsym定义符号表达式定义符号表达式语句格式语句格式 变量变量=sym(=sym(字符串字符串) )symsym函数可以定义符号变量、符号常数和一般的符号表达式函数可以定义符号变量、符号常数和一般的符号表达式。 x=sym(1/3) x=sym(1/3) x = x =1/31/3 y=x+1/2 %y=x+1/2 %符号运算符号运算

2、y = y =5/65/6x=1/3x=1/3x=x= 0.3333 0.3333 y=1/2+1/3 % y=1/2+1/3 %数值计数值计算算y =y = 0.8333 0.8333例1：比较符号运算和数值计算例例2 2：分别利用符号运算和数值运算计算如下的：分别利用符号运算和数值运算计算如下的和和 并比较其计算速度。并比较其计算速度。1000131211S注：注：MatlabMatlab语言中利用语言中利用tic tic 和和toc toc 计算时间，语言格式为计算时间，语言格式为 tic tic 程序块程序块 toc toc显示器显示显示器显示tictic和和toctoc之间的程序块运

3、行时间之间的程序块运行时间符号运算程序：s=sym(0);ticfor k=1:1000 s=s+1/k;endtocs运行结果：elapsed time is 17.471170 seconds. s = 5336291328229478504559104562404298040965247228038426009710134924845626888949710175750609790198503569140908873155046809837844217211788500946430234432656602250210027842563285208140554494121044251014

4、2672770294774712708917963967779610453224692426866468888281582071984897105110796873249319155529397017508931564519976085734473014183284011724412280649074307703736683170055800293659235088589360235285852808160759574737836655413175508131522517/7128865274665093053166384155714272920668358861885893040452001

5、991154324087581111499476444151913871586911717817019575256512980264067621009251465871004305131072686268143200196609974862745937188343705015434452523739745298963145674982128236956232823794011068809262317708861979540791247754558049326475737829923352751796735248042463638051137034331214781746850878453485

6、678021888075373249921995672056932029099390891687487672697950931603520000数值运算程序数值运算程序 s=0;s=0; tic tic for for k=1:1000k=1:1000 s=s+1/k s=s+1/k; ; end end toc toc s s运行结果：运行结果：Elapsed time is 0.000015 secondsElapsed time is 0.000015 seconds. .s s = =7.48557.4855可以看出，数值计算的计算速度远远高于符号运算可以看出，数值计算的计算速度远远高

7、于符号运算。另外，符号运算只能计算非常简单的问题，复杂问题只能利用数值另外，符号运算只能计算非常简单的问题，复杂问题只能利用数值方法。但由于符号运算得到的是解析式，在许多情况下有利于问题方法。但由于符号运算得到的是解析式，在许多情况下有利于问题的进一步分析。的进一步分析。利用语句利用语句 syms syms 定义符号变量定义符号变量语句形式语句形式 syms syms x y x y z zsymssyms语句语句 一次可定义多个符号变量，但不能定义符号常数和表达式一次可定义多个符号变量，但不能定义符号常数和表达式极限运算语句格式 G= limit(F); %表达式F中变量趋向于0的极限 G=

8、 limit(F,a); %表达式F中变量趋向于a的极限 G= limit(F,v,a); %表达式F中变量v趋向于a的极限 G= limit(F,v,a,right); %表达式F中变量趋向于a的右极限。例3：计算 syms xf=(1-x)/(1+x)(2/x) f =(-(x - 1)/(x + 1)(2/x) limit(f,0) ans =1/exp(4)xxx20 x11limit例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/zxzxyxsinsinlim0导数求导数的语句

9、为 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)/xf = sin(x)/x diff(f,2)ans =(2*sin(x)/x3 - sin(x)/x - (2*cos(x)/x2例6:求函数 关于x和y的二阶偏导数。f=(x-y)/(x+y)f = (x - y)/(x + y) dfxy=diff(diff(f,y),x)dfxy =(2*(x - y)/(x + y)3yxyxyxf),(注：如果表达式注：

10、如果表达式F F有多个符号变量，且没有指明对哪一个变量求导，则有多个符号变量，且没有指明对哪一个变量求导，则MatlabMatlab规定对离规定对离x 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)tan()1ln()cos()sin(22xaxxexxdxdxMatlabMatl

11、ab的微积分符号运算都可以对数组进行。的微积分符号运算都可以对数组进行。函数的积分函数的积分积分符号运算的基本语句积分符号运算的基本语句 int(F) int(F)； % %求函数表达式求函数表达式F F的不定积分的不定积分 int(F,v) int(F,v)； % %求函数表达式求函数表达式F F关于变量关于变量v v的不的不定积分定积分 int(F,a,b) int(F,a,b)； % %求函数表达式求函数表达式F F在区间在区间a,ba,b上上的定积分的定积分 int(F,v,a,b); % int(F,v,a,b); %求函数表达式求函数表达式F F在区间在区间a,ba,b上的关于变量

12、上的关于变量v 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)bxdxeaxcos例9：求定积分 syms x I=int(sin(x)15,x,0,pi/2) I =2048/64352/015sindxx对含参数的运算，符号计算的优势在于可以提供一个用于进一对含参数的运算，符号计算的优势在于可以提供一个用于进一步分析的解析表达式。步分析的解析表达式。级数的符号计算级数的符号计算级数符号计算的基本语句级数符号计算的基本语句 S=sym

13、sum(g) S=symsum(g)； % % 求通项求通项g g的无穷和的无穷和 S=symsum(g,n) S=symsum(g,n)； % % 求通项求通项g g的对整变量的对整变量n n的无穷和的无穷和 S=symsum(g,a,b) S=symsum(g,a,b)； % % 求通项求通项g g的下标在的下标在a a和和b b之间的之间的项的和项的和 S=symsum(g,n,a,b) S=symsum(g,n,a,b)； % % 求通项求通项g g的下标变量的下标变量n n在在a a和和b b之间的项的之间的项的和和例10:求级数的和 syms n symsum(1/n2,1,inf

14、) ans = pi2/6 121nn例11: 求幂级数的和 syms x n f=xn/n f = xn/n S=symsum(f,n,1,inf) S = piecewise(1 = x, Inf, abs(x) = 1 and x 1, -log(1 - x)1nnnxTaylorTaylor展式展式TaylorTaylor展式的基本语句展式的基本语句 F=taylor(f); % F=taylor(f); %求函数求函数f f 的马克的马克劳林级数的前劳林级数的前6 6项项 F=taylor(f,n); % F=taylor(f,n); %求函数求函数f f 的马克的马克劳林级数的前劳

15、林级数的前n n项项 F=taylor(f,n,x F=taylor(f,n,x0); %); %求函数求函数f f 的在的在x x0展展开的开的TaylorTaylor级数的前级数的前n n项项 F=taylor(f,v,n,v F=taylor(f,v,n,v0); %); %求函数求函数f f对变量对变量v v在在v v0展展开的开的TaylorTaylor级数的前级数的前n n项项例例12: 12: 求函数求函数 的在的在0点的点的taylortaylor展式展式前前8 8项项 syms x syms x S=taylor(1-x)/(1+x),8) S=taylor(1-x)/(1+

16、x),8) S = S = 2 2* *x6 - 2x6 - 2* *x7 - 2x7 - 2* *x5 + 2x5 + 2* *x4 - 2x4 - 2* *x3 + 2x3 + 2* *x2 - 2x2 - 2* *x + 1x + 1xxxf11)(符号运算有很大的局限性。第一是计算是纯推理过程，效率非常符号运算有很大的局限性。第一是计算是纯推理过程，效率非常低。第二是大量的问题无法利用解析的形式给出，例如我们经常低。第二是大量的问题无法利用解析的形式给出，例如我们经常使用的概率积分使用的概率积分 就无法利用牛顿就无法利用牛顿- -莱布尼兹公式推导得到。莱布尼兹公式推导得到。相比之下，数

17、值计算在实际计算中有符号运算无法比拟的优点。相比之下，数值计算在实际计算中有符号运算无法比拟的优点。一是效率高，对公式稍微复杂一些或规模稍微大一些的问题，只一是效率高，对公式稍微复杂一些或规模稍微大一些的问题，只能利用数值计算。二是我们经常遇到的计算问题，都可以利用数能利用数值计算。二是我们经常遇到的计算问题，都可以利用数值计算得到满意的解，而符号运算只能处理其中一小部分问题。值计算得到满意的解，而符号运算只能处理其中一小部分问题。baxdxe2数值计算需要对连续变量离散求解！数值计算需要对连续变量离散求解！数值计算需要设计数值算法数值计算需要设计数值算法数值计算是在离散节点上的近似计算，通常

18、对连续的区间要首先数值计算是在离散节点上的近似计算，通常对连续的区间要首先进行离散。对计算区间进行离散。对计算区间a.ba.b上的导数上的导数。首先将区间离散成等距的节点首先将区间离散成等距的节点其中步长其中步长 与下标与下标k k无关。无关。 bxxxxan210kkxxnabh1一阶和二阶一阶和二阶导数的计算由公式导数的计算由公式 给出给出, , 误差为误差为 O(h O(h2 2) )。 hxfxfxfkkk2)()()(11211)()(2)()(hxfxfxfxfkkkk MatlabMatlab的数值导数计算可以利用语句的数值导数计算可以利用语句diffdiff实现。实现。diff

19、diff语句对非符语句对非符号计算表示向前差分，即对号计算表示向前差分，即对 函数函数diff(x)diff(x)输出输出diff(x,2)diff(x,2)输出输出 diff(diff(x), diff(diff(x),即即 一般一般, diff(x,k) , diff(x,k) 的计算可依次类推。的计算可依次类推。,210nxxxxx,11201nnxxxxxx2,221012nnnxxxxxx注：注： diff(x) diff(x) 的元素个数比的元素个数比x x少少 1 1！一阶导数的计算一阶导数的计算一阶导数可以利用一阶导数可以利用 计算。由计算。由公公式式dfdf的第的第i i个元

20、素与个元素与 f(x) f(x) 在在 x xi+1/2i+1/2点的值最接近，误差是点的值最接近，误差是O(hO(h2 2) )(/).(xdiffydiffdf hxfxfxfkkk2)()()(11例13：利用diff语句计算函数 y=sin(x) 的导数，并计算误差。参考程序：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)legen

21、d(df,f)legend(error) 结果图：二阶导数的计算二阶导数的计算语句语句 输出二阶向前差分，即输出二阶向前差分，即 可以利用可以利用 来计算来计算y(x)y(x)的二阶导数。的二阶导数。)2 ,(xdiffy iiiixxxy1222, 2 , 1ni 2).(/).2 ,(xdiffydiff计算的是计算的是n-2n-2个内点上的二阶导数，不包含两个端点上的导数！个内点上的二阶导数，不包含两个端点上的导数！例14：计算 y=sin(x)的二阶导数，并与 sin(x)比较。参考程序：x=0:pi/100:2*pi;y=sin(x);dx=diff(x);d2y=diff(y,2)

22、./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上的定积分。 I=quad(fun,a,b,tol) %计算函数fun在a,b上精度tol的定积分。funfun的的 取法：取法： 过去的取法是一个字符串作为所调用的函数名字

23、。如过去的取法是一个字符串作为所调用的函数名字。如 quad(sin, 0, pi/2) quad(sin, 0, pi/2) 在在matlab7matlab7以上的版本中引入了函数句柄，以上的版本中引入了函数句柄，可以在函数可以在函数quadquad中利用中利用 指示一个函数句柄。指示一个函数句柄。例如例如 sin %sin %调用以调用以sinsin为函为函数名的函数子程序数名的函数子程序 (x,y)(x.3+y.3) % (x,y)(x.3+y.3) %建立并调用一个函数句建立并调用一个函数句柄柄 (x)myfun(x,a) % (x)myfun(x,a) %调用含有参数的子程序调用含有

24、参数的子程序myfun myfun 例15：计算积分 k=1,2,3,4,5并与精确值比较。2/0)(sindxxkI=zeros(1,5);I=zeros(1,5);for k=1:5for k=1:5 g=strcat(sin(x).,num2str(k); g=strcat(sin(x).,num2str(k); I(k)=quad(g,0,pi/2); I(k)=quad(g,0,pi/2);endendI0=1,pi/4,2/3,3I0=1,pi/4,2/3,3* *pi/16,8/15;pi/16,8/15;subplot(1,2,1)subplot(1,2,1)plot(I-I0

25、,p)plot(I-I0,p)legend(errors of quad)legend(errors of quad) for k=1:5 for k=1:5 g=strcat(sin(x).,num2str(k); g=strcat(sin(x).,num2str(k); I(k)=quad(g,0,pi/2,1e-12); I(k)=quad(g,0,pi/2,1e-12); end end subplot(1,2,2)subplot(1,2,2)plot(I-I0,p)plot(I-I0,p)legend(errors of quad)legend(errors of quad)积分函数

26、积分函数quadvquadv积分函数积分函数quadvquadv类似于函数类似于函数quadquad，但它可以对含有向量参数的函数，但它可以对含有向量参数的函数积分，而且可以对矩阵函数积分。积分，而且可以对矩阵函数积分。例例1616：计算积分：计算积分k=1,2,3,4,5k=1,2,3,4,52/0)(sindxxk这里这里k k 可以看做参数向量可以看做参数向量! !参考语句： format rat format rat I=quadv(sin(x).1:5,0,pi/2) I=quadv(sin(x).1:5,0,pi/2)I = 1 355/452 2/3 1065/1808 8/15

27、例17：计算矩阵函数 的积分 xxxxxxxexFxsin)1ln(11cos)(221)(dxxF参考语句： format rat format rat I=quadv(x)exp(-x). I=quadv(x)exp(-x).* *cos(x),(1-x)./sqrt(1+x.2);cos(x),(1-x)./sqrt(1+x.2);x.x.* *log(1+x),sin(x)./x,1,2)log(1+x),sin(x)./x,1,2)I = 249/7261 -548/2111 1075/769 689/1045 可用于广义积分的积分函数可用于广义积分的积分函数quadgkquadgk

28、函数函数quadgkquadgk采用采用Gauss-KronrodGauss-Kronrod数值积分方法，可以解决无穷区数值积分方法，可以解决无穷区间、无穷间断点等问题。其使用方法与间、无穷间断点等问题。其使用方法与quadquad类似类似，但积分限可以，但积分限可以取取infinf例例1818：计算：计算 参考语句：参考语句： I=quadgk(x)1./(3 I=quadgk(x)1./(3* *x.2+2),1,inf)x.2+2),1,inf)I =I =0.27950.27951223xdx多重积分多重积分MatlabMatlab的多重积分函数为的多重积分函数为 I=dblquad(fun,xmin,xmax,ymin,ymax) % I=dblquad(fun,xmin,x