3.4.3.1. A skeleton of the user’s main program
The following is a skeleton of the user’s main program (or calling program) for the integration of an ODE IVP using the ERKStep module. Most of the steps are independent of the NVECTOR implementation used. For the steps that are not, refer to §9 for the specific name of the function to be called or macro to be referenced.
Initialize parallel or multi-threaded environment, if appropriate.
For example, call
MPI_Init
to initialize MPI if used, or setnum_threads
, the number of threads to use within the threaded vector functions, if used.Create the SUNDIALS simulation context object.
Call
SUNContext_Create()
to allocate theSUNContext
object.Set problem dimensions, etc.
This generally includes the problem size,
N
, and may include the local vector lengthNlocal
.Note
The variables
N
andNlocal
should be of typesunindextype
.Set vector of initial values
To set the vector
y0
of initial values, use the appropriate functions defined by the particular NVECTOR implementation.For native SUNDIALS vector implementations (except the CUDA and RAJA based ones), use a call of the form
y0 = N_VMake_***(..., ydata);
if the
sunrealtype
arrayydata
containing the initial values of \(y\) already exists. Otherwise, create a new vector by making a call of the formy0 = N_VNew_***(...);
and then set its elements by accessing the underlying data where it is located with a call of the form
ydata = N_VGetArrayPointer_***(y0);
For details on each of SUNDIALS’ provided vector implementations, see the corresponding sections in §9 for details.
Create ERKStep object
Call
arkode_mem = ERKStepCreate(...)
to create the ERKStep memory block.ERKStepCreate()
returns avoid*
pointer to this memory structure. See §3.4.3.2.1 for details.Specify integration tolerances
Call
ERKStepSStolerances()
orERKStepSVtolerances()
to specify either a scalar relative tolerance and scalar absolute tolerance, or a scalar relative tolerance and a vector of absolute tolerances, respectively. Alternatively, callERKStepWFtolerances()
to specify a function which sets directly the weights used in evaluating WRMS vector norms. See §3.4.3.2.2 for details.Set optional inputs
Call
ERKStepSet*
functions to change any optional inputs that control the behavior of ERKStep from their default values. See §3.4.3.2.5 for details.Specify rootfinding problem
Optionally, call
ERKStepRootInit()
to initialize a rootfinding problem to be solved during the integration of the ODE system. See §3.4.3.2.3 for general details, and §3.4.3.2.5 for relevant optional input calls.Advance solution in time
For each point at which output is desired, call
ier = ERKStepEvolve(arkode_mem, tout, yout, &tret, itask);
Here,
itask
specifies the return mode. The vectoryout
(which can be the same as the vectory0
above) will contain \(y(t_\text{out})\). See §3.4.3.2.4 for details.Get optional outputs
Call
ERKStepGet*
functions to obtain optional output. See §3.4.3.2.7 for details.Deallocate memory for solution vector
Upon completion of the integration, deallocate memory for the vector
y
(oryout
) by calling the NVECTOR destructor function:N_VDestroy(y);
Free solver memory
Call
ERKStepFree()
to free the memory allocated for the ERKStep module.Finalize MPI, if used
Call
MPI_Finalize
to terminate MPI.