Skip to content

Commit

Permalink
address more PR comments
Browse files Browse the repository at this point in the history
  • Loading branch information
@balos1
balos1 committed Jun 29, 2023
1 parent c010d8b commit 74383d7
Show file tree
Hide file tree
Showing 5 changed files with 43 additions and 129 deletions.
108 changes: 14 additions & 94 deletions doc/arkode/guide/source/Usage/SPRKStep_c_interface/User_callable.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -281,10 +281,6 @@ Optional inputs for SPRKStep
+-----------------------------------------------------+------------------------------------------+------------------------+
| Set fixed step size (disables time step adaptivity) | :c:func:`SPRKStepSetFixedStep()` | disabled |
+-----------------------------------------------------+------------------------------------------+------------------------+
| Supply an initial step size to attempt | :c:func:`SPRKStepSetInitStep()` | estimated |
+-----------------------------------------------------+------------------------------------------+------------------------+
| Maximum no. of warnings for :math:`t_n+h = t_n` | :c:func:`SPRKStepSetMaxHnilWarns()` | 10 |
+-----------------------------------------------------+------------------------------------------+------------------------+
| Maximum no. of internal steps before *tout* | :c:func:`SPRKStepSetMaxNumSteps()` | 500 |
+-----------------------------------------------------+------------------------------------------+------------------------+
| Set a value for :math:`t_{stop}` | :c:func:`SPRKStepSetStopTime()` | undefined |
Expand Down Expand Up @@ -455,29 +451,6 @@ Optional inputs for SPRKStep
calling :c:func:`SPRKStepEvolve()` to resume integration.
.. c:function:: int SPRKStepSetInitStep(void* arkode_mem, realtype hin)
Specifies the initial time step size SPRKStep should use after
initialization, re-initialization, or resetting.
:param arkode_mem: pointer to the SPRKStep memory block.
:param hin: value of the initial step to be attempted :math:`(\ne 0)`.
:return:
* *ARK_SUCCESS* if successful
* *ARK_MEM_NULL* if the SPRKStep memory is ``NULL``
* *ARK_ILL_INPUT* if an argument has an illegal value
**Notes:**
Pass 0.0 to use the default value.
By default, SPRKStep estimates the initial step size to be
:math:`h = \sqrt{\dfrac{2}{\left\| \ddot{y} \right\|}}`, where
:math:`\ddot{y}` is an estimate of the second derivative of the
solution at :math:`t_0`.
.. c:function:: int SPRKStepSetMaxNumSteps(void* arkode_mem, long int mxsteps)
Specifies the maximum number of steps to be taken by the
Expand All @@ -498,29 +471,6 @@ Optional inputs for SPRKStep
Passing *mxsteps* < 0 disables the test (not recommended).
.. c:function:: int SPRKStepSetMaxHnilWarns(void* arkode_mem, int mxhnil)
Specifies the maximum number of messages issued by the
solver to warn that :math:`t+h=t` on the next internal step, before
SPRKStep will instead return with an error.
**Arguments:**
* *arkode_mem* -- pointer to the SPRKStep memory block.
* *mxhnil* -- maximum allowed number of warning messages :math:`(>0)`.
**Return value:**
* *ARK_SUCCESS* if successful
* *ARK_MEM_NULL* if the SPRKStep memory is ``NULL``
* *ARK_ILL_INPUT* if an argument has an illegal value
**Notes:**
The default value is 10; set *mxhnil* to zero to specify
this default.
A negative value indicates that no warning messages should be issued.
.. c:function:: int SPRKStepSetStopTime(void* arkode_mem, realtype tstop)
Specifies the value of the independent variable
Expand Down Expand Up @@ -590,7 +540,7 @@ Optional inputs for IVP method selection

+-----------------------------+-------------------------------------------+----------+
| Optional input | Function name | Default |
+-----------------------------+-------------------------------------------+----------+
+=============================+===========================================+==========+
| Set integrator method order | :c:func:`SPRKStepSetOrder()` | 4 |
+-----------------------------+-------------------------------------------+----------+
| Set SPRK method | :c:func:`SPRKStepSetMethod()` | internal |
Expand Down Expand Up @@ -692,7 +642,7 @@ described in :numref:`ARKODE.Mathematics.Rootfinding`.
+-----------------------------------------+-------------------------------------------+----------+
| Optional input | Function name | Default |
+-----------------------------------------+-------------------------------------------+----------+
+=========================================+===========================================+==========+
| Direction of zero-crossings to monitor | :c:func:`SPRKStepSetRootDirection()` | both |
+-----------------------------------------+-------------------------------------------+----------+
| Disable inactive root warnings | :c:func:`SPRKStepSetNoInactiveRootWarn()` | enabled |
Expand Down Expand Up @@ -853,7 +803,7 @@ Main solver optional output functions

