利用matlab進(jìn)行微積分的計(jì)算

1、MatlabMatlab定義符號(hào)運(yùn)算和數(shù)值運(yùn)算。如果運(yùn)算過(guò)程中有符號(hào)變量，則定義符號(hào)運(yùn)算和數(shù)值運(yùn)算。如果運(yùn)算過(guò)程中有符號(hào)變量，則運(yùn)算過(guò)程為符號(hào)運(yùn)算運(yùn)算過(guò)程為符號(hào)運(yùn)算。符號(hào)變量可以利用兩種方法定義。符號(hào)變量可以利用兩種方法定義。(1) (1) 利用利用函數(shù)函數(shù)symsym定義符號(hào)表達(dá)式定義符號(hào)表達(dá)式語(yǔ)句格式語(yǔ)句格式 變量變量=sym(=sym(字符串字符串) )symsym函數(shù)可以定義符號(hào)變量、符號(hào)常數(shù)和一般的符號(hào)表達(dá)式函數(shù)可以定義符號(hào)變量、符號(hào)常數(shù)和一般的符號(hào)表達(dá)式。 x=sym(1/3) x=sym(1/3) x = x =1/31/3 y=x+1/2 %y=x+1/2 %符號(hào)運(yùn)算符號(hào)運(yùn)算

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 %數(shù)值計(jì)數(shù)值計(jì)算算y =y = 0.8333 0.8333例1：比較符號(hào)運(yùn)算和數(shù)值計(jì)算例例2 2：分別利用符號(hào)運(yùn)算和數(shù)值運(yùn)算計(jì)算如下的：分別利用符號(hào)運(yùn)算和數(shù)值運(yùn)算計(jì)算如下的和和 并比較其計(jì)算速度。并比較其計(jì)算速度。1000131211S注：注：MatlabMatlab語(yǔ)言中利用語(yǔ)言中利用tic tic 和和toc toc 計(jì)算時(shí)間，語(yǔ)言格式為計(jì)算時(shí)間，語(yǔ)言格式為 tic tic 程序塊程序塊 toc toc顯示器顯示顯示器顯示tictic和和toctoc之間的程序塊運(yùn)

3、行時(shí)間之間的程序塊運(yùn)行時(shí)間符號(hào)運(yùn)算程序：s=sym(0);ticfor k=1:1000 s=s+1/k;endtocs運(yùn)行結(jié)果：elapsed time is 17.471170 seconds. s = 5336291328229478504559104562404298040965247228038426009710134924845626888949710175750609790198503569140908873155046809837844217211788500946430234432656602250210027842563285208140554494121044251014

4、2672770294774712708917963967779610453224692426866468888281582071984897105110796873249319155529397017508931564519976085734473014183284011724412280649074307703736683170055800293659235088589360235285852808160759574737836655413175508131522517/7128865274665093053166384155714272920668358861885893040452001

5、991154324087581111499476444151913871586911717817019575256512980264067621009251465871004305131072686268143200196609974862745937188343705015434452523739745298963145674982128236956232823794011068809262317708861979540791247754558049326475737829923352751796735248042463638051137034331214781746850878453485

6、678021888075373249921995672056932029099390891687487672697950931603520000數(shù)值運(yùn)算程序數(shù)值運(yùn)算程序 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運(yùn)行結(jié)果：運(yùn)行結(jié)果：Elapsed time is 0.000015 secondsElapsed time is 0.000015 seconds. .s s = =7.48557.4855可以看出，數(shù)值計(jì)算的計(jì)算速度遠(yuǎn)遠(yuǎn)高于符號(hào)運(yùn)算可以看出，數(shù)值計(jì)算的計(jì)算速度遠(yuǎn)遠(yuǎn)高

