第三講 (2010)MATLAB在微積分運算中的應(yīng)用

第三講(2010)MATLAB在微積分運算中的應(yīng)用第三講MATLAB在微積分計算上的應(yīng)用MATLAB在微積分計算中的常用命令運算命令求極限limit(F(x),x,a)求導(dǎo)diff(function,'vaiaale',n)求極小點fminbnd(f,a,b)，fminsearch(f(x),x0)求不定積分int(function,vaiaale)求定積分int(function,vaiaale,a,b)求梯度gradient(),求方向?qū)?shù)Jacobian(f)?v級數(shù)求和symsum(function,vaiaale,a,b)泰勒級數(shù)展開taylor(function,n,vaiaale,a)解微分方程dsolve('d_equa1',‘condi',‘var')21.極限(1)對極限定義的認識練習(xí)1用下面的語句觀察數(shù)列的前100項變化情況n=1:100;a=n.^(n.^(-1))為了更清楚地觀察其是否收斂，讀者可將項數(shù)增大一些對于該數(shù)列，我們再用語句：plot(n,a,'.')畫出其散點圖，借助于圖形來觀察它的變化趨勢。3數(shù)列的散點圖從上圖可看出，這個數(shù)列似乎收斂于1.但如何說明它收斂于1，而不是收斂于大于1的某個數(shù)呢（由，若極限存在，則極限必不小于1）？設(shè)該數(shù)列收斂于A=1+u(u

0),我們?nèi)=10-2用程序來檢察接近的程度。4u=10^-2;A=1+u;m=5;n=2;an=sqrt(2);whileabs(an-A)>=10^-mn=n+1;an=n^(1/n);endfprintf('A=%3.2f,n=%3.0f,an=%3.0f,abs(an-A)=%8.7e\n',A,n,an,abs(an-A))結(jié)果為:A=1.01,n=651,an=1,abs(an-A)=1.3098309e-006這說明當(dāng)=651，an=1.01時,an與1+10-2的距離小于10-5。5(2)極限的計算在MATLAB軟件中可以直接用命令limit來求極限，其一般格式是:這里需要說明幾點1.上式求的是符號表達式F(x)當(dāng)x

a時的極限值，若要計算右極限或左極限，可在后指明趨向的方向；練習(xí)7試比較下面語句的區(qū)別.limit(exp(-1/x),x,0)limit(exp(-1/x),x,0,'right')limit(exp(-1/x),x,0,'left')2.在試圖求無窮振蕩點處的極限時，limit語句得到的是函數(shù)振蕩時可能的取值范圍；limit(F(x),x,a)6例1求下列函數(shù)的極限：clearF1=sym('atan(x)/x');F2=sym('((1+x)/(1-x))^(1/x)');F3=sym('(sqrt(1+x^2)-1)/(1-cos(x))');F4=sym('x*log(1+x)/sin(x^2)');F=[F1,F2,F3,F4]limit(F)ans=

[1,exp(2),1,1]7例2求函數(shù)的極限：clearsymsaFxF=(1+a/x)^xlimit(F,'x',inf,'left')ans=

exp(a)8(3)一些數(shù)列的極限的討論設(shè)數(shù)列xn與yn由下式確定：xn與yn的極限存在嗎？9運行該程序可判斷出：xn與yn有極限，且這兩極限值是相等的。（x100=1.456791E+000,y100=1.456791E+000。

