CVode does not return CV_TOO_CLOSE when initial step size is set

[eltix]

Issue description

To illustrate our questions, here is a MWE in the form of a simple one-equation ODE system:

$$\begin{aligned} & \dot y = \mathrm{mod}(t, 1)\\\ & y(t=0) = 0 \\\ & 0 \le t \le 10 \end{aligned}$$

where $\mathrm{mod}(t, 1)$ stands for the fractional part of $t$.
The exact solution is a piece-wise parabola and we are mostly interested in the value at t=10 which is $y(t=10)=5$.

How we solved this problem using CVode

See the example code in the "To replicate" section below.

1. By default, CVode may give incorrect results because the ODE RHS is not continuous and therefore the solver may "step over" the discontinuity should it take a large enough step.

2. This is expected so we added a root function which triggers at each discontinuity of the RHS, that is for $t = k, k \in \mathbb N$:

static int root_fn(sunrealtype t, N_Vector y, sunrealtype *gout, void *user_data) {
  sunrealtype y1, b;

  y1 = Ith(y, 1);
  b = SUN_RCONST(1);

  // a jig-saw function which has roots at t = k, k being any integer
  *gout = fabs( fmod(t-b/2, 2*b) - b) - b/2;
  return 0;
}

But now CVode() returns CV_TOO_CLOSE:

[ERROR][rank 0][/build/sundials-7.2.1/src/cvode/cvode.c:3710][cvHandleFailure] tout too close to t0 to start integration.
At t_out = 1.00e+01 t = 1.000000e+01 y1 = 5.0e+00 tcur = 1.0e+01 step = 0.0e+00 order = 0 err1 = -8.8e-07 
retval = -27, return 1

3. This is somewhat expected since the root time coincides exactly with the asked time t=10. So we special cased CV_TOO_CLOSE. At this point, if you're trying to reproduce the issue, you should uncomment this block at lines 153-159:

    else if (retval == CV_TOO_CLOSE) {
      /* This happens if tout is too close to t. In such a case, no solving was required; just continue as if solving succeeded. 
      However, if the init step is set using CVodeSetInitStep, CVode DOES NOT RETURN CV_TOO_CLOSE, see below
      */
      printf("retval == CV_TOO_CLOSE\n");
      iout++;
    }

4. Up until now, the init step size was not set. Let's set it (uncomment lines 109-110):

  retval = CVodeSetInitStep(cvode_mem, INITSTEP);
  if (check_retval(&retval, "CVodeSetInitStep", 1)) { return (1); }

Now we have an infinite loop because CVode() no longer returns CV_TOO_CLOSE!

5. This is where things get interesting and we are seeking guidance. To prevent the infinite loop, we wrote an additional safeguard (uncomment lines 160-167):

    else if (retval == CV_ROOT_RETURN && t == tout) {
      /* See the above
         if the init step is set using CVodeSetInitStep, CVode DOES NOT RETURN CV_TOO_CLOSE when tout is too close to t
         In such a case, pretend that the root didn't happen, lest we keep handling it forever. 
      */
      printf("t == tout && retval == CV_ROOT_RETURN\n");
      iout++;
    }

This seems to work but the results are actually incorrect since we have y(t=10) = 0.5 instead of 5:

At t_out = 1.00e+01 t = 1.000000e+01 y1 = 5.0e-01 tcur = 1.0e+01 step = 1.0e+00 order = 1 err1 = 0.0e+00

The questions

1. Is it expected that CVode() does not return CV_TOO_CLOSE when init step size is set?
2. In general, how should we handle cases where the root finding function triggers on a time that coincides with a requested solution time?

To replicate

1. Sundials version: 7.2.1

2. example.zip contains an example.c source file which you can compile with

gcc example.c -lsundials_cvode -lsundials_core -lm

and then run with

./a.out

[sundials-users-github-account]

