Azzera filtri
Azzera filtri

The Matlab event function doesn't stop the integration

2 visualizzazioni (ultimi 30 giorni)
Hi,
I have a ode system that I solve with Matlab. I want to find the steady state of the system and to do so, I use the event function as described here .
But some times, the solver doesn't stop even if the criterium is achieved. For example, if you solve the following system with x0 = 10 the solver will stop before 2000 but with x0 = 0.0001 it won't.
The event function (eventfun_t.m)
function [x,isterm,dir] = eventfun_t(t,y)
dy = test_systeme(t,y);
x = norm(dy) - 1e-3;
isterm = 1;
dir = 0; %or -1, doesn't matter
end
The system (test_systeme.m)
function dx = test_systeme(t,x)
v = x./(x+1);
dx = -v;
end
Solve the system
x0 = 10;
eventfonction = @(t,y) eventfun_t(t,y);
optionsode=odeset('Events',eventfonction);
[t x]=ode15s(@(t,x) test_systeme(t,x),[0 2000],x0,optionsode);
I suppose it's because with x0 = 0.0001 norm(dy) is already lower than 1e-3 but in that case, how can I stop the solver without checking by myself ?
  4 Commenti
Jan
Jan il 19 Ott 2018
Modificato: Jan il 19 Ott 2018
This is the first comment to the answer in the thread you've posted a link to.
dy = test_systeme(t,y);
if t~= 0 && x < 1e-3
x = 0
This cannot work, because x is undefined. Even if you mean dy instead of x, this is not a useful event function. Remember, that the event is triggered, when the 1st output of this functions becomes 0. Because it is very unlikely, that 0 is hit exactly, the integrators detect a change of the sign of this value and adjust the step size, to locate the exact zero crossing. Then "if t~= 0 && x < 1e-3, x = 0;" is not useful, because this is not a zero-crossing at t==0 even for tiny dy.
I recommend to read the instructions again: doc: ode-event-location
I did not test it, but Andrew's comment sounds reasonable: The ODE integrators do not detect an event at the initial step. So run the event function manually before starting the integration. This is cheap and easy, so I do not see a reason to implement any further tricks.

Accedi per commentare.

Risposta accettata

Cécile Moulin
Cécile Moulin il 19 Ott 2018
The anwser is in the comments of the question : there is no trick and I have to check before the integration if the conditions are already achieved.

Più risposte (0)

Categorie

Scopri di più su Programming in Help Center e File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by