用MATLAB軟件編出如下程序進行觀察：xn=1;yn=2;forn=2:1:100xN=xn;yN=yn;xn=sqrt(xN*yN);yn=(xN+yN)/2;endfprintf('x100=%E,y100=%E\n',xn,yn)102.導(dǎo)數(shù)與微分在MATLAB中由命令函數(shù)diff()來完成運算,其具體形式為：diff(function,'vaiaale',n)參數(shù)function為需要進行求導(dǎo)運算的函數(shù),vaiaale為求導(dǎo)運算的獨立變量,n為求導(dǎo)的階數(shù)。11(1)求導(dǎo)命令diff在MATLAB軟件中，可用語句diff(f(x),x)計算函數(shù)的導(dǎo)函數(shù)，當(dāng)然在使用前需先將x定義成符號變量。若要求f(x)在x=a處的導(dǎo)數(shù)，可用subs命令，只要將x=a賦給上面的導(dǎo)函數(shù)便可得到。而命令diff(f(x),x,n)求的是函數(shù)f(x)對x的n階導(dǎo)函數(shù)。12例3求下列函數(shù)的導(dǎo)數(shù)clearsymsxy1y2y3y1=(cos(x))^3-cos(3*(x));y2=x*sin(x)*log(x);y3=(x*exp(x)-1)/sin(x);dy1=diff(y1)dy2=diff(y2)dy3=diff(y3)13dy1=

-3*cos(x)^2*sin(x)+3*sin(3*x)dy2=

sin(x)*log(x)+x*cos(x)*log(x)+sin(x)dy3=

(exp(x)+x*exp(x))/sin(x)-(x*exp(x)-1)/sin(x)^2*cos(x)pretty(dy3)exp(x)+xexp(x)(xexp(x)-1)cos(x)---------------------------------------sin(x)sin2(x)

14clearsymsxdiff(2*sqrt(x),x);x=2;f=1/x^(1/2)ans=

1/x^(1/2)f=0.707115clearsymsxy=1/(1-log(x));dy=diff(y,x,6)dy=

720/(1-log(x))^7/x^6-1800/(1-log(x))^6/x^6+2040/(1-log(x))^5/x^6-1350/(1-log(x))^4/x^6+548/(1-log(x))^3/x^6-120/(1-log(x))^2/x^616(2)隱函數(shù)與由參數(shù)方程確定的函數(shù)的導(dǎo)數(shù)求隱函數(shù)的導(dǎo)數(shù)與由參數(shù)方程確定的函數(shù)的導(dǎo)數(shù)，需將求導(dǎo)命令與數(shù)學(xué)公式或方法結(jié)合起來，才能奏效！例5求由方程xy-ex+ey=0確定的函數(shù)y=f(x)的導(dǎo)數(shù)。與筆算的做法一樣，先在方程兩邊對變量求導(dǎo)，再從所得方程中解出即可。clearsymsxyf=(x*y-exp(x)+exp(y))-diff(f,x)/diff(f,y)運行之后可求出y'(x)。這兩個步驟分別可由以下語句完成：ans=

(-y+exp(x))/(x+exp(y))17由下面的語句可得所求的導(dǎo)數(shù)：yans=

sin(t)/(1-cos(t))結(jié)果為:183.導(dǎo)數(shù)的應(yīng)用(1)最值的計算a.直接用MATLAB語句計算MATLAB軟件中提供了求函數(shù)極小點的語句fminbnd(f,a,b)執(zhí)行該語句將得到函數(shù)在區(qū)間[a,b]內(nèi)的極小點.而語句fminsearch(f(x),x0)得到的是離x0最近的極小點.19要求函數(shù)的極大點，可以用命令fminsearch(-f(x),x0)這是因為函數(shù)-f(x)的極小點恰好就是函數(shù)f(x)的極大點．例7求函數(shù)y=2x3-6x2-18x+7的極值。解先作圖了解。clearx=-5:0.1:5;y=2*x.^3-6*x.^2-18*x+7;plot(x,y)20從圖上可知,函數(shù)y在區(qū)間[-5,5]上有極大值與極小值。f=inline('2*x.^3-6*x.^2-18*x+7');