Hello Eliott Tixier, Here are a few comments that occur to me offhand. 1. Since the discontinuities in the RHS function occur at known times, the use of rootfinding is the wrong fix. The easier and more robust fix is to simply stop (using tout) at each discontinuity (t = k in this example), reset the RHS function, and restart CVODE. The RHS reset can be done through user_data, as shown in the example cvRocket_dns. 2. In general, if CVODE is stopped at a root and restarted, the next tout should be beyond the value of t at the root, to avoid the TOO_CLOSE error return. Also, before restarting, the root function should be reset so that the same root is not seen twice. This way there is no possibility of an infinite loop. (cvRocket_dns also illustrates this reset.) 3. A CVODE start with the internal initial step calculation (i.e. not using CVodeSetInitStep) will always give an error return if tout = t0, because the difference tout - t0 is used in that calculation. If the initial step is supplied, that error return is avoided. Perhaps this should be changed.

…

-Alan H

________________________________ From: sundials-users ***@***.***> on behalf of Eliott Tixier ***@***.***> Sent: Friday, June 20, 2025 6:57 AM To: sundials-users ***@***.***> Subject: [sundials-users] [LLNL/sundials] CVode does not return CV_TOO_CLOSE when initial step size is set (Issue #722) [https://avatars.githubusercontent.com/u/28726150?s=20&v=4]eltix created an issue (LLNL/sundials#722)<https://urldefense.us/v3/__https://github.com/LLNL/sundials/issues/722__;!!G2kpM7uM-TzIFchu!1ex3sXtnjl9RU2gwaAfU8UvJnKbgK4OdJDdNkHxbI45GB_kxVBv4lH_ZoebzjLds63xPavHgBJcD9cAOAenkBj1kdGq5xA$> Issue description To illustrate our questions, here is a MWE in the form of a simple one-equation ODE system: $$\begin{aligned} & \dot y = \mathrm{mod}(t, 1)\\\ & y(t=0) = 0 \\\ & 0 \le t \le 10 \end{aligned}$$ where $\mathrm{mod}(t, 1)$ stands for the fractional part of $t$. The exact solution is a piece-wise parabola and we are mostly interested in the value at t=10 which is $y(t=10)=5$. How we solved this problem using CVode See the example code in the "To replicate" section below. 1. By default, CVode may give incorrect results because the ODE RHS is not continuous and therefore the solver may "step over" the discontinuity should it take a large enough step. 2. This is expected so we added a root function which triggers at each discontinuity of the RHS, that is for $t = k, k \in \mathbb N$: static int root_fn(sunrealtype t, N_Vector y, sunrealtype *gout, void *user_data) { sunrealtype y1, b; y1 = Ith(y, 1); b = SUN_RCONST(1); // a jig-saw function which has roots at t = k, k being any integer *gout = fabs( fmod(t-b/2, 2*b) - b) - b/2; return 0; } But now CVode() returns CV_TOO_CLOSE: [ERROR][rank 0][/build/sundials-7.2.1/src/cvode/cvode.c:3710][cvHandleFailure] tout too close to t0 to start integration. At t_out = 1.00e+01 t = 1.000000e+01 y1 = 5.0e+00 tcur = 1.0e+01 step = 0.0e+00 order = 0 err1 = -8.8e-07 retval = -27, return 1 1. This is somewhat expected since the root time coincides exactly with the asked time t=10. So we special cased CV_TOO_CLOSE. At this point, if you're trying to reproduce the issue, you should uncomment this block at lines 153-159: else if (retval == CV_TOO_CLOSE) { /* This happens if tout is too close to t. In such a case, no solving was required; just continue as if solving succeeded. However, if the init step is set using CVodeSetInitStep, CVode DOES NOT RETURN CV_TOO_CLOSE, see below */ printf("retval == CV_TOO_CLOSE\n"); iout++; } 1. Up until now, the init step size was not set. Let's set it (uncomment lines 109-110): retval = CVodeSetInitStep(cvode_mem, INITSTEP); if (check_retval(&retval, "CVodeSetInitStep", 1)) { return (1); } Now we have an infinite loop because CVode() no longer returns CV_TOO_CLOSE! 1. This is where things get interesting and we are seeking guidance. To prevent the infinite loop, we wrote an additional safeguard (uncomment lines 160-167): else if (retval == CV_ROOT_RETURN && t == tout) { /* See the above if the init step is set using CVodeSetInitStep, CVode DOES NOT RETURN CV_TOO_CLOSE when tout is too close to t In such a case, pretend that the root didn't happen, lest we keep handling it forever. */ printf("t == tout && retval == CV_ROOT_RETURN\n"); iout++; } This seems to work but the results are actually incorrect since we have y(t=10) = 0.5 instead of 5: At t_out = 1.00e+01 t = 1.000000e+01 y1 = 5.0e-01 tcur = 1.0e+01 step = 1.0e+00 order = 1 err1 = 0.0e+00 The questions 1. Is it expected that CVode() does not return CV_TOO_CLOSE when init step size is set? 2. In general, how should handle cases where the root finding function triggers on a time that coincides with a requested solution time? To replicate 1. Sundials version: 7.2.1 2. example.zip<https://urldefense.us/v3/__https://github.com/user-attachments/files/20837371/example.zip__;!!G2kpM7uM-TzIFchu!1ex3sXtnjl9RU2gwaAfU8UvJnKbgK4OdJDdNkHxbI45GB_kxVBv4lH_ZoebzjLds63xPavHgBJcD9cAOAenkBj1-R5Z70A$> contains an example.c source file which you can compile with gcc example.c -lsundials_cvode -lsundials_core -lm and then run with ./a.out — Reply to this email directly, view it on GitHub<https://urldefense.us/v3/__https://github.com/LLNL/sundials/issues/722__;!!G2kpM7uM-TzIFchu!1ex3sXtnjl9RU2gwaAfU8UvJnKbgK4OdJDdNkHxbI45GB_kxVBv4lH_ZoebzjLds63xPavHgBJcD9cAOAenkBj1kdGq5xA$>, or unsubscribe<https://urldefense.us/v3/__https://github.com/notifications/unsubscribe-auth/AXXIGEVJMONTTH7GDTRK4HL3EQHNTAVCNFSM6AAAAAB7YMUAJOVHI2DSMVQWIX3LMV43ASLTON2WKOZTGE3DGNBUGEYTONI__;!!G2kpM7uM-TzIFchu!1ex3sXtnjl9RU2gwaAfU8UvJnKbgK4OdJDdNkHxbI45GB_kxVBv4lH_ZoebzjLds63xPavHgBJcD9cAOAenkBj1j6cSNvw$>. You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

________________________________ To unsubscribe from the SUNDIALS-USERS list: write to: ***@***.***

[eltix]

@sundials-users-github-account

Hello Alan, thank you very much for your answer.

1. Since the discontinuities in the RHS function occur at known times, the use of rootfinding is the wrong fix

The above was just an example, in general we do not know when the roots trigger, hence the need to use root finding. For fixed-time events (where we know in advance at what time they occur), we simply stop and re-init the solver as you explained.

2. before restarting, the root function should be reset so that the same root is not seen twice

In the cvRocket example, CVodeRootInit is called mid-solving to change the number of root functions from 2 (running out of gas + max height) to 1 (zero velocity). While this works in the rocket example (because you can't run out of gas, then have gas again then run out of gas once more), in general we cannot do that for a given root function may trigger multiple times during the simulation.

3. If the initial step is supplied, that error return is avoided. Perhaps this should be changed.

Since this was a surprise from a user perspective, I would say that indeed the same error should be returned when the initial step is supplied such that the behavior is consistent.

[sundials-users-github-account]

Hi Elliot, I suspected that your real problem does not have known stop-times, but wanted to make the distinction clear. And for root-defined stops, there is very little chance that the stop time will coincide exactly with tout. I will recommend that the TOO_CLOSE error test be added when the initial step is supplied. But it is always better for the user program logic to avoid having the initial tout equal to t0.

…

-Alan H

________________________________ From: sundials-users ***@***.***> on behalf of Eliott Tixier ***@***.***> Sent: Thursday, June 26, 2025 12:47 AM To: sundials-users ***@***.***> Subject: Re: [sundials-users] [LLNL/sundials] CVode does not return CV_TOO_CLOSE when initial step size is set (Issue #722) [https://avatars.githubusercontent.com/u/28726150?s=20&v=4]eltix left a comment (LLNL/sundials#722)<https://urldefense.us/v3/__https://github.com/LLNL/sundials/issues/722*issuecomment-3007525949__;Iw!!G2kpM7uM-TzIFchu!w1mIgGmzoSlMK1oZ5R2c_NJcRxtKDr_t6OZmPhP3r530EZLaYyjC8gzuyjHJc7h9BOdQsi87Oom2dxHHl4r_w5NLADsE2Q$> @sundials-users-github-account<https://urldefense.us/v3/__https://github.com/sundials-users-github-account__;!!G2kpM7uM-TzIFchu!w1mIgGmzoSlMK1oZ5R2c_NJcRxtKDr_t6OZmPhP3r530EZLaYyjC8gzuyjHJc7h9BOdQsi87Oom2dxHHl4r_w5N3L2G1gA$> Dear Alan, thank you for your answer 1. Since the discontinuities in the RHS function occur at known times, the use of rootfinding is the wrong fix The above was just an example, in general we do not know when the roots trigger, hence the need to use root finding. For fixed-time events (where we know in advance at what time they occur), we simply stop and re-init the solver as you explained. 1. before restarting, the root function should be reset so that the same root is not seen twice This is the part we were missing, thank you very much 1. If the initial step is supplied, that error return is avoided. Perhaps this should be changed. Since this was a surprise from a user perspective, I would say that indeed the same error should be returned when the initial step is supplied such that the behavior is consistent. — Reply to this email directly, view it on GitHub<https://urldefense.us/v3/__https://github.com/LLNL/sundials/issues/722*issuecomment-3007525949__;Iw!!G2kpM7uM-TzIFchu!w1mIgGmzoSlMK1oZ5R2c_NJcRxtKDr_t6OZmPhP3r530EZLaYyjC8gzuyjHJc7h9BOdQsi87Oom2dxHHl4r_w5NLADsE2Q$>, or unsubscribe<https://urldefense.us/v3/__https://github.com/notifications/unsubscribe-auth/AXXIGESGBH6EHNP7SJYVLLD3FOQPTAVCNFSM6AAAAAB7YMUAJOVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZTAMBXGUZDKOJUHE__;!!G2kpM7uM-TzIFchu!w1mIgGmzoSlMK1oZ5R2c_NJcRxtKDr_t6OZmPhP3r530EZLaYyjC8gzuyjHJc7h9BOdQsi87Oom2dxHHl4r_w5Mvj-TD6A$>. You are receiving this because you were mentioned.Message ID: ***@***.***>

________________________________ To unsubscribe from the SUNDIALS-USERS list: write to: ***@***.***

[eltix]

@sundials-users-github-account

Hi Alan,

Thank you again for your answer.

And for root-defined stops, there is very little chance that the stop time will coincide exactly with tout

However unlikely, this did happen to us, hence the admittedly hackish if (retval == CV_TOO_CLOSE) check we added.

I should have added more context as to how we use Sundials. We are primarily interested in solving Quantitative Systems Pharmacology models which we programatically convert to IVP problems before solving them with CVODE. We often deal with models we did not create ourselves, for instance SBML models from public databases which can for instance define events that trigger when the condition t >= te is satisifed, te being also part of the requested solution times.
This is just an example and if we want something robust, we must be completely agnostic as to what the root conditions and requested solution times look like.

I will recommend that the TOO_CLOSE error test be added when
the initial step is supplied. But it is always better for the
user program logic to avoid having the initial tout equal to t0

Agreed. Ultimately it would be ideal if the following logic that evaluates whether tout "is too close to" t0

  /* If tout is too close to tn, give up */

  if ((tdiff = tout - cv_mem->cv_tn) == ZERO) { return (CV_TOO_CLOSE); }

  sign   = (tdiff > ZERO) ? 1 : -1;
  tdist  = SUNRabs(tdiff);
  tround = cv_mem->cv_uround * SUNMAX(SUNRabs(cv_mem->cv_tn), SUNRabs(tout));

  if (tdist < TWO * tround) { return (CV_TOO_CLOSE); }

were extracted out of cvHin() such that:

1. it can be called in CVode() regardless of whether the initial step is supplied
2. it is exposed to users such that we can test whether tout is equal to t0 in the user program