ARKODE run as a fully implicit ESDIRK method (singly
diagonally implicit Runge-Kutta with an explicit first stage) with
Newton iteration and a dense linear solver. Its
sundialr_arkode_* entry points mirror the
sundialr_cvode_* ones one for one — the same callback
types, handle lifecycle, tolerance and step-size setters,
set_udata, step counting, and the same no-throw,
no-allocation-on-reuse and one-handle-per-thread contracts — so a host
loop written against the CVODE API can switch integrator by changing
which entry points it looks up. The method defaults to
ARKODE_ESDIRK436L2SA_6_3_4 (order 4);
sundialr_arkode_set_order() picks ARKODE’s default ESDIRK
of order 2 to 5 and sundialr_arkode_set_table_name() any
ARKODE DIRK table by name. Status codes are ARKODE’s, mirrored as
SUNDIALR_ARK_*; they coincide with the CVODE codes for the
common failures but not everywhere.
sundialr_arkode_reinit() is warm: it uses
ARKodeReset, which keeps the step size and the error
controller’s history from the previous segment, the natural restart
after a dose within one trajectory.
sundialr_arkode_reinit_cold() is the full re-initialisation
(ARKStepReInit), with CVodeReInit’s semantics:
it leaves a used handle exactly as a new one would be (pinned bit for
bit by a test), so a handle reused across unrelated problems, such as
the subjects of a population, gives results that do not depend on what
it solved before; the step count keeps running across it. The addition
is binary-compatible, so SUNDIALR_ABI_VERSION stays 1. To
support it the bundled build now compiles ARKODE, with
CRAN-compliance patches for the stdout references in its
parameter-printing pathsSUNNonlinearSolver_Auto stiffness switch no longer moves to
Newton iteration prematurely, and IDA no longer writes
duplicate keys in its logging output (logging is disabled in this
build). The sundialr interface and the C API are unchanged,
and the existing CRAN-compliance patches apply to the new sources
without modification. The installed headers change accordingly: the
N_Vector operations table gains an
nvsetdevicearraypointer entry, so code compiled against
these headers must also be built against 7.9.0 sourcescmake call now uses the current option names
SUNDIALS_ENABLE_C_EXAMPLES and
SUNDIALS_ENABLE_CXX_EXAMPLES, silencing two deprecation
warnings during installationopen module. The step
that normalises the line endings of the two cmake-generated
headers used perl -Mopen=IO,:raw, and open is
not part of the minimal perl-base, so installing from
source on a minimal Debian/Ubuntu image (for example
rocker/r-ver) failed with “Can’t locate open.pm in @INC”. The step now uses POSIX
tr, which behaves the same on Linux, macOS and Windows
(Rtools)sundialr_cvode_set_udata() replaces the callback-data
pointer on an existing handle, so one handle can serve many parameter
sets — a population of subjects — instead of paying a create/free cycle
(SUNContext, cvode_mem, matrix and linear
solver allocations) per set. sundialr_cvode_reset_stats()
zeroes the accumulated step count so
sundialr_cvode_get_num_steps() can report per-subject
counts from a shared handle. Both are additive and binary-compatible, so
SUNDIALR_ABI_VERSION stays 1sundialr_cvode_get_num_steps() under the
reinit-per-segment pattern: CVodeReInit zeroes
CVODE’s internal step counter, so the count previously
restarted at every sundialr_cvode_reinit() and never
reflected more than the current segment, despite being documented as a
running total. The handle now accumulates the count across reinits,
counting from create() or the last
reset_stats()sundialr_capi.h and pinned by tests: after
the first reinit+solve cycle, further
reinit/solve calls on a handle perform no heap
allocation (a tight-loop test drives 10^4 segments and checks the state
buffer is never reallocated), and distinct handles — each owning its own
SUNContext and error record, with no global mutable state
and no R API on any entry point — may be driven from concurrent threads,
one handle per thread (pinned by a two-thread test)ARKODE and
CVODE from compiled code”, showing how another package uses
the C API: it looks up the sundialr_arkode_* entry points
with R_GetCCallable, simulates repeated doses of a
one-compartment model with ARKODE, restarting the
integrator at each dose, and checks the result against the closed-form
solution. It also covers method selection, warm and cold
reinitialisation and switching to CVODE. The Introduction
and Summary no longer describe ARKODE as a future addition;
it is available through the C API, though not yet as an R
functioncvode(),
cvodes(), ida() and cvsolve()
accept an optional jacobian argument, an R
function giving the Jacobian of the system analytically. When it is not
supplied, which remains the default, SUNDIALS approximates
the Jacobian by finite differences as before, so existing code is
unaffected. For cvode(), cvodes() and
cvsolve() the function has the same signature as the system
itself, function(t, y, p), and returns an n-by-n matrix
whose [i, j] entry is d(ydot[i])/d(y[j]). For
ida() it is function(t, y, ydot, cj, p) and
returns dF/dy + cj * dF/dydot, cj being a
scalar supplied by the solver. Supplying the Jacobian is usually faster
and more accurate on stiff systemsinst/include/sundialr_capi.h and registered as
R_RegisterCCallable entry points under the package name
sundialr. A consumer reaches it with
R_GetCCallable("sundialr", ...) after
Imports: sundialr; the header exposes only
double/int and an opaque handle, so no
SUNDIALS headers or LinkingTo are required. The handle owns
a persistent state vector and is reused across many segments via
sundialr_cvode_reinit(), the pattern a host integrator loop
needs; tolerances may be scalar or per-equation (the latter realised
with a custom error-weight function, so vector relative tolerances are
supported), an analytic Jacobian is optional, and every entry point
returns a CVODE status code and records a retrievable
message rather than throwing across the foreign call boundary. It lets
CVODE be embedded as an integrator inside another package’s
own solve loop; the existing R interface is unchangedcvodes() accepts an
optional sensitivity argument, an R function
giving the sensitivity right-hand side analytically. When it is not
supplied, which remains the default, SUNDIALS approximates
the sensitivity equations by finite differences of the state right-hand
side as before, so existing code is unaffected. The function has
signature function(t, y, ydot, iS, yS, p) and returns
d(yS_iS)/dt = J %*% yS_iS + df/dp_iS as a numeric vector of
length(y), where iS is the 1-based index of
the parameter. It is independent of the jacobian argument —
the caller embeds the Jacobian inside their own sensitivity
function — so supplying one does not require supplying the other.
Providing the sensitivity right-hand side is usually faster and more
accurate than the finite-difference default, which evaluates the state
right-hand side an extra length(p) times per internal
stepsrc/scripts/cran_patches.sh) which now verifies that no
flagged call survives patchingtools/strip_sundials_tarball.sh to reproducibly
strip the bundled SUNDIALS tarball of examples, docs, tests and
benchmarksEvents record at
time 0 in cvsolve() now adds its value to the initial
condition instead of replacing it, so the IC supplied by
the caller is always the starting state and t = 0 behaves like every
other event time. Several events at t = 0 on one state now accumulate
rather than the last one winning. Results change for any call that doses
at time 0: with IC = c(1) and a value of 100 at t = 0 the
state now starts at 101 rather than 100, and
inst/examples/cvsolve_1D.r starts at 101 accordingly. To
reproduce the old numbers, subtract the initial condition from the value
of any event at time 0cvsolve() treating the initial time as 0 rather
than as the first element of time_vector. An event at the
initial time is now applied to the initial conditions as intended;
previously, if the time vector began anywhere other than 0, such an
event was silently dropped and its value lost with no error or
warningcvsolve() now rejects an event time outside the range
of time_vector instead of accepting it. An event before the
first output time cannot be applied, since the solve starts there; it
previously produced a row of zeros in the returned matrix and discarded
the event. Event times that are NA are rejected for the
same reasoncvsolve() no longer returns rows of zeros. A record
that does not advance the solution — one at the initial time, or one
sharing a time with the record before it — was skipped without storing
anything, leaving that row filled with zeros instead of the stateida() and cvsolve() where a
scalar abstolerance (the default) allocated a length-1
tolerance vector instead of one entry per state. On systems with more
than one state this read and wrote past the end of that vector, and the
solvers failed immediately with a corrector convergence error at the
initial time. Passing abstolerance as a vector of the same
length as IC was unaffected, as were cvode()
and cvodes()cvode(), cvodes(), ida() and
cvsolve() no longer leak memory when a call ends in an
error. The SUNDIALS objects behind each solve are now released by a
scope guard instead of by free calls at the end of the function, so they
are also freed when an error unwinds the stack. Previously a failed call
leaked everything allocated up to the point of the error, which
accumulated when solvers were called repeatedly (for example from an
optimiser or a fitting loop). This includes errors raised by SUNDIALS
itself, such as a failure to converge: these previously aborted the call
in a way that skipped all cleanup, and were the most common source of
the leak, since a solver that fails to converge is exactly the case a
user is likely to retry in a loop. Separately, cvsolve()
never freed its constraints vector on any path, including successful
onescvsolve() now rejects an invalid state index in the
first column of Events. The index is converted to a
position in the state vector and used to subscript it directly, so a
value outside 1:length(IC) wrote outside that vector: an
index of 0, the natural mistake when counting from zero,
corrupted the heap and could abort the R session.
Fractional and NA indices are rejected for the same reason,
a fractional index having previously been truncated to a different state
without warningR function is
now returned to the solver as a failure code rather than being thrown
out through SUNDIALS’ own C code, which lets SUNDIALS unwind its
internal state instead of being abandoned part-way through a step. The
R error itself is unchanged: the original condition still
reaches the caller, so tryCatch() on a specific condition
class continues to workR function returning the wrong number
of values is now reported instead of being read past the end. Previously
a right-hand side or residual function returning fewer values than there
are states, or a Jacobian function returning a matrix of the wrong
shape, was copied out element by element regardless and the solve
continued on uninitialised memoryida() and
cvsolve(), which previously had none, including a check
that a scalar and an equivalent vector abstolerance give
identical results, and for the callback error handling aboveCVodeSetLinearSolver was reported
as CVDlsSetLinearSolver (a name removed in SUNDIALS 4) in
cvode() and cvsolve(),
CVodeWFtolerances as CVodeSetEwtFn and
CVodeSensInit1 as CVodeSensInit in
cvodes(), and IDASetUserData and
IDASolve as their cvode() equivalents in
ida()/tools/cmake_call.sh file to remove warnings
with respect to stdout/stderr/fprintf etc in upstream C libraryUpdated build system to use cmake
Bug fixes: * Fixed the bug in assigning absolute tolerances
Bug fixes: * Fixed the linking problem due to multiple defined symbols
New features: * added the cvsolve function to solve ODEs
with multiple discontinuities
Bug fixes: * None
API change: * NONE is existing functions * A new function
cvsolve is added
Internal changes: * Now uses SUNDIALS version 5.2.0 (released March 2020) at the back end