pmin=fminbnd(f,-5,5)pmin=-5g=inline('-2*x.^3+6*x.^2+18*x-7');pmax=fminbnd(g,-5,5)pmax=5fprintf('%g,%g,%g,%g\n',pmin,f(pmin),pmax,f(pmax))-5,-303,5,17因5與-5皆在區(qū)間[-53，5]內(nèi)，故所求的最大值為17，最小值為-303。21b.利用導(dǎo)數(shù)計算由高等數(shù)學(xué)的知識可知，函數(shù)y=f(x)的駐點與一階不可導(dǎo)的孤立點可能是函數(shù)的極值點。若這些點的個數(shù)有限，我們只要比較這些點與區(qū)間端點、處的函數(shù)值，便能求出函數(shù)在[a,b]上的最大值與最小值了。clearsymsxy=2*x.^3-6*x.^2-18*x+7;

dy=diff(y)dy=

6*x^2-12*x-1822Px=solve(dy)Px=

3-1ezplot(y)234.不定積分在MATLAB中由命令函數(shù)int()來完成積分運算,其具體形式為：int(function,vaiaale)參數(shù)function為需要進行求導(dǎo)運算的函數(shù),vaiaale為積分變量。24clearsymsxy=x^5+x^3-sqrt(x)/4y=

x^5+x^3-1/4*x^(1/2)int(y)ans=

1/6*x^6+1/4*x^4-1/6*x^(3/2)pretty(ans)pretty(ans)

643/21/6x+1/4x-1/6x25clearsymsx;y=sin(2*x)*sin(3*x)*sin(4*x)ans=

-1/20*cos(5*x)-1/12*cos(3*x)+1/36*cos(9*x)-1/4*cos(x)int(y)pretty(ans)-1/20cos(5x)-1/12cos(3x)+1/36cos(9x)-1/4cos(x)y=

sin(2*x)*sin(3*x)*sin(4*x)265.定積分在MATLAB中由命令函數(shù)int()來完成積分運算,其具體形式為：int(function,vaiaale,a,b)參數(shù)function為需要進行求導(dǎo)運算的函數(shù),vaiaale為積分變量,a、b分別為積分下、上限。(1)用命令函數(shù)計算定積分27clearsymsx;y=(x*exp(x))/(1+x)^2ans=

1/2*exp(1)-1int(y,0,1)y=

x*exp(x)/(1+x)^228clearsymsx;y=1/(x^2+2*x+3);int(y,-inf,inf)ans=

1/2*pi*2^(1/2)29(2)由定義計算定積分在定積分的定義中，劃分積分區(qū)間的方法與在每個小區(qū)間上取的點都是任意的，求積分和極限時要求每個小區(qū)間長度的最大值趨向于0，這些都給我們直接由定義來驗證一個定積分是否存在帶來了很大的困難?，F(xiàn)在我們借助于計算機，按照定義的要求，對積分和的極限作近似計算，根據(jù)結(jié)果對定積分是否存在作出判斷。下面以積分為例來說明這個問題：a.首先在區(qū)間[0,1]中插入n-1個分點，為使分點任意，可用能產(chǎn)生隨機數(shù)的函數(shù)rand()。b.為保證分割加細時，各小區(qū)間的長度趨于0，在取分點時，讓相鄰兩分點間的距離小于2/n。30在下列給出的程序中分點取為因0≤ui≤1,故c.其次在每個小區(qū)間[xi-1,xi](i=1,2,…n)上任取一點

i，為使

i具有任意性，我們同樣利用函數(shù)rand()來實現(xiàn)，31為了提高精確度，我們讓分點不斷增多反復(fù)進行計算。計算程序如下：clearallf=inline('x^2');a=0;b=1;n=20;%n為分割成的小區(qū)間個數(shù),初始值取為20x=[];x(1)=a;fork=1:6x(n+1)=b;s=0;fori=1:n-1x(i+1)=(i+rand())*(b-a)/n;%取區(qū)間的分割點end32fori=1:ndxi=x(i+1)-x(i);%計算第i個區(qū)間的長度c=x(i)+dxi*rand();%在第i個區(qū)間上任取一點s=s+f(c)*dxi;%逐步求積分和endfprintf('n=%g,s=%g\n',n,s);n=n*2;end程序中分割小區(qū)間的個數(shù)n的初值取為20，然后每循環(huán)一次放大一倍，共放大了5次。程序運行結(jié)果為：這樣做的目的是為了盡快獲得結(jié)果，當(dāng)然我們也可取其它的值作為的初值，也可以用其它的方式讓增大。33n=20,s=0.340088n=40,s=0.335021n=80,s=0.333965n=160,s=0.333266n=320,s=0.333431n=640,s=0.333373由分割的任意性及

