Infinite loop in CVodeF

[aseyboldt]

Under certain conditions the forward solver of cvodes enters an infinite loop if called using CV_NORMAL_STEP. If the step-size decreases until t + h = t numerically, then CVode will print a warning, but still return a 0 return value. The loop here in CVodeF will continue forever.

I do not have a reproducing example in C, but here is one that uses the sunode package to wrap sundials:

We start integrating the ode dy/dt = 1 / (t - 1) at t0 = 1 + 1e-15:

import sunode
import sunode.symode
import numpy as np

params = {'p': ()}
states = {'a': ()}

def rhs(t, y, p):
    return {
        'a': 1 / (t - 1)
    }

problem = sunode.SympyProblem(
    params=params,
    states=states,
    derivative_params=[('p',)],
    rhs_sympy=rhs,
)

solver = sunode.solver.AdjointSolver(problem)
solver.set_params_dict({'p': 0.})

h = 1e-15
tvals = np.linspace(1 + h, 1 + 1000 * h)
out, *_ = solver.make_output_buffers(tvals)

y0 = np.array(0.)
solver.solve_forward(1 + h, tvals, y0, out)

The solver prints error messages repeatedly until the maximum number of warnings in exhausted, then silently continues:

[CVODES WARNING]  CVode
  Internal t = 1 and h = 4.96507e-17 are such that t + h = t on the next step. The solver will continue anyway.


[CVODES WARNING]  CVode
  Internal t = 1 and h = 4.96507e-17 are such that t + h = t on the next step. The solver will continue anyway.


[CVODES WARNING]  CVode
  Internal t = 1 and h = 4.96507e-17 are such that t + h = t on the next step. The solver will continue anyway.


[CVODES WARNING]  CVode
  Internal t = 1 and h = 4.96507e-17 are such that t + h = t on the next step. The solver will continue anyway.


[CVODES WARNING]  CVode
  The above warning has been issued mxhnil times and will not be issued again for this problem.

CC @astoeriko

[balos1]

What version of sundials does sunode interface to?

[aseyboldt]

This is using the latest version on conda-forge, so 5.3.0.

[aseyboldt]

@balos1 Is there anything I can do to help fix this issue and #31?
Both make it quite hard for me to actually use sundials. Weird parameter combinations happen a lot when sampling posteriors.

[balos1]

@aseyboldt A reproducing example in C would be very helpful. I am going to try and take a look at this soon.

[aseyboldt]

Thats great, thanks!
I'll port the example from the description tomorrow.

[aseyboldt]

I modified one of the examples so that it shows this infinite loop problem:
https://gist.github.com/aseyboldt/714633c5c3b6e8f91e1ec629d85aba0e

[balos1]

@aseyboldt I was able to reproduce the problem with that example, and we are working on fixing the problem.

[aseyboldt]

That's great.
Just to make sure we don't do the same thing more than one time unnecessarily: I think the patch in #45 solves the issue.

[aseyboldt]

Any news on this?

[balos1]

@aseyboldt This will be fixed in the next release.

[aseyboldt]

Thanks!