matlab解常微分方程,Matlab中解常微分方程的ode45
ode是專門用于解微分方程的功能函數,他有ode23,ode45,ode23s等等,采用的是Runge-Kutta算法。ode45表示采用四階,五階runge-kutta單步算法,截斷誤差為(Δx)^3。解決的是Nonstiff(非剛性)的常微分方程.是解決數值解問題的首選方法,若長時間沒結果,應該就是剛性的,換用ode23來解.其他幾個也是類似的用法
使用方法
[T,Y] = ode45(odefun,tspan,y0)
odefun 是函數句柄,可以是函數文件名,匿名函數句柄或內聯函數名
tspan 是區間
[t0 tf] 或者一系列散點[t0,t1,...,tf]
y0 是初始值向量
T 返回列向量的時間點
Y 返回對應T的求解列向量
[T,Y] = ode45(odefun,tspan,y0,options)
options 是求解參數設置,可以用odeset在計算前設定誤差,輸出參數,事件等
[T,Y,TE,YE,IE] =ode45(odefun,tspan,y0,options)
每組(t,Y)之產生稱為事件函數。每次均會檢查是否函數等于零。并決定是否在零時終止運算。這可以在函數中之特性上設定。例如以events 或@events產生一函數。[value, isterminal,direction]=events(t,y)其中,value(i)為函數之值,isterminal(i)=1時運算在等于零時停止,=0時繼續;direction(i)=0時所有零時均需計算(默認值), +1在事件函數增加時等于零, -1在事件函數減少時等于零等狀況。此外,TE, YE, IE則分別為事件發生之時間,事件發生時之答案及事件函數消失時之指針i。
sol =ode45(odefun,[t0 tf],y0...)
sol 結構體輸出結果
應用舉例
1 求解一階常微分方程
程序:
odefun=@(t,y) (y+3*t)/t^2; %定義函數
tspan=[1 4]; %求解區間
y0=-2; %初值
[t,y]=ode45(odefun,tspan,y0);
plot(t,y) %作圖
title('t^2y''=y+3t,y(1)=-2,1
legend('t^2y''=y+3t') xlabel('t')
ylabel('y') % 精確解
% dsolve('t^2*Dy=y+3*t','y(1)=-2')
% ans =
% (3*Ei(1) - 2*exp(1))/exp(1/t) - (3*Ei(1/t))/exp(1/t)
2 求解高階常微分方程
關鍵是將高階轉為一階,odefun的書寫.
F(y,y',y''...y(n-1),t)=0用變量替換,y1=y,y2=y'...注意odefun方程定義為列向量
dxdy=[y(1),y(2)....]
程序:
function Testode45
tspan=[3.9 4.0]; %求解區間
y0=[2 8]; %初值
[t,x]=ode45(@odefun,tspan,y0);
plot(t,x(:,1),'-o',t,x(:,2),'-*')
legend('y1','y2')
title('y'' ''=-t*y + e^t*y'' +3sin2t')
xlabel('t') ylabel('y')
function y=odefun(t,x)
y=zeros(2,1); % 列向量
y(1)=x(2);
y(2)=-t*x(1)+exp(t)*x(2)+3*sin(2*t);
end
end
總結
以上是生活随笔為你收集整理的matlab解常微分方程,Matlab中解常微分方程的ode45的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: php array product,PH
- 下一篇: ubuntu如何修改php文件夹权限,L