ODE45的意义

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

matlab ode45 解微分方程

在用odesolver(ode45, ode15s, …)来解微分方程的时候,最基本的用法是:

[t, y] = odesolver(odefun, tspan, y0);

这里的odefun是待求的微分方程。那么odefun中一般会含有多个系统参数,通常要通过改变参数来观察系统动态的变化。那么如何在调用odesolver的时候传递参数呢?

以前,我都是用全局变量的写法,将参数在主函数和子函数中分别都定义为global,这样做有一个弱点:针对系统不同,参数的表达与数量有变化的时候,程序通常要做变化,通用性不强。那么最好是在调用的时候进行传递,方法如下:

实际上很简单,就是将一切其他的参数都写在括号中就可以了!但是要注意的是:odesolver的第四个参数一定是options,也就是对微分方程添加补充功能的参数(类型为structure,要用odeset来定义),那么其他系统参数就只能从第五个参数写起。也就是说,第四个参数不可以为空,一定要定义某种option加进去,或者用使用空白矩阵(placeholder)。这样调用的时候格式就是:

[t, y] = odesolver(odefun, tspan, y0, options, parameter1, parameter2);

或者

[t, y] = odesolver(odefun, tspan, y0, [], parameter1, parameter2);

然后定义微分方程的时候也要有参数的地方:

function dydt = odefun(t, y, parameter1, parameter2)

dydt = [ eqn-1; eqn-2; …];

就OK了。

另一种用法:

[T,Y] = ode45(@(t,y) rigid(t,y,Lambda,Lami,Cs,ionization,Ed),[0 100],[PHI_0,E_0,M_0,Gma_0],options); Lambda,Lami,Cs,ionization,Ed是参数之前有声明t y是变量

,[PHI_0,E_0,M_0,Gma_0],是初始值

--------------------------------------------------

function dy = rigid(t,y,Lambda,Lami,Cs,ionization,Ed)

dy = zeros(4,1); % a column vector

dy(1) = y(2);

dy(2) = exp(y(1))-y(3)/y(4);

dy(3)=Lambda*ionization/Cs*exp(y(1));

dy(4) = -y(4)*Lambda/Lami-Lambda*ionization/Cs-y(2)/y(4)-Ed/y(4);

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

[t,H]=ode45('solitiontry1',tspan,h0,[],m1,epsinon) ——————————————————————————————

function solfi=solitiontry1(t,H,flag,m1,epsinon)

%m1=1.15;

y=(1-epsinon)*(1-exp(H))+epsinon*m1^2*(1-sqrt(1+2*H*m1^-2))+m1^2*(1-sqrt(1-2*H*m1^-2));

solfi=-sqrt(-2*y);

;;注意前面两个函数之间的区别,第二个多了一股flag参数。应该是针对options为空的那个参数的,但是在这里,如果options不为空的话,还不知道有什么区别。

[原创]Matlab中解常微分方程的ode45

2010-09-24 15:18

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);

相关文档
最新文档