+-----------------------------------------------------+--------------------------------------------+
| Optional output | Function name |
+-----------------------------------------------------+--------------------------------------------+
+=====================================================+============================================+
| Cumulative number of internal steps | :c:func:`SPRKStepGetNumSteps()` |
+-----------------------------------------------------+--------------------------------------------+
| Step size used for the last successful step | :c:func:`SPRKStepGetLastStep()` |
Expand Down Expand Up @@ -1077,7 +1027,7 @@ Rootfinding optional output functions
+--------------------------------------------------+---------------------------------+
| Optional output | Function name |
+--------------------------------------------------+---------------------------------+
+==================================================+=================================+
| Array showing roots found | :c:func:`SPRKStepGetRootInfo()` |
+--------------------------------------------------+---------------------------------+
| No. of calls to user root function | :c:func:`SPRKStepGetNumGEvals()`|
Expand Down Expand Up @@ -1195,23 +1145,14 @@ Runge--Kutta stages, denoted by *s*, be no larger for the new problem than
for the previous problem. This condition is automatically fulfilled
if the method order *q* is left unchanged.
One important use of the :c:func:`SPRKStepReInit()` function is in the
treating of jump discontinuities in the RHS function. Except in cases
of fairly small jumps, it is usually more efficient to stop at each
point of discontinuity and restart the integrator with a readjusted
ODE model, using a call to this routine. To stop when the location
of the discontinuity is known, simply make that location a value of
``tout``. To stop when the location of the discontinuity is
determined by the solution, use the rootfinding feature. In either
case, it is critical that the RHS function *not* incorporate the
discontinuity, but rather have a smooth extension over the
discontinuity, so that the step across it (and subsequent rootfinding,
if used) can be done efficiently. Then use a switch within the RHS
function (communicated through ``user_data``) that can be flipped
between the stopping of the integration and the restart, so that the
restarted problem uses the new values (which have jumped). Similar
comments apply if there is to be a jump in the dependent variable
vector.
One potential use of the :c:func:`SPRKStepReInit()` function is in the
treating of jump discontinuities in the RHS function :cite:p:`Tao:22`.
In lieu of including if statements within the RHS function to handle
discontinuities, it may be more computationally efficient to stop at each
point of discontinuity (e.g., through use of tout or the rootfinding feature)
and restart the integrator with a readjusted ODE model, using a call to
this routine. We note that for the solution to retain temporal accuracy,
the RHS function should not incorporate the discontinuity.
.. c:function:: int SPRKStepReInit(void* arkode_mem, ARKRhsFn f1, ARKRhsFn f2, realtype t0, N_Vector y0)
Expand Down Expand Up @@ -1260,26 +1201,7 @@ via a call to :c:func:`SPRKStepSetStopTime()`. Following a successful call to
:c:func:`SPRKStepReset()`, call :c:func:`SPRKStepEvolve()` again to continue
solving the problem. By default the next call to :c:func:`SPRKStepEvolve()` will
use the step size computed by SPRKStep prior to calling :c:func:`SPRKStepReset()`.
To set a different step size or have SPRKStep estimate a new step size use
:c:func:`SPRKStepSetInitStep()`.
One important use of the :c:func:`SPRKStepReset()` function is in the
treating of jump discontinuities in the RHS functions. Except in cases
of fairly small jumps, it is usually more efficient to stop at each
point of discontinuity and restart the integrator with a readjusted
ODE model, using a call to :c:func:`SPRKStepReset()`. To stop when
the location of the discontinuity is known, simply make that location
a value of ``tout``. To stop when the location of the discontinuity
is determined by the solution, use the rootfinding feature. In either
case, it is critical that the RHS functions *not* incorporate the
discontinuity, but rather have a smooth extension over the
discontinuity, so that the step across it (and subsequent rootfinding,
if used) can be done efficiently. Then use a switch within the RHS
functions (communicated through ``user_data``) that can be flipped
between the stopping of the integration and the restart, so that the
restarted problem uses the new values (which have jumped). Similar
comments apply if there is to be a jump in the dependent variable
vector.
.. c:function:: int SPRKStepReset(void* arkode_mem, realtype tR, N_Vector yR)
Expand All @@ -1298,9 +1220,7 @@ vector.
**Notes:**
By default the next call to :c:func:`SPRKStepEvolve()` will use the step size
computed by SPRKStep prior to calling :c:func:`SPRKStepReset()`. To set a
different step size or have SPRKStep estimate a new step size use
:c:func:`SPRKStepSetInitStep()`.
computed by SPRKStep prior to calling :c:func:`SPRKStepReset()`.
All previously set options are retained but may be updated by calling
the appropriate "Set" functions.
Expand Down
43 changes: 28 additions & 15 deletions doc/shared/sundials.bib
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2032,24 +2032,37 @@ @article{Ver:78
}

