ode45 Event function does not trigger correctly with nested function

Hello guys,

I tried to simplify my code to show the underlying problem:

I try to solve an ode with ode45:

[ t , y ] = ode45 @(t ,y) myOde 

This is the simplified ode:

function dydt = myOde
    dydt = func(y)
end

The problem is: func (y) is only determined for positive y and otherwise there will be an error.

For the actual code the solution for y will decrease to zero. I wrote an event function that will stop when the defined precision is reached.

function [position,isterminal,direction] = EventMyOde( ~ , y )
    precision = 1e-4 ;
    position =  y - precision; 
    isterminal = 1;  
    direction = 0;   
end

ode45 calculates y for different timesteps, gets closer to 0 and eventually crosses the precision defined in the Eventfunction.

BUT: To find the zero-crossing in the event function exaclty, ode45 keeps on calculating for some additional steps. In some cases it returns a negative y. So for the next evaluation func(y) (now undetermined, with a negative y) will cause an error.

Is there a way to stop this behaviour of ode45, or do you have some clever ideas for workarounds?

I had three Ideas to solve this issue but they did not really work:

• With odeset 'RelTol' and 'AbsTol' I can arrange the stepsize and stop the calculation so it won't cause this error. But I have to solve this kind of problem for a lot of parameterised functions, so this is not really an option I think. ( and does not solve the underlying problem.)
• I am not quite sure if I can use the option 'NonNegative', 'on'. In my testruns it did not work for me. The error is always triggered.
• Is there a way to build a try - catch block around this issue?

Thank you in advance!

Hans

Answer 1

[Raag]

Hi Hans,

The behaviour you are experiencing with ‘ode45’ is due to the solver's method of bracketing zero-crossings, which can cause slight overshoots into negative values. This is problematic when your function ‘func(y)’ is undefined for negative ‘y’.

To address this, you can adjust your event function to detect only downward crossings, which will help prevent negative evaluations:

function [position, isterminal, direction] = EventMyOde(~, y)
   precision = 1e-4;
   position = y - precision;
   isterminal = 1;  
   direction = -1;  % Detect when y crosses downward through 'precision'
end

By setting ‘direction = -1’, you ensure that the event is triggered only when y decreases through the threshold, which helps avoid negative values. Additionally, I would suggest you to consider setting 'NonNegative' to ‘1’ in your ‘odeset’ options.

For a better understanding of the above solution, refer to the following MATLAB documentations:

1. ODE Event Location: https://www.mathworks.com/help/releases/R2021a/matlab/math/ode-event-location.html

2. ode45: https://www.mathworks.com/help/releases/R2021a/matlab/ref/ode45.html

3. odeset: https://www.mathworks.com/help/releases/R2021a/matlab/ref/odeset.html

Answer 2

[Walter Roberson]

Use odeset() to set NonNegative for the component; https://www.mathworks.com/help/matlab/math/nonnegative-ode-solution.html