Non-negative ODE solution

Is there any Julia option equivalent to “Nonnegative” in “odeset” in Matlab?

Here is the documentation about the function I’m talking about:

https://www.mathworks.com/help/matlab/math/nonnegative-ode-solution.html

I think it’s a very useful function and it’d be nice to have something similar in Julia.

Thanks in advance!

Hey,
There are two ways to do this. One is to use isoutofdomain. It’s documented with the other solver options:

http://docs.juliadiffeq.org/latest/basics/common_solver_opts.html#Miscellaneous-1

isoutofdomain = (u,p,t)->any(x->x<0,u) makes any step with a negative value get rejected and use a smaller dt. This is a method that is always safe, but may get a few more rejections.

The other method is to use PositiveDomain() as a callback on the problem. This uses continuous extensions in the event handling to perform the pullback to the positive time.

http://docs.juliadiffeq.org/latest/features/callback_library.html#PositiveDomain-1

This is more like MATLAB’s version. This is less safe than isoutofdomain since it requires that the function is differentiable around zero in order for the interpolation to make sense, but it can be (a tiny bit) faster.

4 Likes

Hi Chris,

it seems these links have broken sometime in the past year.

I tried the isoutofdomain approach and I got the following warning and the simulation would stop early.

Is there an example available to see where the line PositiveDomain() goes? Thanks!

Best,
Brad

Usually that happens if the true solution is negative. Are you sure it’s not?

How horribly embaressing… I thought I had correctly addressed the local vs. global scope situation with regards to using the eval function inside of the function to solve (e.g. latka_volterra(du,u,p,t)), but it seems I didn’t fully. I wasn’t updating states correctly that depending on the values of other states. It’s working beautifully now and I’ll begin fitting soon. :grin:

function tau_agg!(du, u_local, p, t)
global u = u_local
du[:] = eval.(doyle_eqs[:,5])
end

This appears to work now. :slight_smile:

Shall I delete my comments to prevent confusion for future onlookers?

Sincerest apologies for taking your time on a personal fumble. I’ll try to make it the last.

No worries. I think I need to just keep dropping that nugget around that isoutofdomain isn’t magic: it makes things not numerically leave the domain by forcing solver rejections whenever it does, but if the solution truly leaves the domain, you’ll get dtmin issues as it tries to make the calculation more and more accurate to not leave the domain, but dt->0 as it accurately keeps leaving the domain.

2 Likes