ode求解器的事件(Event)属性

转载 2015年07月07日 09:58:18

检测事件

matlab微分方程如何设定变量的范围 

dy1=y2
dy2=y1+1
其中y1的范围为0

【解】

M文件:
function [value,isterminal,direction] events1(t,y)
value y(1)-4;
isterminal= 1;
direction 0;

命令窗口:
dy @(t,y) [y(2);y(1) 1];
options=odeset('events',@events1);
[t,y] ode45(dy,[0 12],[0 1],options);
plot(t,y)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%function [value,isterminal,direction]=events(t,x)

事件检查函数,此时需要做的是过零点检测
ode45函数自动检查当value=0是否成立
如果我们要求检测Y=0的点,设置value=Y
这里我们要检测Y=4,那么就设置value=Y-4

isterminal检测到指定条件时,是否终止ode45函数的运行
1表示终止,0表示继续
在我们这个问题上,我们只要检测到零点时就停止程序



direction:value过零点检测的方向
-1表示由正到负,+1表示由负到正
对于我们这个问题,当然是由正到负

%.........................

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


一个例子和说明:

这是一个用到简单微分方程的物理情景
      一个质量m=100kg的物体从高处竖直落下,加速度会受到空气阻力的影响,这里简单的认为重力加速度g=9.8不变,空气阻力f=k*v^2  ,简单起见,k=1。只考虑竖直方向速度v,且速度,加速度,竖直位移都以向下为正。初速度,初位移都为0;那么有以下微分方程:
dy/dt=v                   
dv/dt=9.8-1*v^2/m
m=100,v0=y0=0

然后我用MATLAB的ode45函数求这个微分方程的数值解
先编写函数
function dx=fun(t,x)
% x(1)表示下落的距离y(向下为正),x(2)表示下落速度v(向下为正)
k=1;     % k=1为表示空气阻力的一个常量,这里简化空气阻力f=k*v^2
m=100;% m为质量=100kg
dx=zeros(2,1);
dx(1)=x(2);             % 下落距离对时间的导数=速度
a=9.8-(k*x(2)^2)/m;% a加速度(向下为正)=重力加速度 - 空气阻力产生的加速度
dx(2)=a;                 % 速度对时间的导数=加速度
end

现在我想要得到t=15s时的位移和速度
那么输入
[T,X]=ode45('fun',[0,15],[0 0]); 

返回的X中的最后一列就是我想要的值;

但假如我想知道当竖直向下的位移刚好=100米时的时间和速度,那该怎么办?现在我的做法是先将解一个充分大的时间,然后在里面找位移在100两侧的时间和速度,再通过插值得到位移刚好=100时的时间和速度。但这样很麻烦,也不见得准确,MATLAB有什么自带的语句能实现这个功能吗?或是有什么更好的方法?


在不知道结果时间的时候是需要先设定一个比较大的时间范围计算的
但是并不需要将整个范围的结果都算出来再插值
这个时候可以设定触发事件函数在一定条件下停止计算
用odeset可以为ode45求解器设定触发事件的函数
详细的用法要仔细查看matlab的帮助文件,这里我以你的题目,举个例子

微分方程还是用你的函数fun
在用ode45解方程之前,再写一个函数:事件触发函数eventfun,
它的格式固定要返回三个值,这三个值的意思是
当第一个值vaule的值到达0时,时间会触发
而根据第二个值isterminal设置,触发时间会否终止求解
第三个值是设置触发过0的方向

function [value,isterminal,direction] = eventfun(t,x)
value=x(1)-100; %触发时间,当其值为0的时候,时间会触发
isterminal=1;   %设为1时会,触发时间会停止求解器,设0时触发不影响工作
direction=1;    %触发方向设1时是上升触发,设-1是下降触发,设0是双向触发
end

写好fun和eventfun之后,你就可以调用ode45解方程了
但用ode45之前记得先用odeset,将触发函数加入哦

op=odeset('Events',@eventfun);
[T,X,Tend,Xend,evennum]=ode45(@fun,[0,15],[0 0],op); 

这样你看看,到达100米时,求解器就停住了
细心的你注意,ode45多返回了Tend,Xend,evennum三个参数
第一个Tend是触发事件发生的时间
第二个Xend是触发时间发生时刻的X
第三个evenum是标识触发事件的编号
由于这里只设置了一个触发事件,所以编号肯定是1

其实你也可以将该时间isterminal设为0
那么求解器不会因为触发事件而停止15秒内的解都会算出
你仍旧可以根据Tend 和 Xend得到到达100米时的时间和状态

额外提示,有时候你不知道到底取多大的时间范围才能够等到你要的触发时间
那么你可以用循环判断的方法,先设一个时间范围,然后求解
到最后都没有等到触发事件(Tend和Xend都是空矩阵)
那么适当延长求解时间区间,将上次的最后时刻和状态作为初始状态
再一次求解时间范围更大的解,如此直到找到触发事件才停止

引用:

1.http://zhidao.baidu.com/question/1046647733744158819.html?oldq=1

2.http://blog.sina.com.cn/s/blog_3ecbcc0701012u89.html


相关文章推荐

Javascript中Event事件的属性大全

1.altKey 描述: 检查alt键的状态。  语法: event.altKey  可能的值: 当alt键按下时,值为 TRUE ,否则为 FALSE 。只读。  2.button 描述:...

JavaScript Event学习第七章:事件属性

为了理解Event属性,我将在这里给出一些示例代码。在这个范畴内有非常严重的浏览器兼容性问题。 当我们想去读一读关于Event的一些资料时,常常会湮没在大量的属性里面,这些属性其中的大多数不能良...

HTML5中的新标签和新事件属性(Tags and Event Attributes)

Defines an article Defines content aside from the page content Defines sound content Isolates ...

兼容IE和FF的手动触发事件,并自定义Event的属性值

兼容IE和FF的手动触发事件,并自定义Event的属性值

javascript学习笔记—event属性和方法以及window事件

鼠标 / 键盘属性 属性 描述 IE F O W3C altKey 返回当事件被触发时,"ALT" 是否...

Web前端教程js-event(事件对象及其属性、方法)

js event

求解微分方程组的ODE算法

4阶5级龙格库塔法用于解一阶微分方程(组),对于高阶微分方程,可以将其转换为一阶微分方程组求解。原程序由John.H.Mathews编写(数值方法matlab版),但只能解微分方程,不能解微分方程组。...

JS-事件类型(鼠标事件中event对象的独有属性)

ClientX与ClientY返回客户区坐标位置 var bd=document.getElementById("bd");//body EventUtil.addHand...
内容举报
返回顶部
收藏助手
不良信息举报
您举报文章:深度学习:神经网络中的前向传播和反向传播算法推导
举报原因:
原因补充:

(最多只允许输入30个字)