7、于符號(hào)運(yùn)算。另外，符號(hào)運(yùn)算只能計(jì)算非常簡(jiǎn)單的問(wèn)題，復(fù)雜問(wèn)題只能利用數(shù)值另外，符號(hào)運(yùn)算只能計(jì)算非常簡(jiǎn)單的問(wèn)題，復(fù)雜問(wèn)題只能利用數(shù)值方法。但由于符號(hào)運(yùn)算得到的是解析式，在許多情況下有利于問(wèn)題方法。但由于符號(hào)運(yùn)算得到的是解析式，在許多情況下有利于問(wèn)題的進(jìn)一步分析。的進(jìn)一步分析。利用語(yǔ)句利用語(yǔ)句 syms syms 定義符號(hào)變量定義符號(hào)變量語(yǔ)句形式語(yǔ)句形式 syms syms x y x y z zsymssyms語(yǔ)句語(yǔ)句 一次可定義多個(gè)符號(hào)變量，但不能定義符號(hào)常數(shù)和表達(dá)式一次可定義多個(gè)符號(hào)變量，但不能定義符號(hào)常數(shù)和表達(dá)式極限運(yùn)算語(yǔ)句格式 G= limit(F); %表達(dá)式F中變量趨向于0的極限 G=

8、 limit(F,a); %表達(dá)式F中變量趨向于a的極限 G= limit(F,v,a); %表達(dá)式F中變量v趨向于a的極限 G= limit(F,v,a,right); %表達(dá)式F中變量趨向于a的右極限。例3：計(jì)算 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：計(jì)算 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導(dǎo)數(shù)求導(dǎo)數(shù)的語(yǔ)句