i的任意性，我們有理由認為,即使n固定，每次運行該程序所得的結(jié)果也很可能是不同的。事實上，在實驗時得到了下表34第1次第2次第3次第4次第5次200.3318340.3310790.3346160.3270940.337772400.333490.3328930.3333720.3327080.333732800.334510.3333370.3334970.3332840.3326611600.3334730.3330610.3333940.333430.3331243200.3333460.3334620.3332340.3334690.3333786400.3333270.3333650.3333430.3333590.333333表中任何兩個數(shù)據(jù)都不完全相同，但可看出它們間的差異不是很大，特別是最后一行當(dāng)n=640時，5次運行的結(jié)果前4位有效數(shù)字是一樣的，故我們猜測：35(3)從圖形觀察積分和與定積分的關(guān)系從圖形上來觀察隨著分割點的增多，積分和是否越來越接近定積分的值。下面以積分為例來說明：36clearxab;f=inline('sin(x)');a=0;b=pi/2;n=0;symsx;axismanual;set(gca,'nextplot','replacechildren');forj=3:20:103n=j;t(1)=a;t(n+1)=b; fori=1:n-1t(i+1)=(i+rand())*(b-a)/n;end ezplot(f(x),[a,b]);holdon

37fori=1:nc=t(i)+(t(i+1)-t(i))*rand(); bar([(t(i)+t(i+1))/2,(t(i+1)+t(i))/2],[f(c),f(c)],t(i+1)-t(i)); endtext(1,1,[num2str(j),'個分割點']);M(:,(j-3)/20+1)=getframe;holdoff;endmovie(M,10)38選擇特殊的n值，將上面程序略加修改后，可得到如下圖形：積分和39圖中陰影部分為積分和，從圖中可看出在分割點為20時，陰影部分的上邊界還很粗糙；當(dāng)分割點為80時該邊界已比較光滑，若再增加分割點的個數(shù)；當(dāng)分割點為640時，已很難看出該邊界與曲線有什么地方不一樣了，這說明此時用積分和近似定積分產(chǎn)生的誤差已非常小。40當(dāng)定積分存在時，所有任取的積分和在

趨向于0時的極限都相同，此時可以選擇較簡單的劃分與簡單的

i，一般的做法是將區(qū)間等分，且讓小區(qū)間的某端點作為