@article{Ruth:93,
title={A canonical integration technique},
author={Ruth, Ronald D},
journal={IEEE Trans. Nucl. Sci.},
volume={30},
number={CERN-LEP-TH-83-14},
pages={2669--2671},
year={1983}
title = {A canonical integration technique},
author = {Ruth, Ronald D},
journal = {IEEE Trans. Nucl. Sci.},
volume = {30},
number = {CERN-LEP-TH-83-14},
pages = {2669--2671},
year = {1983}
}

@article{Sofroniou:05,
title={Derivation of symmetric composition constants for symmetric integrators},
author={Sofroniou, Mark and Spaletta, Giulia},
journal={Optimization Methods and Software},
volume={20},
number={4-5},
pages={597--613},
year={2005},
publisher={Taylor \& Francis}
title = {Derivation of symmetric composition constants for symmetric integrators},
author = {Sofroniou, Mark and Spaletta, Giulia},
journal = {Optimization Methods and Software},
volume = {20},
number = {4-5},
pages = {597--613},
year = {2005},
publisher = {Taylor \& Francis}
}

@article{Tao:22,
title = {Accurate and efficient simulations of Hamiltonian mechanical systems with discontinuous potentials},
journal = {Journal of Computational Physics},
volume = {450},
pages = {110846},
year = {2022},
issn = {0021-9991},
doi = {https://doi.org/10.1016/j.jcp.2021.110846},
url = {https://www.sciencedirect.com/science/article/pii/S0021999121007415},
author = {Molei Tao and Shi Jin},
keywords = {Symplectic integrator, Time-reversible / symmetric integrator, Hamiltonian with discontinuous potential, Contact and impact at interface, Reflection and refraction, Sauteed mushroom}
}

@article{KnWo:98,
Expand Down
2 changes: 1 addition & 1 deletion examples/arkode/C_serial/ark_kepler.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -63,7 +63,7 @@

#include <arkode/arkode_arkstep.h> /* prototypes for ARKStep fcts., consts */
#include <arkode/arkode_sprk.h>
#include <arkode/arkode_sprkstep.h> /* prototypes for MRIStep fcts., consts */
#include <arkode/arkode_sprkstep.h> /* prototypes for SPRKStep fcts., consts */
#include <math.h>
#include <nvector/nvector_serial.h> /* serial N_Vector type, fcts., macros */
#include <stdio.h>
Expand Down
4 changes: 0 additions & 4 deletions include/arkode/arkode_sprkstep.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -68,11 +68,7 @@ SUNDIALS_EXPORT int SPRKStepSetInterpolantType(void* arkode_mem, int itype);
SUNDIALS_EXPORT int SPRKStepSetInterpolantDegree(void* arkode_mem, int degree);
/* TODO(CJB): should we remove this from the initial release and wait for the OO
* adaptivity? */
SUNDIALS_EXPORT int SPRKStepSetAdaptivityFn(void* arkode_mem, ARKAdaptFn hfun,
void* h_data);
SUNDIALS_EXPORT int SPRKStepSetMaxNumSteps(void* arkode_mem, long int mxsteps);
SUNDIALS_EXPORT int SPRKStepSetMaxHnilWarns(void* arkode_mem, int mxhnil);
SUNDIALS_EXPORT int SPRKStepSetInitStep(void* arkode_mem, realtype hin);
SUNDIALS_EXPORT int SPRKStepSetStopTime(void* arkode_mem, realtype tstop);
SUNDIALS_EXPORT int SPRKStepSetFixedStep(void* arkode_mem, realtype hfixed);
SUNDIALS_EXPORT int SPRKStepSetErrHandlerFn(void* arkode_mem,
Expand Down
15 changes: 0 additions & 15 deletions src/arkode/arkode_sprkstep_io.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -58,16 +58,6 @@ int SPRKStepSetMaxNumSteps(void* arkode_mem, long int mxsteps)
return (arkSetMaxNumSteps(arkode_mem, mxsteps));
}

int SPRKStepSetMaxHnilWarns(void* arkode_mem, int mxhnil)
{
return (arkSetMaxHnilWarns(arkode_mem, mxhnil));
}

int SPRKStepSetInitStep(void* arkode_mem, realtype hin)
{
return (arkSetInitStep(arkode_mem, hin));
}

int SPRKStepSetStopTime(void* arkode_mem, realtype tstop)
{
return (arkSetStopTime(arkode_mem, tstop));
Expand Down Expand Up @@ -98,11 +88,6 @@ int SPRKStepSetPostprocessStageFn(void* arkode_mem, ARKPostProcessFn ProcessStag
return (arkSetPostprocessStageFn(arkode_mem, ProcessStage));
}

int SPRKStepSetAdaptivityFn(void* arkode_mem, ARKAdaptFn hfun, void* h_data)
{
return (arkSetAdaptivityFn(arkode_mem, hfun, h_data));
}

int SPRKStepSetFixedStep(void* arkode_mem, realtype hfixed)
{
return (arkSetFixedStep(arkode_mem, hfixed));
Expand Down

0 comments on commit 74383d7

Please sign in to comment.