9、為 df=diff(F) %求表達(dá)式F的一階導(dǎo)數(shù) df=diff(F,n) %求表達(dá)式F的n階導(dǎo)數(shù) df=diff(F,u,n) %求表達(dá)式F的關(guān)于變量u的n階導(dǎo)數(shù)例5：求函數(shù) f=sin(x)/x 的二階導(dǎo)數(shù)。 f=sin(x)/xf = sin(x)/x diff(f,2)ans =(2*sin(x)/x3 - sin(x)/x - (2*cos(x)/x2例6:求函數(shù) 關(guān)于x和y的二階偏導(dǎo)數(shù)。f=(x-y)/(x+y)f = (x - y)/(x + y) dfxy=diff(diff(f,y),x)dfxy =(2*(x - y)/(x + y)3yxyxyxf),(注：如果表達(dá)式注：

10、如果表達(dá)式F F有多個(gè)符號(hào)變量，且沒(méi)有指明對(duì)哪一個(gè)變量求導(dǎo)，則有多個(gè)符號(hào)變量，且沒(méi)有指明對(duì)哪一個(gè)變量求導(dǎo)，則MatlabMatlab規(guī)定對(duì)離規(guī)定對(duì)離x x最近的變量求導(dǎo)。其他運(yùn)算也有類似規(guī)定。最近的變量求導(dǎo)。其他運(yùn)算也有類似規(guī)定。例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的微積分符號(hào)運(yùn)算都可以對(duì)數(shù)組進(jìn)行。的微積分符號(hào)運(yùn)算都可以對(duì)數(shù)組進(jìn)行。函數(shù)的積分函數(shù)的積分積分符號(hào)運(yùn)算的基本語(yǔ)句積分符號(hào)運(yùn)算的基本語(yǔ)句 int(F) int(F)； % %求函數(shù)表達(dá)式求函數(shù)表達(dá)式F F的不定積分的不定積分 int(F,v) int(F,v)； % %求函數(shù)表達(dá)式求函數(shù)表達(dá)式F F關(guān)于變量關(guān)于變量v v的不的不定積分定積分 int(F,a,b) int(F,a,b)； % %求函數(shù)表達(dá)式求函數(shù)表達(dá)式F F在區(qū)間在區(qū)間a,ba,b上上的定積分的定積分 int(F,v,a,b); % int(F,v,a,b); %求函數(shù)表達(dá)式求函數(shù)表達(dá)式F F在區(qū)間在區(qū)間a,ba,b上的關(guān)于變量

12、上的關(guān)于變量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對(duì)含參數(shù)的運(yùn)算，符號(hào)計(jì)算的優(yōu)勢(shì)在于可以提供一個(gè)用于進(jìn)一對(duì)含參數(shù)的運(yùn)算，符號(hào)計(jì)算的優(yōu)勢(shì)在于可以提供一個(gè)用于進(jìn)一步分析的解析表達(dá)式。步分析的解析表達(dá)式。級(jí)數(shù)的符號(hào)計(jì)算級(jí)數(shù)的符號(hào)計(jì)算級(jí)數(shù)符號(hào)計(jì)算的基本語(yǔ)句級(jí)數(shù)符號(hào)計(jì)算的基本語(yǔ)句 S=sym

13、sum(g) S=symsum(g)； % % 求通項(xiàng)求通項(xiàng)g g的無(wú)窮和的無(wú)窮和 S=symsum(g,n) S=symsum(g,n)； % % 求通項(xiàng)求通項(xiàng)g g的對(duì)整變量的對(duì)整變量n n的無(wú)窮和的無(wú)窮和 S=symsum(g,a,b) S=symsum(g,a,b)； % % 求通項(xiàng)求通項(xiàng)g g的下標(biāo)在的下標(biāo)在a a和和b b之間的之間的項(xiàng)的和項(xiàng)的和 S=symsum(g,n,a,b) S=symsum(g,n,a,b)； % % 求通項(xiàng)求通項(xiàng)g g的下標(biāo)變量的下標(biāo)變量n n在在a a和和b b之間的項(xiàng)的之間的項(xiàng)的和和例10:求級(jí)數(shù)的和 syms n symsum(1/n2,1,inf

14、) ans = pi2/6 121nn例11: 求冪級(jí)數(shù)的和 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展式的基本語(yǔ)句展式的基本語(yǔ)句 F=taylor(f); % F=taylor(f); %求函數(shù)求函數(shù)f f 的馬克的馬克勞林級(jí)數(shù)的前勞林級(jí)數(shù)的前6 6項(xiàng)項(xiàng) F=taylor(f,n); % F=taylor(f,n); %求函數(shù)求函數(shù)f f 的馬克的馬克勞林級(jí)數(shù)的前勞

15、林級(jí)數(shù)的前n n項(xiàng)項(xiàng) F=taylor(f,n,x F=taylor(f,n,x0); %); %求函數(shù)求函數(shù)f f 的在的在x x0展展開(kāi)的開(kāi)的TaylorTaylor級(jí)數(shù)的前級(jí)數(shù)的前n n項(xiàng)項(xiàng) F=taylor(f,v,n,v F=taylor(f,v,n,v0); %); %求函數(shù)求函數(shù)f f對(duì)變量對(duì)變量v v在在v v0展展開(kāi)的開(kāi)的TaylorTaylor級(jí)數(shù)的前級(jí)數(shù)的前n n項(xiàng)項(xiàng)例例12: 12: 求函數(shù)求函數(shù) 的在的在0點(diǎn)的點(diǎn)的taylortaylor展式展式前前8 8項(xiàng)項(xiàng) 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)(符號(hào)運(yùn)算有很大的局限性。第一是計(jì)算是純推理過(guò)程，效率非常符號(hào)運(yùn)算有很大的局限性。第一是計(jì)算是純推理過(guò)程，效率非常低。第二是大量的問(wèn)題無(wú)法利用解析的形式給出，例如我們經(jīng)常低。第二是大量的問(wèn)題無(wú)法利用解析的形式給出，例如我們經(jīng)常使用的概率積分使用的概率積分 就無(wú)法利用牛頓就無(wú)法利用牛頓- -萊布尼茲公式推導(dǎo)得到。萊布尼茲公式推導(dǎo)得到。相比之下，數(shù)

17、值計(jì)算在實(shí)際計(jì)算中有符號(hào)運(yùn)算無(wú)法比擬的優(yōu)點(diǎn)。相比之下，數(shù)值計(jì)算在實(shí)際計(jì)算中有符號(hào)運(yùn)算無(wú)法比擬的優(yōu)點(diǎn)。一是效率高，對(duì)公式稍微復(fù)雜一些或規(guī)模稍微大一些的問(wèn)題，只一是效率高，對(duì)公式稍微復(fù)雜一些或規(guī)模稍微大一些的問(wèn)題，只能利用數(shù)值計(jì)算。二是我們經(jīng)常遇到的計(jì)算問(wèn)題，都可以利用數(shù)能利用數(shù)值計(jì)算。二是我們經(jīng)常遇到的計(jì)算問(wèn)題，都可以利用數(shù)值計(jì)算得到滿意的解，而符號(hào)運(yùn)算只能處理其中一小部分問(wèn)題。值計(jì)算得到滿意的解，而符號(hào)運(yùn)算只能處理其中一小部分問(wèn)題。baxdxe2數(shù)值計(jì)算需要對(duì)連續(xù)變量離散求解！數(shù)值計(jì)算需要對(duì)連續(xù)變量離散求解！數(shù)值計(jì)算需要設(shè)計(jì)數(shù)值算法數(shù)值計(jì)算需要設(shè)計(jì)數(shù)值算法數(shù)值計(jì)算是在離散節(jié)點(diǎn)上的近似計(jì)算，通常

18、對(duì)連續(xù)的區(qū)間要首先數(shù)值計(jì)算是在離散節(jié)點(diǎn)上的近似計(jì)算，通常對(duì)連續(xù)的區(qū)間要首先進(jìn)行離散。對(duì)計(jì)算區(qū)間進(jìn)行離散。對(duì)計(jì)算區(qū)間a.ba.b上的導(dǎo)數(shù)上的導(dǎo)數(shù)。首先將區(qū)間離散成等距的節(jié)點(diǎn)首先將區(qū)間離散成等距的節(jié)點(diǎn)其中步長(zhǎng)其中步長(zhǎng) 與下標(biāo)與下標(biāo)k k無(wú)關(guān)。無(wú)關(guān)。 bxxxxan210kkxxnabh1一階和二階一階和二階導(dǎo)數(shù)的計(jì)算由公式導(dǎo)數(shù)的計(jì)算由公式 給出給出, , 誤差為誤差為 O(h O(h2 2) )。 hxfxfxfkkk2)()()(11211)()(2)()(hxfxfxfxfkkkk MatlabMatlab的數(shù)值導(dǎo)數(shù)計(jì)算可以利用語(yǔ)句的數(shù)值導(dǎo)數(shù)計(jì)算可以利用語(yǔ)句diffdiff實(shí)現(xiàn)。實(shí)現(xiàn)。diff

19、diff語(yǔ)句對(duì)非符語(yǔ)句對(duì)非符號(hào)計(jì)算表示向前差分，即對(duì)號(hào)計(jì)算表示向前差分，即對(duì) 函數(shù)函數(shù)diff(x)diff(x)輸出輸出diff(x,2)diff(x,2)輸出輸出 diff(diff(x), diff(diff(x),即即 一般一般, diff(x,k) , diff(x,k) 的計(jì)算可依次類推。的計(jì)算可依次類推。,210nxxxxx,11201nnxxxxxx2,221012nnnxxxxxx注：注： diff(x) diff(x) 的元素個(gè)數(shù)比的元素個(gè)數(shù)比x x少少 1 1！一階導(dǎo)數(shù)的計(jì)算一階導(dǎo)數(shù)的計(jì)算一階導(dǎo)數(shù)可以利用一階導(dǎo)數(shù)可以利用 計(jì)算。由計(jì)算。由公公式式dfdf的第的第i i個(gè)元

20、素與個(gè)元素與 f(x) f(x) 在在 x xi+1/2i+1/2點(diǎn)的值最接近，誤差是點(diǎn)的值最接近，誤差是O(hO(h2 2) )(/).(xdiffydiffdf hxfxfxfkkk2)()()(11例13：利用diff語(yǔ)句計(jì)算函數(shù) y=sin(x) 的導(dǎo)數(shù)，并計(jì)算誤差。參考程序：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) 結(jié)果圖：二階導(dǎo)數(shù)的計(jì)算二階導(dǎo)數(shù)的計(jì)算語(yǔ)句語(yǔ)句 輸出二階向前差分，即輸出二階向前差分，即 可以利用可以利用 來(lái)計(jì)算來(lái)計(jì)算y(x)y(x)的二階導(dǎo)數(shù)。的二階導(dǎo)數(shù)。)2 ,(xdiffy iiiixxxy1222, 2 , 1ni 2).(/).2 ,(xdiffydiff計(jì)算的是計(jì)算的是n-2n-2個(gè)內(nèi)點(diǎn)上的二階導(dǎo)數(shù)，不包含兩個(gè)端點(diǎn)上的導(dǎo)數(shù)！個(gè)內(nèi)點(diǎn)上的二階導(dǎo)數(shù)，不包含兩個(gè)端點(diǎn)上的導(dǎo)數(shù)！例14：計(jì)算 y=sin(x)的二階導(dǎo)數(shù)，并與 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)數(shù)值積分一元數(shù)值積分Matlab中的數(shù)值積分函數(shù)quad采用自適應(yīng)復(fù)化Simpson公式，語(yǔ)句格式為 I=quad(fun,a,b) %計(jì)算函數(shù)fun在a,b上的定積分。 I=quad(fun,a,b,tol) %計(jì)算函數(shù)fun在a,b上精度tol的定積分。funfun的的 取法：取法： 過(guò)去的取法是一個(gè)字符串作為所調(diào)用的函數(shù)名字

23、。如過(guò)去的取法是一個(gè)字符串作為所調(diào)用的函數(shù)名字。如 quad(sin, 0, pi/2) quad(sin, 0, pi/2) 在在matlab7matlab7以上的版本中引入了函數(shù)句柄，以上的版本中引入了函數(shù)句柄，可以在函數(shù)可以在函數(shù)quadquad中利用中利用 指示一個(gè)函數(shù)句柄。指示一個(gè)函數(shù)句柄。例如例如 sin %sin %調(diào)用以調(diào)用以sinsin為函為函數(shù)名的函數(shù)子程序數(shù)名的函數(shù)子程序 (x,y)(x.3+y.3) % (x,y)(x.3+y.3) %建立并調(diào)用一個(gè)函數(shù)句建立并調(diào)用一個(gè)函數(shù)句柄柄 (x)myfun(x,a) % (x)myfun(x,a) %調(diào)用含有參數(shù)的子程序調(diào)用含有

24、參數(shù)的子程序myfun myfun 例15：計(jì)算積分 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)積分函數(shù)

26、積分函數(shù)quadvquadv積分函數(shù)積分函數(shù)quadvquadv類似于函數(shù)類似于函數(shù)quadquad，但它可以對(duì)含有向量參數(shù)的函數(shù)，但它可以對(duì)含有向量參數(shù)的函數(shù)積分，而且可以對(duì)矩陣函數(shù)積分。積分，而且可以對(duì)矩陣函數(shù)積分。例例1616：計(jì)算積分：計(jì)算積分k=1,2,3,4,5k=1,2,3,4,52/0)(sindxxk這里這里k k 可以看做參數(shù)向量可以看做參數(shù)向量! !參考語(yǔ)句： 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：計(jì)算矩陣函數(shù) 的積分 xxxxxxxexFxsin)1ln(11cos)(221)(dxxF參考語(yǔ)句： 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 可用于廣義積分的積分函數(shù)可用于廣義積分的積分函數(shù)quadgkquadgk

28、函數(shù)函數(shù)quadgkquadgk采用采用Gauss-KronrodGauss-Kronrod數(shù)值積分方法，可以解決無(wú)窮區(qū)數(shù)值積分方法，可以解決無(wú)窮區(qū)間、無(wú)窮間斷點(diǎn)等問(wèn)題。其使用方法與間、無(wú)窮間斷點(diǎn)等問(wèn)題。其使用方法與quadquad類似類似，但積分限可以，但積分限可以取取infinf例例1818：計(jì)算：計(jì)算 參考語(yǔ)句：參考語(yǔ)句： 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的多重積分函數(shù)為的多重積分函數(shù)為 I=dblquad(fun,xmin,xmax,ymin,ymax) % I=dblquad(fun,xmin,x