i(i=1,2,…n),這樣積分和便成為或?qū)τ谶B續(xù)函數(shù)的定積分，用這兩個式子來近似計算是比較簡單的。在數(shù)值計算中，稱用上兩式近似定積分的方法為矩形法，這是因為這兩個式子在幾何上表示一些矩形面積的和。41(4)定積分近似計算的梯形法在近似計算中也常用式子近似計算定積分。該式是上兩式的的平均。它在幾何上表示一些梯形面積的和，故它被稱為近似計算定積分的梯形法。42例13用梯形法近似計算定積分解我們知道f(x)=ex在[0,1]上連續(xù),所以定積分存在?，F(xiàn)將區(qū)間[0,1]均分為等分，由梯形公式得：利用此式，編程如下：43f=inline('exp(x)');a=0;b=1;s0=1;s1=0;n=20;m=6;whileabs(s0-s1)>10^-ms1=s0;i=1:n;s0=sum((f(a+(i-1)*(b-a)/n)*(b-a)/n+f(a+i*(b-a)/n)*(b-a)/n)/2);n=n*2;endfprintf('%s%s%g\n','exp(x)','在[0,1]上的積分約為',s0)運行結(jié)果為：exp(x)在[0,1]上的積分約為1.71828。上面的數(shù)據(jù)1.71828除了整數(shù)位外，小數(shù)部分的前5位與的前5位完全一樣。44(5)定積分近似計算的MonteCarlo方法設(shè)函數(shù)f(x)定義在區(qū)間[a,b]上，當(dāng)axb時有0f(x)H，其中H是某個非負數(shù)（如圖）。今向圖中的矩形內(nèi)隨機投點，對于位于曲線y=f(x)下方的圖形，其面積的一種合理估計應(yīng)該是矩形的面積乘以落在該圖形內(nèi)的隨機點個數(shù)占總隨機點數(shù)的百分比，即：其中，為m隨機點總數(shù)，s是落在位于曲線y=f(x)下方的圖形中的隨機點個數(shù)。上面的結(jié)論意味著下式成立：45作為例子，我們用上面的方法來求定積分為此，編程如下：a=0;b=1;m=1000;s=0;H=exp(1);%s設(shè)置為落在曲邊梯形內(nèi)的點數(shù)fori=1:mxi=rand();yi=H*rand();ifyi<exp(xi)s=s+1;end%如果隨機點落在曲邊梯形內(nèi),s增加1endfprintf('%s%g\n','exp(x)在[0,1]上的積分約等于',H*(b-a)*s/m)exp(x)在[0,1]上的積分約等于1.75601。466.多元函數(shù)的極限與求導(dǎo)在MATLAB中多元函數(shù)求極限、求導(dǎo)仍分別由命令函數(shù)limit()、diff()來完成運算。還有命令函數(shù)Jacobian()求偏導(dǎo)。47clearsymsxy;f=(1-x*y)/(x^2+y^2);fx=limit(f,'x',0);fxy=limit(fx,'y',1)fxy=

148clearsymsxy;F=atan(y/x)-log(sqrt(x^2+y^2));diff(F,'x')ans=

-y/x^2/(1+y^2/x^2)-1/(x^2+y^2)*xans=

1/x/(1+y^2/x^2)-1/(x^2+y^2)*ydiff(F,'y')49此例使用Jacobian求偏導(dǎo)可能更方便。clearsymsxy;z=atan(y/x)-log(sqrt(x^2+y^2));dxy=jacobian(z)dxy=

[-y/x^2/(1+y^2/x^2)-1/(x^2+y^2)*x,1/x/(1+y^2/x^2)-1/(x^2+y^2)*y]50clearsymsxyz;F=x^2+y^2+z^2-4*zF=

x^2+y^2+z^2-4*zFx=diff(F,'x')Fx=

2*x51Fz=diff(F,'z')Fz=

2*z-4G=-Fx/FzG=

-2*x/(2*z-4)52在MATLAB的符號運算中，用命令Jacobian(f,v)來計算梯度。7.梯度計算與方向?qū)?shù)例17已知函數(shù)求梯度gradf（(1,-1,2)）。clearsymsxyzf=x^2+y^2+z^2;gr=jacobian(f)gf=subs(subs(subs(gr,1),-1),2)gf=2-2453數(shù)值運算中，我們?nèi)匀豢梢杂妹頳iff()對數(shù)據(jù)矩陣進行梯度向量的計算，但是我們常用的命令gradient(),其格式如下：Fx=gradient(F)[Fx,Fy]=gradient(F)[Fx,Fy,Fz]=gradient(F)分別用來求一維、二維與多維向量F的數(shù)值梯度。54clearsymsxyz=x^2+y^2;dxy=jacobian(z)dxy=

[2*x,2*y]55clearx=[-1:0.1:1];[X,Y]=meshgrid(x);%生成三維圖形的x,y坐標矩陣Z=X.^2+Y.^2;%計算z坐標矩陣[DX,DY]=gradient(Z);%y計算矩陣Z的梯度向量subplot(2,1,1);%生成第一個繪圖區(qū)間mesh(Z);%繪制z=x2+y2的三維網(wǎng)線圖subplot(2,1,2);%生成第二個繪圖區(qū)間Quiver(DX,DY);%繪出矩陣Z的二維向量場56方向?qū)?shù)顯示了函數(shù)在某點沿著某個方向上的變化率，而模取得最大值的方向?qū)?shù)就是該點的梯度。函數(shù)在某點沿著某個方向上的方向?qū)?shù)等于在這點的梯度和沿該方向的單位向量的數(shù)量積。設(shè)某個方向上的單位向量為v,則方向?qū)?shù)為：Jacobian(f)?v57clearsymsxy;z=x*exp(2*y);v=[sqrt(2)/2,-sqrt(2)/2]；dxy=jacobian(z)；%計算梯度dxy=

