Controlling output limits in Scipy's ODEINT
<p>I have the following simplified pseudocode setup which involves a system of ODEs, which I attempt to solve with scipy's odeint.</p>

<pre><code>from scipy.integrate import odeint

def diff_func(y, time, parms):

    # Do stuff with parms that depends upon y and t.
    new_parms = other_funcs(y,time, parms)

    # Now calculate derivatives
    dy1_dt = dy1_func(y, new_parms)
    dy2_dt = dy2_func(y, new_parms)

# Setup up initial conditions
y_0 = [y_1_0, y_2_0]
time = np.linspace(0, 1000, 1000)
parms = list_o_constants

# Solve diffeq
yout, info = odeint(diff_func, y_0, time, args=(parms,), full_output=True)
</code></pre>

<p>Depending upon the inputs I give it, I can lead to nonsensical results - or crashes. I have narrowed down the problem to the <code>y</code> inputs go negative occasionally, since one of them is temperature - it leads to all kinds of issues. The reason why it is going negative in the first place is unclear, and might be a stability issue.</p>

<p>I have attempted to change the diff function to include a check like:</p>

<pre><code>def diff_func(y, time, parms):
    if y[1] &lt;= 0.1:
        y[1] = 0.1
</code></pre>

<p>This seems to work okay, I am still testing it. But, I wanted to ask if: this is going to cause problems down the line, and if there is a better way to force odeint to only allow positive outputs (i.e., change step size if a negative is found).</p>

<p>Thanks!</p>