布仁滿都拉
(赤峰學(xué)院數(shù)學(xué)與計(jì)算機(jī)科學(xué)學(xué)院,赤峰024000)
在物理、化學(xué)、金融和醫(yī)學(xué)等很多領(lǐng)域中經(jīng)常建立數(shù)學(xué)模型研究相關(guān)問(wèn)題,其中常微分方程是經(jīng)常遇到的微分方程。因此,常微分方程是理科、工科和文科的很多專(zhuān)業(yè)的必修課。結(jié)合常微分方程和MATLAB 解算的教學(xué)方法比傳統(tǒng)的理論講授教法,更能提高學(xué)生的學(xué)習(xí)興趣和知識(shí)的混合應(yīng)用能力,下面用實(shí)例介紹常微分方程和MATLAB 的結(jié)合用法。
求方程y''-6y'+9y=4e3x通解。
解:我們先理論上算出方程組的解析解,然后給出它的MATLAB 解算。先求齊次方程的通解,特征方程為:
特征根為:
λ1,2=3
因此,齊次方程的通解為:
y=C1e3x+C2xe3x
由于3 是二重特征根,非齊次方程有形如:
y*=Ax2e3x
的特解。將它代入非齊次方程,比較x 的同次冪系數(shù),得A=2,所以:
y*=2x2e3x
所求通解為:
y=2x2e3x+C1e3x+C2xe3x
下面用MATLAB 符號(hào)法求微分方程的解。
鍵入:
y=dsolve('D2y-6*Dy+9*y=4*exp(3*x)','x')
回車(chē)得出:
y=2*x^2*exp(3*x)+C5*exp(3*x)+C6*x*exp(3*x)
求解方程組
解:我們先理論上算出方程組的解析解,然后給出它的MATLAB 解算。
系數(shù)矩陣為:
特征方程為:
特征根為:λ1=3,λ1=2,λ1=1。
先求λ1=3 對(duì)應(yīng)的特征向量:
a,b,c 滿足方程組:
即:
可得a=c,b=0,取一組非零解,例如令c=1,就有a=1,即:
下面求λ2=2 對(duì)應(yīng)的特征向量:
a,b,c 滿足方程組:
即:
可得a=c,b=c,取一組非零解,例如令c=1,就有a=1,b=1,即:
最后求λ3=1對(duì)應(yīng)的特征向量:
a,b,c 滿足方程組:
即:
可得a=0,b=c,取一組非零解,例如令c=1,就有b=1,即:
故方程組的通解是:
下面用MATLAB 符號(hào)求上述方程組的解析解:
鍵入:
回車(chē)得出:
在很多問(wèn)題中遇到的常微分方程的解析解是很難算出來(lái)的,這時(shí),我們可以用數(shù)值方法求近似解。
求解2y'+4xy=1,y(0)=0。
解:先用分離變量法和常數(shù)變易試算齊次方程的通解:
2y'+4xy=0
分離變量,得:
兩端積分,得:
解出y,得:
y=Ce-x2
由常數(shù)變易法,令:
y=C(x)e-x2
為非齊次方程的解,代入后得:
由于ex2的原函數(shù)不是初等函數(shù),積分:
的計(jì)算無(wú)法進(jìn)行。
下面用MATLAB 求該微分方程的數(shù)值解。
先建立如下M-函數(shù)文件
function y1=CDE(x,y)
y1=0.5-2*x*y;
在指令窗口中鍵入:
[x,y]=ode23(@CDE,[0 3],[0])
回車(chē)得出:
鍵入:
回車(chē)得出:
因此,在區(qū)間[0,3]上把x 分成了21 節(jié)點(diǎn),對(duì)應(yīng)地得出y 的21 個(gè)取值,即得到了微分方程的數(shù)值解。
進(jìn)一步,鍵入:
plot(x,y)
回車(chē)畫(huà)出數(shù)值解的圖形,見(jiàn)圖1。
圖1 數(shù)值解的圖形