[exp(2*y),2*x*exp(2*y)]58T=dot(dxy,v)%計算數(shù)量積T=

1/2*exp(2*conj(y))*2^(1/2)-conj(x*exp(2*y))*2^(1/2)subs(subs(T,'x',1),'y',0)ans=-0.7071598.多元函數(shù)的極值要使用MATLAB解決多元函數(shù)的極值問題，首先需要明確“充分條件定理”。設(shè)z=f(x,y)在點（x0,y0）的某個鄰域內(nèi)連續(xù)且有一階二階連續(xù)偏導(dǎo)數(shù)，又fx(x0,y0)=0,fy(x0,y0)=0,令A(yù)=fxx(x0,y0),B=fxy(x0,y0),C=fyy(x0,y0),則有：。（1）當(dāng)AC-B2>0時有極值，且A<0時有極大值，A>0時有極小值（2）當(dāng)AC-B2<0時沒有極值（3）當(dāng)AC-B2=0時極值不確定，需要具體討論。60clearsymsxyz;f=x^3-y^3+3*x^2+3*y^2-9*x;F=jacobian(f);%計算梯度F=

[3*x^2+6*x-9,-3*y^2+6*y][X,Y]=solve(F(1),x,F(2),y);%求駐點61X=

1-31-3

Y=

0022[X,Y]=solve(F(1),x,F(2),y);%求駐點62dxx=diff(F(1));%計算Fxx(x,y)dxx=

6*x+6dyy=diff(F(2));%計算Fyy(x,y)dyy=

-6*y+6dxy=diff(F(1),y);%計算Fxy(x,y)dxy=

063例題分析：一共四個駐點：(-3,0),(-3,2),(1,0),(1,2)（1）在駐點(-3,0)處，A=-12，B=0，C=6，AC-B2=-72<0，所以此點不是極值點。（2）在駐點(-3,2)處，A=-12，B=0，C=-6，AC-B2=72<0，所以此點是極大值點。（3）在駐點(1,0)處，A=12，B=0，C=6，AC-B2=72<0，所以此點是極小值點。（4）在駐點(1,2)處，A=12，B=0，C=-6，AC-B2=-72<0，所以此點不是極值點。649.重積分65clearsymsxyz;f=z;A1=int(int(int(f,z,x^2+y^2,2),y,0,sqrt(2-x^2)),x,0,sqrt(2));A1=

2/3*piA2=int(int(int(f,z,x^2+y^2,1),y,0,sqrt(1-x^2)),x,0,1);A2=

1/12*piA1-A2ans=

7/12*pi6610.曲線積分與曲面積分67clearsymsabxyt;x=a*cos(t);y=b*sin(t);dx=diff(x,t);dy=diff(y,t);dx=

-a*sin(t)dy=

b*cos(t)68ds=sqrt(dx^2+dy^2);f=x*y;ds=

(a^2*sin(t)^2+b^2*cos(t)^2)^(1/2)I=int(f*ds,t,0,pi/2)f=

a*cos(t)*b*sin(t)I=

1/3*a*b*((a^2)^(1/2)*a^2-(b^2)^(1/2)*b^2)/(a^2-b^2)

6970symsthr;symsapositive;y=(a*r)/(a^2-r^2);I=int(int(y,r,0,sqrt(a^2-(a/2)^2)),th,0,2*pi)I=

2*a*log(2)*pi71clearsymsxyz;z=sqrt(1-x^2-y^2);f=x*y*z;z=

