實(shí)驗(yàn)六常微分方程

函數(shù)表達(dá)Matlab表達(dá)函數(shù)表達(dá)Matlab表達(dá)

y=sinxy=sin(x)y=exy=exp(x)

y=cosxy=cos(x)y=Inxy=iog(x)

y=x2y=xA2y=arctanxy=atan(x)

sinx

y=x2sinxy=x2*sin(x)y=xy=sin(x)/x

實(shí)驗(yàn)六常微分方程

dsolve]dx=-a*x')%輸出ans=Cl*exp(-(2*x)

1=/%{'小=一〃*工；比(0)=1'，'5)%輸出x=exp(-〃*s)

w=dsolve(D3w=-w7w(0)=1,Dw(0)=0,D2w(0)=0,)%輸出略

[7,g]=^Zve('zy=/+g,Dg=-/+g','/(o)=r,'g(o)=2')%輸出略

一、實(shí)驗(yàn)內(nèi)容

求微分方程的解析解

2,

例1求微分方程y+2.二比一"的通解。

輸入：

dsolve(*dy+2*x*y=x*exp(-xA2)',*x*)

便得到微分方程的通解：

(l/2*xA2+Cl)*exp(-x-'2)

其中Cl是任意常數(shù)。

例2求微分方程=0在初始條件y(l)=2e下的特解。

輸入：

dsolve(*x*dy+y-exp(-x)=0*z*y(1)=2*exp(1)*z*x*)

便得到微分方程的特解：

(-exp(-x)+exp(-1)+2*exp(1))/x

實(shí)驗(yàn)七空間圖形的畫法

(一)、三維曲線的繪制

(1)plot3命令

繪制三維曲線，基本使用形式是：

plot3(x,y,z,'s')

例如，一條空間螺旋線的參數(shù)方程是：x=cosf,y=sinf,z=看(0W8%)。

輸入：

t=0:0.1:8*pi;

x=cos(t);

y=sin(t);

z=t/10;

plot3(x,y,z)

xlabel(1x1);

ylabel(1y1);

zlabel(1z1);

則輸出了一條螺旋線(見圖6-1)o

圖6-1

(2)ezplot3命令

同繪制二維曲線相似，也有簡潔的繪制空間曲線命令ez〃"3,具體使用方法

同相似，形式是：

ezpI3o(fx(t)*,*y(t)*,*z

其中x(t),y(t),z(t)是曲線的參數(shù)方程的表達(dá)式。tl,t2是作圖時(shí)參數(shù)t的范圍。

例如上述螺線也可輸入下面命令得到：

/p，“3('cos(t)；'sin(t)7t/10\[0,8*pi\)

輸出圖形略。

(二)三維曲面網(wǎng)線圖與曲面圖的繪制

Meshgrid：此函數(shù)用來產(chǎn)生三維繪圖時(shí)的陣列(網(wǎng)格).

例：

x=l:3;

y=10:14;

[X,Y]=meshgrid(x,y)

X=

123

123

123

123

123

101010

111111

121212

131313

141414

mesh:繪制3d網(wǎng)狀立體圖

surf:繪制3d彩色表面圖

title:標(biāo)識(shí)文字，title('text',…)

xlabel：對(duì)x軸的坐標(biāo)

ylabel：對(duì)y軸的坐標(biāo)

gridon:在圖表中交網(wǎng)格線

holdon:保持目前的圖表，后續(xù)可附加

subplot(m,n,p)：多個(gè)窗口同時(shí)顯示，子窗口圖形個(gè)數(shù)為m*n矩陣

MATLAB軟件繪制曲面圖基本形式是：

(1)[X,Y]=meshgrid(x,y)

(2)Z=/(X,Y)

(3)mesh(X,Y,Z)%繪制網(wǎng)線圖

(4)5wrf(X,Y,Z)%繪制曲面圖

例如，畫出曲面Z=%2+y2的圖形。輸入：

x=-2:0.1:2;

y=-2:0.1:2;

[X,Y]=meshgrid(x,y);

Z=X.A2+Y.A2;

Mesh(X,Y,Z)

surf(X,Y,Z)

得到曲面Z=x2+y2(見圖6-2)。

圖6-2

又如，參數(shù)方程：x=2sin℃ose,y=2sin°sine,z=2cos0是以原點(diǎn)為圓心,

半徑為2的球面，其中0<0<萬萬，因此只要輸入：

t=0:0.1:pi;

r=0:0.1:2*pi;

[R,T]=meshgrid(r,t);

x=2.*sin(T).*cos(R);

y=2.*sin(T).*sin(R);

z=2,*cos(T);

surf(x,y,z)

便作出了方程為4+y2+z2=4的球面(見圖6?3方

圖6-3

繪制二元函數(shù)圖形也可用簡捷繪制的ezs〃^命令，它的使用格式為:

ezsurf(f(x,y),[a,b,u,v])

即可繪制函數(shù)/(x,y)在區(qū)域［a,勿上的圖形。當(dāng)省略區(qū)域時(shí)，默認(rèn)區(qū)域

是[-2肛21]x[-2乃,21]o

一、實(shí)驗(yàn)內(nèi)容

(一)空間曲線

例1作出空間曲線x=/cos/,y=/sin/,z=2/(0W/W6萬)的圖形。

輸入：

ezplot3('t*cos(t)'；t*sin(t)','2*t',[0,6*pi])

輸出如圖6-4所示。輸出如圖6-4所示。

（二）空間二次曲面

2,2

yz，

例2作橢球面工+—+—=1的圖形。

491

該曲面的參數(shù)方程是x=2sin^>cos0,y=3sin^sin0,z=cos(p,其中

G<(p<7i.O<0<o輸入：

t=0:0.1:pi;

r=0:0.1:2*pi;

[R,T]=meshgrid(r,t);

x=2,*sin(T),*cos(R);

y=3.*sin(T).*sin(R);

z=cos(T);

surf(x,y,z)

輸出如圖6-5所示。

圖6-5

222

例3作單葉雙曲面三+匕—二=1的圖形。

149

該曲面的參數(shù)方程是％=se&?svny=,2?secvo,其中

jrjr

——<w<—,Q^v<TQ.O輸入：

22

t=-pi/4:0.1:pi/4;

r=0:0.1:2*pi;

[R,T]=meshgrid(r,t);

x=sin(R).*sec(T);

y=cos(R).*sec(T);

z=3*tan(T);

surf(x,y,z)

輸出如圖6-6所示。

圖6-6

222

例4作雙葉雙曲面三+三一三=一1的圖形。

該曲面的參數(shù)方程是x=1.5cotw-cosv,=1.4cotw-sinv,z=1.3cscw,其中參數(shù)

rr.TT

0<M?—乃4vW萬對(duì)應(yīng)雙葉雙曲面的一葉，參數(shù)——〈”<0,—乃乃對(duì)應(yīng)雙

22

葉雙曲面的另一葉。輸入：

t=pi/1000:0.1:pi/2;

r=-pi:0.1:pi;

[R,T]=meshgrid(r,t);

xl=l.5.*cos(R),*cot(T);

yl=l.4.*sin(R).*cot(T);

zl=l.3.*csc(T);

mesh(xl,yl,zl)

holdon

mesh(-xlr-yl,-zl)

輸出如圖6-7所示。

例3求微分方程產(chǎn)2)/+5y=e，cos2x的通解。(simplify)化簡結(jié)果

輸入：

simplify(dsolve(*D2y-2dy+5*y=exp(x)*cos(2*x)*x*))

便得到微分方程的通解：

sin(5A(1/2)*x)*C2+cos(5A(1/2)*x)*Cl+2/5+l/10*exp(x)*cos(2*x)

+1/5*exp(x)*sin(2*x)

其中Cl,C2是任意常數(shù)。

解微分方程組的命令格式為：

[x,y]=dsolve^dx=/(x,y)'6=g(x,y)1)

dxct

---bx+2y=e黑))::下的特解。

例4求微分方程組出在初始條件＜

電-x-y=0

〔dt

輸入:

[xzy]=dsolve('dx=-x-2*y+exp(t)','dy=x+y','x(0)=l','y(0)=0')

便得到微分方程的特解:

x=cos(t)

y=l/2*sin(t)-1/2*cos(t)+1/2*exp(t)

實(shí)驗(yàn)八

(-)求偏導(dǎo)數(shù)命令

命令d療既可以用于求一元函數(shù)的導(dǎo)數(shù)，也可以用于求多元函數(shù)的偏導(dǎo)數(shù)。

用于求偏導(dǎo)數(shù)時(shí)，可根據(jù)需要分別采用如下幾種形式：

y,z),x)%求/(%,y/)對(duì)刀的偏導(dǎo)數(shù)

或「(/(X,Xz),y)%求于(x,y,z)對(duì)y的偏導(dǎo)數(shù)

d療(/(x,y,z),無,2)或diff(diff(f(x,y,z),x),x)%求/(x,y,z)對(duì)％的

二階偏導(dǎo)數(shù)

diff{diff(f(x,y,z),x),y%求/(%y,z)對(duì)x、y的混合偏導(dǎo)數(shù)

其余類推。

一、實(shí)驗(yàn)內(nèi)容

(-)求多元函數(shù)的偏導(dǎo)數(shù)與全微分

例1iSz=sin(xy)+cos2(xy),求當(dāng)

oxoyoxoxoy

輸入：

symsxy

z=*sin(x*y+(cos(x*y))A2)*;

diff(z,x)

diff(z,y)

diff(z,x,2)

diff(diff(z,x),y)

便依次得到函數(shù)表達(dá)式及所求的四個(gè)偏導(dǎo)數(shù)的結(jié)果：

ans=cos(cos(x*y)A2+x*y)*(y-2*y*cos(x*y)*sin(x*y))

ans=cos(cos(x*y)A2+x*y)*(x-2*x*cos(x*y)*sin(x*y))

ans=-sin(cos(x*y)A2+x*y)*(y-2*y*cos(x*y)*sin(x*y))A2-

cos(cos(x*y)A2+x*y)*(2*yA2*cos(x*y)A2-2*yA2*sin(x*y)A2)

ans=-cos(cos(x*y)A2+x*y)*(2*x*y*cos(x*y)A2+

2*cos(x*y)*sin(x*y)-2*x*y*sin(x*y)A2-1)-sin(cos(x*y)A2+x*y)*(x

-2*x*cos(x*y)*sin(x*y))*(y-2*y*cos(x*y)*sin(x*y))

例2設(shè)z=(a+砂)"求二，一。

dxdy

輸入：

symsxy

z=*(a+x*y)Ay*;

diff(z,x)

diff(z,y)

則有輸出：

ans=yA2*(a+x*y)A(y-1)

ans=log(a+x*y)*(a+x*y)Ay+x*y*(a+x*y)A(y-1)

(二)多元函數(shù)的極值

例3求函數(shù)/(%?)=%3-/+3_￥2+3,2一9%的極值。

輸入：

symsxy

f=,xA3-yA3+3*xA2+3*yA2-9*x,;

fx=diff(f,x)

fy=diff(f,y)

輸出為：

fx=3*xA2+6*x-9

fy=6*y-3*yA2

再輸入：

[x0,y0]=solve('3*x/'2+6*x-9=0','-3*y/'2+6*y','x','y')

輸出為如下四個(gè)駐點(diǎn)：

xO=

1

-3

1

-3

yO=

0

0

2

2

再輸入：

fxx=diff(f,x,2);

fyy=diff(f,y,2);

fxy=diff(fx,y);

a=(fxx)*(fyy)-(fxy)A2;

x=-3;

y=0；

al=eval(a)

bl=eval(fxx)

cl=eval(f)

x=-3;

y=2;

a2=eval(a)

b2=eval(fxx)

c2=eval(f)

x=l;

y=0;

a3=eval(a)

b3=eval(fxx)

c3=eval(f)

x=l;

y=2;

a4=eval(a)

b4=eval(fxx)

c4=eval(f)

我們得到了四個(gè)駐點(diǎn)處的判別式函數(shù)A(x,y)=%4-抬，兒與/的值。歸

納以后用表格形式列出。

y九A(x,y)f

-30-12-7227

-32-127231

101272-5

1212-72-1

從中可見：

%=_3,y=2時(shí)，判別式A(x,y)=72,九=-12,因此函數(shù)有極大值31.

x=l,y=O時(shí)，判別式4(x,y)=72,九=12,因此函數(shù)有極小值-5.

龍=—3,y=0和無=l,y=2時(shí)，判別式A(x,y)=—72,函數(shù)在這些點(diǎn)不取極值.

二、實(shí)驗(yàn)作業(yè)

“、兒2_p.&dz

1.設(shè)z=e”,求——,—o

dxdy

-(x2+j;2)//822

2.z=^(cosx+siny),求包2,。％o

dxdydxdy

3.求f(x,y)=-120x3-30x4+18x5+5x6+30xy2的極值。

實(shí)驗(yàn)八

實(shí)驗(yàn)?zāi)康?/p>

掌握用MATLAB計(jì)算二重積分的方法。提高應(yīng)用重積分解決各種應(yīng)用問題的

能力。

學(xué)習(xí)MATLAB命令

(-)重積分命令

命令int也可用于計(jì)算重積分，常用格式如下：

int(int(/(x,y),y,%(x),y2(x),x,a,b)

用來求二重積分力力y,其中。是X型區(qū)域aWxWA,%(x)WyW%(x)。

D

注：y型區(qū)域上的二重積分的格式也是類似的。

例如要計(jì)算xy2dydx,

輸入：

symsxy

int(int(x*yA2,0,x),x,O,l)

得到：

ans=1/15

二、實(shí)驗(yàn)內(nèi)容

(一)計(jì)算重積分

例1計(jì)算口劃其中。為由x+y=2,x=6，y=2所圍成的有界

D

區(qū)域。

先做出區(qū)域。的草圖，手工就可以確定積分限。應(yīng)先對(duì)X積分，

輸入：

symsxy

A

int(int(x*y2zx,2-y,sqrt(y)),y,l,2)

輸出為：

ans=193/120

例2計(jì)算口"(,+廿)辦dy,其中。為Y+V41。

D

如果用直角坐標(biāo)計(jì)算，輸入：

symsxyreal

f=exp(-(xA2+yA2));

int(int(f,y,-sqrt(l-xA2),sqrt(l-xA2)),x,-l,l)

輸出為:

ans=int((piA(1/2)*erf((1-xA2)A(1/2)))/exp(xA2),x=-1..1)

積分遇到了困難。改用極坐標(biāo)，也用手工確定積分限，輸入:

symsrs

f=exp(-(r人2))*r;

int(int(f,r,0,1),s,0,2*pi)

輸出為:

ans=-pi*(1/exp(1)-1)

（二）重積分的應(yīng)用

例3求曲面z=4-x?-y2在My平面上的面積s。

輸入：

symsxyz

ezsurf(14-xA2-yA21)

輸出如圖8-1所示。

圖8-1

觀察曲面的圖形，可見是一個(gè)旋轉(zhuǎn)拋物面。計(jì)算曲面面積公式是

S=JJyjl+z[+z^dxdy。

D

輸入：

symsxy

z=*4-xA2-yA2,;

f=sqrt(1+diff(z,x)A2+diff(z,y)人2)

輸出為：

f=(4*x-2+4*yb+1)A(1/2)

因此用極坐標(biāo)計(jì)算。

輸入：

symsrt

f=*sqrt(l+4*rA2)*r*;

s=int(int(f,r,0,2),t,0,2*pi)

輸出為：

s=(pi*(17*17A(1/2)-1))/6

三、實(shí)驗(yàn)作業(yè)

1.計(jì)算［%jj（ysin九一xsin丁）為協(xié)。

2.求積分Jsin（exy）dxdy的近似值。

用極坐標(biāo)計(jì)算￡￡2：2dydxo

實(shí)驗(yàn)九二重積分

（二）重積分命令

命令int也可用于計(jì)算重積分，常用格式如下:

int(int(/(x,y),y,%(x),y2(x),x,a,b)

用來求二重積分U/(%,y)dxdy,其中。是X型區(qū)域。4%<九%(%)<%(》)。

D

注：F型區(qū)域上的二重積分的格式也是類似的。

例如要計(jì)算xy2dydx,

symsxy

int(int(x*yA2,y,0,x),x,O,l)

得到：ans=1/15

(三)二元函數(shù)的數(shù)值積分

函數(shù)加的功能是求矩形區(qū)域上的二元函數(shù)的數(shù)值積分。其格式如下:

q=dblquad{jun,xmin,xmax,ymin,ymax,tol)

最后一個(gè)參數(shù)的意義是用指定的精度⑹代替默認(rèn)精度KF',再進(jìn)行計(jì)算。

例如輸入：

f=inline(*sqrt(max(1-(x.A2+y.A2),0)),);

Q=dblquad(f,-1,1,-lf1)

得到：

Q=

2.0944

四、實(shí)驗(yàn)內(nèi)容

計(jì)算重積分

例1計(jì)算孫2dmy,其中。為由x+y=2,x=6，y=2所圍成的有界區(qū)域。

D

先做出區(qū)域。的草圖，手工就可以確定積分限。應(yīng)先對(duì)X積分，輸入：

symsxy

A

int(int(x*y2zx,2-y,sqrt(y)),y,l,2)

輸出為：

ans=193/120

例2計(jì)算JJ如力，其中。為好+J?<i。

D

如果用直角坐標(biāo)計(jì)算，輸入:

symsxyreal

f=exp(-(xA2+yA2));

int(int(f,y,-sqrt(l-xA2),sqrt(l-xA2))

輸出為：

ans=int((piA(1/2)*erf((1-xA2)A(1/2)))/exp(xA2),x=-1..1)

積分遇到了困難。

(選學(xué)內(nèi)容一一極坐標(biāo))解決方法:

改用極坐標(biāo),也用手工確定積分限,輸入:

symsrs

f=exp(-(r^2))★r;

int(int(ffrfOfl)fs,0,2*pi)

輸出為:

ans=-pi★(1/exp(1)-1)

實(shí)驗(yàn)十無窮級(jí)數(shù)

符號(hào)表達(dá)式求和函數(shù)

symsum{S{ky)%返回符號(hào)表達(dá)式S(左)中的符號(hào)變量左從0到左-1的和值

symsum(S(k),v)%返回符號(hào)表達(dá)式S(左)中指定的符號(hào)變量左由v代替，再

對(duì)v從0至―-1的和值

symsum(S(k),a,b)%返回符號(hào)表達(dá)式S(左)中的符號(hào)變量左從a到。的和

值

symsum(S(^),v,a,b)%返回符號(hào)表達(dá)式S伏)中的符號(hào)變量左由v代替，再

對(duì)從。到。的和值

(二)符號(hào)函數(shù)的泰勒級(jí)數(shù)展開式函數(shù)

taylor(f(x))或taylor(y,x)%求函數(shù)y=/(x)的5階麥克勞林展開式

taylor(f(x),n)%求函數(shù)/(%)的〃-1階麥克勞林展開式

taylor(f(x),n,a)%求函數(shù)/(%)在％=a處的〃-1階泰勒展開式

(三)在符號(hào)表達(dá)式中進(jìn)行符號(hào)替換的函數(shù)

subs(S,old,new)%將符號(hào)表達(dá)式S中的符號(hào)變量old用new代替

(四)符號(hào)表達(dá)式的化簡函數(shù)

simple(expr)或simplify(expr)%用于化簡符號(hào)表達(dá)式expr

一、實(shí)驗(yàn)內(nèi)容

(-)級(jí)數(shù)求和

當(dāng)符號(hào)變量的和存在時(shí)，可以用symw/m命令來求無窮級(jí)數(shù)的和。

001

例1求—。

念+8〃+3

輸入：

symsn

sl=symsum(1/(4*nA2+8*n+3),n,1,inf)

得到該級(jí)數(shù)的和為：

si=1/6

(二)求寨級(jí)數(shù)的收斂域

例2求￡4"'白:3)"的收斂域與和函數(shù)。

〃=0〃+1

輸入：

symsnx;

al=4A(2*n)*(x-3)An/(n+1);

a2=subs(al,n,n+1);

a=simplify(a2/al)

p=limit(a,n,inf)

輸出為:

p=16*x-48

注意，這里對(duì)的和田都沒有加絕對(duì)值。因此上式的絕對(duì)值小于1時(shí)，幕級(jí)數(shù)

收斂，大于1時(shí)發(fā)散。求出收斂區(qū)間的端點(diǎn)，輸入：

Xl=solve(*48-x*l6=11)

X2=solve(,x*16-48=l,)

輸出為：

XI=47/16

X2=49/16

由止匕可知竺<x(竺時(shí)收斂，X<”或x>竺時(shí)發(fā)散。為了判斷端點(diǎn)的斂散

16161616

性，輸入：

simplify(subs(al,'x',49/16))

得到X為右端點(diǎn)時(shí)塞級(jí)數(shù)的一般項(xiàng)為：

ans=1/(n+1)

由于級(jí)數(shù)￡一~7發(fā)散，故當(dāng)冗=”時(shí)發(fā)散。再輸入：

+l16

simplify(subs(al,*x1,47/16))

輸出為：

ans=(-1)An/(n+1)

由于級(jí)數(shù)El收斂，因此當(dāng)X=%時(shí)，級(jí)數(shù)收斂。

也可以在收斂域內(nèi)求得這個(gè)級(jí)數(shù)的和函數(shù)。輸入：

symsnx

s3=symsum(4A(2*n)*(x-3)An/(n+l),n,0,inf)

輸出為：

piecewise([49/16<=x.Inf]z[x<>49/16andabs(16*x-48)<=

1,-log(49-16*x)/(16*x-48)])

(三)將函數(shù)展開為嘉級(jí)數(shù)

例4求arctanx的5階麥克勞林