12.6 一階微分方程之解法

12.6一階微分方程之解法

連續性函數之微分可以利用多種指令為之，而基本之細段微分法則是最早使用的概念，亦可利用MATLAB運算指令直接求解。

尤拉法(Euler method)

尤拉法為數值積分法最基本之一種。它是利用時間細段微分，利用初值及一階微分方程進行求解。以下列方程式為例：

 y'=dy/dt = f(t,y)

式中t為時間，y為因變數。y', dy/dt 均為其一階微分的表示式。m為t之函數。根據微分之定義，可以改寫如下：

  y'=dy/dt = limΔ-0 [y(t+Δt)-y(t)]/dt = f(t, y)

展開，可以求得經過Δt之時間後，y(t+Δt)之結果為：

 y(t+Δt)=y(t)+Δtf(t,y)

若設Δt為一均勻之區段時間，或稱為步進值，則可以給予y(t)一個初值，然後利用增加步進值之個數逐步累積至最終的區間，最後得到y之答案值。其求解之型式如下：

 y(tk+1) = y(tk) + Δtf(tk,y(tk))

故只要給定時間範圍及f(t0,y0)等初值，即可利用迴圈的方式得到最終答案。其中，Δt設定值之大小則會影響答案之準確度，Δt值愈大，準確度會差，Δt值愈小則運算時間加長。

例：設f(t,y)=-10y，初值t=0, y0=2。且y'=-10y 之真實解為y(t)=2e^-10t。設Δt＝0.02，以此利用MATLAB求解。

function y=euler(y0,time,delta)
% Using Euler method to solve differential eqs.
% Inputs:
%    y0:initial value
%    time:time limit
%    delta:step size
% Example: y=euler(2,0.5,0.02)
r=-10;k=0;
t=0:delta:time;y=zeros(size(t));
y(1)=y0;
for i=2:length(t),y(i)=y(i-1)+r*y(i-1)*delta;end
y_true=2*exp(-10*t);
plot(t,y,'o',t,y_true),xlabel('t'),ylabel('y');

執行結果：

>> y=euler(2,0.5,0.02); 

同樣的問題亦可使用ode45指令計算，則先設定其函數為fun：

fun=@(t,y)(-10*y);
[t,y]=ode45(fun,[0 0.5],2);
plot(t,y,'ro',t,2*exp(-10*t)),xlabel('t'),ylabel('y'); 

利用MATLAB之ODE45之指令所得之結果比前面所用的方法準確度高。以下我們將介紹數種微分方程之解法及相關指令。