(1-x^2-y^2)^(1/2)f=

x*y*(1-x^2-y^2)^(1/2)72I=int(int(f,y,0,sqrt(1-x^2)),x,0,1)I=

1/1573clearsymsaxy;y=sqrt(a*x-x^2);dyx=diff(y,x);f=x^2+y^2+4*x*y*dyx;dyx=

1/2/(a*x-x^2)^(1/2)*(a-2*x)f=

a*x+2*x*(a-2*x)74I=int(f,x,0,a)I=

1/6*a^37511.數(shù)項級數(shù)的求和與收斂展開MATLAB符號工具箱給我們提供了功能強大的級數(shù)求和的命令函數(shù)symsum(),其格式如下：symsum(function,vaiaale,a,b)Function為級數(shù)的通項表達式,vaiaale用來聲明通項中的求和變量,a和b分別為求和變量的起始點與終止點。76clearsymsn;f1=(2*n-1)/2^n;f2=1/(n*(2*n+1));I1=symsum(f1,n,1,inf);I2=symsum(f2,n,1,inf)I1=

3I2=

2-2*log(2)77clearsymsnx;f1=1/n;f2=sqrt(n+1)-sqrt(n);I2=symsum(f2,n,1,inf);MATLAB計算不出結(jié)果，因此無法斷定其斂散性。I1=

InfI1=symsum(f1,n,1,inf);I2=

sum((n+1)^(1/2)-n^(1/2),n=1..Inf)級數(shù)發(fā)散。78clearsymsnx;f1=sin(x)/n^2;f2=(-1)^(n-1)*x^n/n;f3=1/n;I2=symsum(f2,n,1,inf);I2=

log(x+1)I1=

1/6*sin(x)*pi^2I1=symsum(f1,n,1,inf);79然而，正如我們看到的那樣，用symsum語句得不到冪級數(shù)的收斂半徑和收斂區(qū)間。要求一個冪級數(shù)的收斂半徑和收斂區(qū)間，我們可以象手算那樣根據(jù)有關(guān)結(jié)論去求，也可以用計算機做些數(shù)值實驗來判斷。8012.泰勒級數(shù)展開對冪級數(shù)來說，除了要會求冪級數(shù)的和函數(shù)外，還存在著將函數(shù)展開為冪級數(shù)的問題，其目的就是求一個冪級數(shù)收斂到已知的函數(shù)，這等價于用一系列多項式逼近已知函數(shù)。在MATLAB中泰勒級數(shù)展開由命令函數(shù)taylor(),其格式如下：taylor(function,n,vaiaale,a)Function為待展開函數(shù)表達式,n為展開階數(shù)，缺省是6階。Vaiaale為聲明function中的變量,a為變量求導(dǎo)的取值點，缺省為0,即麥克勞林級數(shù)展開。81clearsymsx;f=1/(1+x^2);taylor(f);taylor(f,20)%20階麥克勞林級數(shù)ans=

1-x^2+x^4-x^6+x^8-x^10+x^12-x^14+x^16-x^18ans=

1-x^2+x^482clearsymsx;f=1/(x^2+4*x+3);taylor(f,10,x,1);ans=

7/32-3/32*x+7/128*(x-1)^2-15/512*(x-1)^3+31/2048*(x-1)^4-63/8192*(x-1)^5+127/32768*(x-1)^6-255/131072*(x-1)^7+511/524288*(x-1)^8-1023/2097152*(x-1)^983pretty(ans)

21533147/32-3/32x+7/128(x-1)----(x-1)+----(x-1)5122048635127625575118-----(x-1)+-----(x-1)-------(x-1)+------(x-1)81923276813107252428810239--------(x-1)2097152pretty(ans)84下面我們從圖形上來觀察多項式逼近函數(shù)的過程，以函數(shù)f(x)=sinx為例：13.函數(shù)逼近symsx;f=inline('sin(x)');g=0;fori=0:10s=diff(f(x),x,i);s1=subs(s,x,0);g=g+s1*x^i/prod(1:i);ifi==1|i==4|i==7|i==10subplot(2,2,1+(i-1)/3);ezplot(f(x),[0,2*pi]);holdon;85ezplot(g,[0,2*pi]);axis([0,6,-4,4]);title(i);endendholdoff86圖中變動的曲線依次為的1、5、8、10階Taylor多項式的圖形，我們看到：階數(shù)越高，Taylor多項式與其函數(shù)之間的差異越小。對于連續(xù)函數(shù)f(x)，在某個閉區(qū)間內(nèi)，能不能找到一個近似程度既好、階數(shù)又不是很高的多項式呢？在這里，我們作一些討論：如何尋找函數(shù)的最佳一次逼近多項式。在這里，我們作一些討論：如何尋找函數(shù)的最佳一次逼近多項式。這說明，在Taylor級數(shù)的收斂區(qū)間內(nèi)，要取到函數(shù)的比較滿意的近似值，選取的部分和階數(shù)必須足夠的高。87最佳一次逼近多項式用多項式逼近函數(shù)，當(dāng)次數(shù)固定時，存在著最佳逼近多項式，所謂最佳，就是產(chǎn)生的誤差最小，而誤差通常是指最大偏差。具體如下：設(shè)f(x)C(a,b)(C(a,b)表示[a,b]區(qū)間上連續(xù)函數(shù)的集合）,則在[a,b]上以多項式Pn(x)代替函數(shù)f(x)產(chǎn)生的最大偏差為。在上面這種意義下的最佳逼近多項式稱為最佳一致逼近多項式，一般而言，Taylor多項式非最佳一致逼近多項式。88在理論上，若函數(shù)f(x)在閉區(qū)間[a,b]上連續(xù)，則在該區(qū)間內(nèi)存在唯一一個次最佳一致逼近多項式，但是，要求出的次最佳一致逼近多項式卻非常困難．下面我們僅討論當(dāng)f''(x)在(a,b)內(nèi)不變號時，如何確定其最佳一次逼近多項式。我們知道，當(dāng)f''(x)在(a,b)內(nèi)時，在[a,b]上的圖形是一段凹弧或一段凸弧，不妨假設(shè)其為凹弧。現(xiàn)連接兩個端點M、N，得線段MN，它所在的直線方程為：89因弧段y=f(x)在[a,b])上光滑，而f(x)的二階導(dǎo)數(shù)不變號，故在該弧段內(nèi)存在唯一的一點C(x2,f(x2)，此點到線段MN的距離最遠。所以，此x2一定也是函數(shù)在[a,b]上的最大值點，因而有這說明弧段在C處的切線平行于MN。90今讓直線M1N1平行MN及l(fā),并位于MN與l的正中間，則該直線位于區(qū)間[a,b]內(nèi)的部分M1N1就是所要求的最佳一次逼近多項式函數(shù)的圖形。對于位于[a,b]區(qū)間內(nèi)的任意其他的直線段M2N2，下面就M2N2的不同位置來說明上述結(jié)論。911.若M2N2平行M1N1(如圖)，顯然有最佳一次逼近(之一)92最佳一次逼近（之二）2.若M2N2與M1N1相交(如圖6.4)，且交點在(a,b)內(nèi)，則其兩個端點M2、N2中，總有一個位于對應(yīng)端點M1,N1的下方，故也有93最佳一次逼近（之三）3.若M2N2與M1N1相交(如圖6.4)，且交點不在(a,b)內(nèi)，可直接看出上式是成立的。94現(xiàn)在，我們來求M1N1所在的直線方程，由前面的討論知，直線的方程為而MN的方程(6.5)即為95因M2N2平行M1N1與l(如圖)，且位于這兩條平行直線的正中間，故它的方程為將此式整理，得f(x)在[a,b]上的最佳一次逼近多項式為：我們還可進一步求出，用最佳一次逼近多項式P1(x)近似函數(shù)f(x)所產(chǎn)生的最大偏差，