8.1. Description of the NVECTOR Modules
SUNDIALS solvers are written in a data-independent manner. They all operate on
generic vectors (of type N_Vector
) through a set of operations defined
by, and specific to, the particular vector implementation. Users can provide a
custom vector implementation or use one provided with SUNDIALS. The generic
operations are described below. In the sections following, the implementations
provided with SUNDIALS are described.
An N_Vector
is a pointer to the _generic_N_Vector
structure:
-
typedef struct _generic_N_Vector *N_Vector
-
struct _generic_N_Vector
The structure defining the SUNDIALS vector class.
-
void *content
Pointer to vector-specific member data.
-
N_Vector_Ops ops
A virtual table of vector operations provided by a specific implementation.
-
SUNContext sunctx
The SUNDIALS simulation context
-
void *content
The virtual table structure is defined as
-
typedef _generic_N_Vector_Ops *N_Vector_Ops
-
struct _generic_N_Vector_Ops
The structure defining
N_Vector
operations.-
N_Vector_ID (*nvgetvectorid)(N_Vector)
The function implementing
N_VGetVectorID()
-
N_Vector (*nvclone)(N_Vector)
The function implementing
N_VClone()
-
N_Vector (*nvcloneempty)(N_Vector)
The function implementing
N_VCloneEmpty()
-
void (*nvdestroy)(N_Vector)
The function implementing
N_VDestroy()
-
void (*nvspace)(N_Vector, sunindextype*, sunindextype*)
The function implementing
N_VSpace()
-
sunrealtype *(*nvgetarraypointer)(N_Vector)
The function implementing
N_VGetArrayPointer()
-
sunrealtype *(*nvgetdevicearraypointer)(N_Vector)
The function implementing
N_VGetDeviceArrayPointer()
-
void (*nvsetarraypointer)(sunrealtype*, N_Vector)
The function implementing
N_VSetArrayPointer()
-
SUNComm (*nvgetcommunicator)(N_Vector)
The function implementing
N_VGetCommunicator()
-
sunindextype (*nvgetlength)(N_Vector)
The function implementing
N_VGetLength()
-
sunindextype (*nvgetlocallength)(N_Vector)
The function implementing
N_VGetLocalLength()
-
void (*nvlinearsum)(sunrealtype, N_Vector, sunrealtype, N_Vector, N_Vector)
The function implementing
N_VLinearSum()
-
void (*nvconst)(sunrealtype, N_Vector)
The function implementing
N_VConst()
-
void (*nvscale)(sunrealtype, N_Vector, N_Vector)
The function implementing
N_VScale()
-
void (*nvaddconst)(N_Vector, sunrealtype, N_Vector)
The function implementing
N_VAddConst()
-
sunrealtype (*nvdotprod)(N_Vector, N_Vector)
The function implementing
N_VDotProd()
-
sunrealtype (*nvmaxnorm)(N_Vector)
The function implementing
N_VMaxNorm()
-
sunrealtype (*nvwrmsnorm)(N_Vector, N_Vector)
The function implementing
N_VWrmsNorm()
-
sunrealtype (*nvwrmsnormmask)(N_Vector, N_Vector, N_Vector)
The function implementing
N_VWrmsNormMask()
-
sunrealtype (*nvmin)(N_Vector)
The function implementing
N_VMin()
-
sunrealtype (*nvwl2norm)(N_Vector, N_Vector)
The function implementing
N_VWL2Norm()
-
sunrealtype (*nvl1norm)(N_Vector)
The function implementing
N_VL1Norm()
-
void (*nvcompare)(sunrealtype, N_Vector, N_Vector)
The function implementing
N_VCompare()
-
sunbooleantype (*nvinvtest)(N_Vector, N_Vector)
The function implementing
N_VInvTest()
-
sunbooleantype (*nvconstrmask)(N_Vector, N_Vector, N_Vector)
The function implementing
N_VConstrMask()
-
sunrealtype (*nvminquotient)(N_Vector, N_Vector)
The function implementing
N_VMinQuotient()
-
SUNErrCode (*nvlinearcombination)(int, sunrealtype*, N_Vector*, N_Vector)
The function implementing
N_VLinearCombination()
-
SUNErrCode (*nvscaleaddmulti)(int, sunrealtype*, N_Vector, N_Vector*, N_Vector*)
The function implementing
N_VScaleAddMulti()
-
SUNErrCode (*nvdotprodmulti)(int, N_Vector, N_Vector*, sunrealtype*)
The function implementing
N_VDotProdMulti()
-
SUNErrCode (*nvlinearsumvectorarray)(int, sunrealtype, N_Vector*, sunrealtype, N_Vector*, N_Vector*)
The function implementing
N_VLinearSumVectorArray()
-
SUNErrCode (*nvscalevectorarray)(int, sunrealtype*, N_Vector*, N_Vector*)
The function implementing
N_VScaleVectorArray()
-
SUNErrCode (*nvconstvectorarray)(int, sunrealtype, N_Vector*)
The function implementing
N_VConstVectorArray()
-
SUNErrCode (*nvwrmsnormvectorarray)(int, N_Vector*, N_Vector*, sunrealtype*)
The function implementing
N_VWrmsNormVectorArray()
-
SUNErrCode (*nvwrmsnormmaskvectorarray)(int, N_Vector*, N_Vector*, N_Vector, sunrealtype*)
The function implementing
N_VWrmsNormMaskVectorArray()
-
SUNErrCode (*nvscaleaddmultivectorarray)(int, int, sunrealtype*, N_Vector*, N_Vector**, N_Vector**)
The function implementing
N_VScaleAddMultiVectorArray()
-
SUNErrCode (*nvlinearcombinationvectorarray)(int, int, sunrealtype*, N_Vector**, N_Vector*)
The function implementing
N_VLinearCombinationVectorArray()
-
sunrealtype (*nvdotprodlocal)(N_Vector, N_Vector)
The function implementing
N_VDotProdLocal()
-
sunrealtype (*nvmaxnormlocal)(N_Vector)
The function implementing
N_VMaxNormLocal()
-
sunrealtype (*nvminlocal)(N_Vector)
The function implementing
N_VMinLocal()
-
sunrealtype (*nvl1normlocal)(N_Vector)
The function implementing
N_VL1NormLocal()
-
sunbooleantype (*nvinvtestlocal)(N_Vector, N_Vector)
The function implementing
N_VInvTestLocal()
-
sunbooleantype (*nvconstrmasklocal)(N_Vector, N_Vector, N_Vector)
The function implementing
N_VConstrMaskLocal()
-
sunrealtype (*nvminquotientlocal)(N_Vector, N_Vector)
The function implementing
N_VMinQuotientLocal()
-
sunrealtype (*nvwsqrsumlocal)(N_Vector, N_Vector)
The function implementing
N_VWSqrSumLocal()
-
sunrealtype (*nvwsqrsummasklocal)(N_Vector, N_Vector, N_Vector)
The function implementing
N_VWSqrSumMaskLocal()
-
SUNErrCode (*nvdotprodmultilocal)(int, N_Vector, N_Vector*, sunrealtype*)
The function implementing
N_VDotProdMultiLocal()
-
SUNErrCode (*nvdotprodmultiallreduce)(int, N_Vector, sunrealtype*)
The function implementing
N_VDotProdMultiAllReduce()
-
SUNErrCode (*nvbufsize)(N_Vector, sunindextype*)
The function implementing
N_VBufSize()
-
SUNErrCode (*nvbufpack)(N_Vector, void*)
The function implementing
N_VBufPack()
-
SUNErrCode (*nvbufunpack)(N_Vector, void*)
The function implementing
N_VBufUnpack()
-
void (*nvprint)(N_Vector)
The function implementing
N_VPrint()
-
void (*nvprintfile)(N_Vector, FILE*)
The function implementing
N_VPrintFile()
-
N_Vector_ID (*nvgetvectorid)(N_Vector)
The generic NVECTOR module defines and implements the vector
operations acting on a N_Vector
. These routines are nothing but
wrappers for the vector operations defined by a particular NVECTOR
implementation, which are accessed through the ops field of the
N_Vector
structure. To illustrate this point we show below the
implementation of a typical vector operation from the generic NVECTOR
module, namely N_VScale
, which performs the operation \(z\gets cx\)
for vectors \(x\) and \(z\) and a scalar \(c\):
void N_VScale(sunrealtype c, N_Vector x, N_Vector z) {
z->ops->nvscale(c, x, z);
}
§8.2 contains a complete list of all standard vector operations defined by the generic NVECTOR module. §8.2.2, §8.2.3, §8.2.4, §8.2.5, and §8.2.6 list optional fused, vector array, local reduction, single buffer reduction, and exchange operations, respectively.
Fused and vector array operations (see §8.2.2 and
§8.2.3) are intended to increase data reuse, reduce
parallel communication on distributed memory systems, and lower the number of
kernel launches on systems with accelerators. If a particular NVECTOR
implementation defines a fused or vector array operation as NULL
, the
generic NVECTOR module will automatically call standard vector operations as
necessary to complete the desired operation. In all SUNDIALS-provided
NVECTOR implementations, all fused and vector array operations are
disabled by default. However, these implementations provide
additional user-callable functions to enable/disable any or all of the
fused and vector array operations. See the following sections
for the implementation specific functions to enable/disable operations.
Local reduction operations (see §8.2.4) are
similarly intended to reduce parallel
communication on distributed memory systems, particularly when
NVECTOR objects are combined together within an NVECTOR_MANYVECTOR
object (see §8.22). If a
particular NVECTOR implementation defines a local reduction
operation as NULL
, the NVECTOR_MANYVECTOR module will
automatically call standard vector reduction operations as necessary
to complete the desired operation. All SUNDIALS-provided NVECTOR
implementations include these local reduction operations, which may be
used as templates for user-defined implementations.
The single buffer reduction operations
(§8.2.5) are used in low-synchronization
methods to combine separate reductions into one MPI_Allreduce
call.
The exchange operations (see §8.2.6) are intended only for use with the XBraid library for parallel-in-time integration (accessible from ARKODE) and are otherwise unused by SUNDIALS packages.
8.1.1. NVECTOR Utility Functions
The generic NVECTOR module also defines several utility functions to aid in
creation and management of arrays of N_Vector
objects – these functions
are particularly useful for Fortran users to utilize the NVECTOR_MANYVECTOR
or SUNDIALS’ sensitivity-enabled packages CVODES and IDAS.
The functions N_VCloneVectorArray()
and
N_VCloneVectorArrayEmpty()
create (by cloning) an array of count
variables of type N_Vector
, each of the same type as an existing
N_Vector
input:
-
N_Vector *N_VCloneVectorArray(int count, N_Vector w)
Clones an array of
count
N_Vector
objects, allocating their data arrays (similar toN_VClone()
).- Arguments:
count
– number ofN_Vector
objects to create.w
– templateN_Vector
to clone.
- Return value:
pointer to a new
N_Vector
array on success.NULL
pointer on failure.
-
N_Vector *N_VCloneVectorArrayEmpty(int count, N_Vector w)
Clones an array of
count
N_Vector
objects, leaving their data arrays unallocated (similar toN_VCloneEmpty()
).- Arguments:
count
– number ofN_Vector
objects to create.w
– templateN_Vector
to clone.
- Return value:
pointer to a new
N_Vector
array on success.NULL
pointer on failure.
An array of variables of type N_Vector
can be destroyed
by calling N_VDestroyVectorArray()
:
-
void N_VDestroyVectorArray(N_Vector *vs, int count)
Destroys an array of
count
N_Vector
objects.- Arguments:
vs
–N_Vector
array to destroy.count
– number ofN_Vector
objects invs
array.
- Notes:
This routine will internally call the
N_Vector
implementation-specificN_VDestroy()
operation.If
vs
was allocated usingN_VCloneVectorArray()
then the data arrays for eachN_Vector
object will be freed; ifvs
was allocated usingN_VCloneVectorArrayEmpty()
then it is the user’s responsibility to free the data for eachN_Vector
object.
Finally, we note that users of the Fortran 2003 interface may be interested in
the additional utility functions N_VNewVectorArray()
,
N_VGetVecAtIndexVectorArray()
, and N_VSetVecAtIndexVectorArray()
,
that are wrapped as FN_NewVectorArray
, FN_VGetVecAtIndexVectorArray
, and
FN_VSetVecAtIndexVectorArray
, respectively. These functions allow a Fortran
2003 user to create an empty vector array, access a vector from this array, and
set a vector within this array:
-
N_Vector *N_VNewVectorArray(int count, SUNContext sunctx)
Creates an array of
count
N_Vector
objects, the pointers to each are initialized asNULL
.- Arguments:
count
– length of desiredN_Vector
array.sunctx
– aSUNContext
object
- Return value:
pointer to a new
N_Vector
array on success.NULL
pointer on failure.
Changed in version 7.0.0: The function signature was updated to add the
SUNContext
argument.
-
N_Vector *N_VGetVecAtIndexVectorArray(N_Vector *vs, int index)
Accesses the
N_Vector
at the locationindex
within theN_Vector
arrayvs
.- Arguments:
vs
–N_Vector
array.index
– desiredN_Vector
to access from withinvs
.
- Return value:
pointer to the indexed
N_Vector
on success.NULL
pointer on failure (index < 0
orvs == NULL
).
- Notes:
This routine does not verify that
index
is within the extent ofvs
, sincevs
is a simpleN_Vector
array that does not internally store its allocated length.
-
void N_VSetVecAtIndexVectorArray(N_Vector *vs, int index, N_Vector w)
Sets a pointer to
w
at the locationindex
within the vector arrayvs
.- Arguments:
vs
–N_Vector
array.index
– desired location to place the pointer tow
withinvs
.w
–N_Vector
to set withinvs
.
- Notes:
This routine does not verify that
index
is within the extent ofvs
, sincevs
is a simpleN_Vector
array that does not internally store its allocated length.
8.1.2. Implementing a custom NVECTOR
A particular implementation of the NVECTOR module must:
Specify the content field of the
N_Vector
structure.Define and implement the vector operations. Note that the names of these routines should be unique to that implementation in order to permit using more than one NVECTOR module (each with different
N_Vector
internal data representations) in the same code.Define and implement user-callable constructor and destructor routines to create and free an
N_Vector
with the new content field and with ops pointing to the new vector operations.Optionally, define and implement additional user-callable routines acting on the newly-defined
N_Vector
(e.g., a routine to print the content for debugging purposes).Optionally, provide accessor macros as needed for that particular implementation to be used to access different parts in the content field of the newly-defined
N_Vector
.
To aid in the creation of custom NVECTOR modules, the generic NVECTOR module
provides two utility functions N_VNewEmpty()
and
N_VCopyOps()
. When used in custom NVECTOR constructors and clone
routines these functions will ease the introduction of any new optional vector
operations to the NVECTOR API by ensuring that only required operations need
to be set, and that all operations are copied when cloning a vector.
-
N_Vector N_VNewEmpty(SUNContext sunctx)
This allocates a new generic
N_Vector
object and initializes its content pointer and the function pointers in the operations structure toNULL
.Return value: If successful, this function returns an
N_Vector
object. If an error occurs when allocating the object, then this routine will returnNULL
.
-
void N_VFreeEmpty(N_Vector v)
This routine frees the generic
N_Vector
object, under the assumption that any implementation-specific data that was allocated within the underlying content structure has already been freed. It will additionally test whether the ops pointer isNULL
, and, if it is not, it will free it as well.- Arguments:
v – an N_Vector object
-
SUNErrCode N_VCopyOps(N_Vector w, N_Vector v)
This function copies the function pointers in the
ops
structure ofw
into theops
structure ofv
.- Arguments:
w – the vector to copy operations from
v – the vector to copy operations to
Return value: Returns a
SUNErrCode
.
-
enum N_Vector_ID
Each
N_Vector
implementation included in SUNDIALS has a unique identifier specified in enumeration and shown in Table 8.1. It is recommended that a user supplied NVECTOR implementation use theSUNDIALS_NVEC_CUSTOM
identifier.
Vector ID |
Vector type |
ID Value |
---|---|---|
SUNDIALS_NVEC_SERIAL |
Serial |
0 |
SUNDIALS_NVEC_PARALLEL |
Distributed memory parallel (MPI) |
1 |
SUNDIALS_NVEC_OPENMP |
OpenMP shared memory parallel |
2 |
SUNDIALS_NVEC_PTHREADS |
PThreads shared memory parallel |
3 |
SUNDIALS_NVEC_PARHYP |
hypre ParHyp parallel vector |
4 |
SUNDIALS_NVEC_PETSC |
PETSc parallel vector |
5 |
SUNDIALS_NVEC_CUDA |
CUDA vector |
6 |
SUNDIALS_NVEC_HIP |
HIP vector |
7 |
SUNDIALS_NVEC_SYCL |
SYCL vector |
8 |
SUNDIALS_NVEC_RAJA |
RAJA vector |
9 |
SUNDIALS_NVEC_OPENMPDEV |
OpenMP vector with device offloading |
10 |
SUNDIALS_NVEC_TRILINOS |
Trilinos Tpetra vector |
11 |
SUNDIALS_NVEC_MANYVECTOR |
“ManyVector” vector |
12 |
SUNDIALS_NVEC_MPIMANYVECTOR |
MPI-enabled “ManyVector” vector |
13 |
SUNDIALS_NVEC_MPIPLUSX |
MPI+X vector |
14 |
SUNDIALS_NVEC_CUSTOM |
User-provided custom vector |
15 |
8.1.3. Support for complex-valued vectors
While SUNDIALS itself is written under an assumption of real-valued data, it does provide limited support for complex-valued problems. However, since none of the built-in NVECTOR modules supports complex-valued data, users must provide a custom NVECTOR implementation for this task. Many of the NVECTOR routines described in the subsection §8.2 naturally extend to complex-valued vectors; however, some do not. To this end, we provide the following guidance:
N_VMin()
andN_VMinLocal()
should return the minimum of all real components of the vector, i.e., \(m = \displaystyle \min_{0\le i< n} \operatorname{real}(x_i)\).N_VConst()
(and similarlyN_VConstVectorArray()
) should set the real components of the vector to the input constant, and set all imaginary components to zero, i.e., \(z_i = c + 0 j\) for \(0\le i<n\).N_VAddConst()
should only update the real components of the vector with the input constant, leaving all imaginary components unchanged.N_VWrmsNorm()
,N_VWrmsNormMask()
,N_VWSqrSumLocal()
andN_VWSqrSumMaskLocal()
should assume that all entries of the weight vectorw
and the mask vectorid
are real-valued.N_VDotProd()
should mathematically return a complex number for complex-valued vectors; as this is not possible with SUNDIALS’ currentsunrealtype
, this routine should be set toNULL
in the custom NVECTOR implementation.N_VCompare()
,N_VConstrMask()
,N_VMinQuotient()
,N_VConstrMaskLocal()
andN_VMinQuotientLocal()
are ill-defined due to the lack of a clear ordering in the complex plane. These routines should be set toNULL
in the custom NVECTOR implementation.
While many SUNDIALS solver modules may be utilized on complex-valued data,
others cannot. Specifically, although each package’s linear solver
interface (e.g., ARKLS or CVLS) may be used on complex-valued problems,
none of the built-in SUNMatrix or SUNLinearSolver modules will work (all
of the direct linear solvers must store complex-valued data, and all of
the iterative linear solvers require N_VDotProd()
). Hence a
complex-valued user must provide custom linear solver modules for their
problem. At a minimum this will consist of a custom SUNLinearSolver
implementation (see §10.1.8), and optionally a
custom SUNMatrix as well. The user should then attach these modules as
normal to the package’s linear solver interface.
Similarly, although both the
SUNNonlinearSolver_Newton and
SUNNonlinearSolver_FixedPoint modules
may be used with any of the IVP solvers (CVODE(S), IDA(S) and ARKODE) for
complex-valued problems, the Anderson-acceleration option with
SUNNonlinearSolver_FixedPoint cannot be used due to its reliance on
N_VDotProd()
. By this same logic, the Anderson acceleration
feature within KINSOL will also not work with complex-valued vectors.
Finally, constraint-handling features of each package cannot be used
for complex-valued data, due to the issue of
ordering in the complex plane discussed above with
N_VCompare()
, N_VConstrMask()
,
N_VMinQuotient()
, N_VConstrMaskLocal()
and
N_VMinQuotientLocal()
.
We provide a simple example of a complex-valued example problem,
including a custom complex-valued Fortran 2003 NVECTOR module, in the
files examples/arkode/F2003_custom/ark_analytic_complex_f2003.f90
,
examples/arkode/F2003_custom/fnvector_complex_mod.f90
, and
examples/arkode/F2003_custom/test_fnvector_complex_mod.f90
.
8.2. Description of the NVECTOR operations
8.2.1. Standard vector operations
The standard vector operations defined by the generic N_Vector
module are defined as follows. For each of these operations, we give
the name, usage of the function, and a description of its mathematical
operations below.
-
N_Vector_ID N_VGetVectorID(N_Vector w)
Returns the vector type identifier for the vector
w
. It is used to determine the vector implementation type (e.g. serial, parallel, …) from the abstractN_Vector
interface. Returned values are given in Table 8.1.Usage:
id = N_VGetVectorID(w);
-
N_Vector N_VClone(N_Vector w)
Creates a new
N_Vector
of the same type as an existing vector w and sets the ops field. It does not copy the vector, but rather allocates storage for the new vector.Usage:
v = N_VClone(w);
-
N_Vector N_VCloneEmpty(N_Vector w)
Creates a new
N_Vector
of the same type as an existing vector w and sets the ops field. It does not allocate storage for the new vector’s data.Usage:
v = N VCloneEmpty(w);
-
void N_VDestroy(N_Vector v)
Destroys the
N_Vector
v and frees memory allocated for its internal data.Usage:
N_VDestroy(v);
-
void N_VSpace(N_Vector v, sunindextype *lrw, sunindextype *liw)
Returns storage requirements for the
N_Vector
v:lrw contains the number of
sunrealtype
wordsliw contains the number of integer words.
This function is advisory only, for use in determining a user’s total space requirements; it could be a dummy function in a user-supplied NVECTOR module if that information is not of interest.
Usage:
N_VSpace(nvSpec, &lrw, &liw);
-
sunrealtype *N_VGetArrayPointer(N_Vector v)
Returns a pointer to a
sunrealtype
array from theN_Vector
v. Note that this assumes that the internal data in theN_Vector
is a contiguous array ofsunrealtype
and is accessible from the CPU.This routine is only used in the solver-specific interfaces to the dense and banded (serial) linear solvers, and in the interfaces to the banded (serial) and band-block-diagonal (parallel) preconditioner modules provided with SUNDIALS.
Usage:
vdata = N_VGetArrayPointer(v);
-
sunrealtype *N_VGetDeviceArrayPointer(N_Vector v)
Returns a device pointer to a
sunrealtype
array from theN_Vector
v
. Note that this assumes that the internal data inN_Vector
is a contiguous array ofsunrealtype
and is accessible from the device (e.g., GPU).This operation is optional except when using the GPU-enabled direct linear solvers.
Usage:
vdata = N_VGetArrayPointer(v);
-
void N_VSetArrayPointer(sunrealtype *vdata, N_Vector v)
Replaces the data array pointer in an
N_Vector
with a given array ofsunrealtype
. Note that this assumes that the internal data in theN_Vector
is a contiguous array ofsunrealtype
. This routine is only used in the interfaces to the dense (serial) linear solver, hence need not exist in a user-supplied NVECTOR module.Usage:
N_VSetArrayPointer(vdata,v);
-
SUNComm N_VGetCommunicator(N_Vector v)
Returns the
SUNComm
(which is just anMPI_Comm
when SUNDIALS is built with MPI, otherwise it is anint
) associated with the vector (if applicable). For MPI-unaware vector implementations, this should returnSUN_COMM_NULL
.Usage:
MPI_Comm comm = N_VGetCommunicator(v); // Works if MPI is enabled int comm = N_VGetCommunicator(v); // Works if MPI is disabled SUNComm comm = N_VGetCommunicator(v); // Works with or without MPI
-
sunindextype N_VGetLength(N_Vector v)
Returns the global length (number of “active” entries) in the NVECTOR v. This value should be cumulative across all processes if the vector is used in a parallel environment. If v contains additional storage, e.g., for parallel communication, those entries should not be included.
Usage:
global_length = N_VGetLength(v);
-
sunindextype N_VGetLocalLength(N_Vector v)
Returns the local length (number of “active” entries) in the NVECTOR v. This value should be the length of the array returned by
N_VGetArrayPointer()
orN_VGetDeviceArrayPointer()
.Usage:
local_length = N_VGetLocalLength(v);
-
void N_VLinearSum(sunrealtype a, N_Vector x, sunrealtype b, N_Vector y, N_Vector z)
Performs the operation z = ax + by, where a and b are
sunrealtype
scalars and x and y are of typeN_Vector
:\[z_i = a x_i + b y_i, \quad i=0,\ldots,n-1.\]The output vector z can be the same as either of the input vectors (x or y).
Usage:
N_VLinearSum(a, x, b, y, z);
-
void N_VConst(sunrealtype c, N_Vector z)
Sets all components of the
N_Vector
z tosunrealtype
c:\[z_i = c, \quad i=0,\ldots,n-1.\]Usage:
N_VConst(c, z);
-
void N_VProd(N_Vector x, N_Vector y, N_Vector z)
Sets the
N_Vector
z to be the component-wise product of theN_Vector
inputs x and y:\[z_i = x_i y_i, \quad i=0,\ldots,n-1.\]Usage:
N_VProd(x, y, z);
-
void N_VDiv(N_Vector x, N_Vector y, N_Vector z)
Sets the
N_Vector
z to be the component-wise ratio of theN_Vector
inputs x and y:\[z_i = \frac{x_i}{y_i}, \quad i=0,\ldots,n-1.\]The \(y_i\) may not be tested for 0 values. It should only be called with a y that is guaranteed to have all nonzero components.
Usage:
N_VDiv(x, y, z);
-
void N_VScale(sunrealtype c, N_Vector x, N_Vector z)
Scales the
N_Vector
x by thesunrealtype
scalar c and returns the result in z:\[z_i = c x_i, \quad i=0,\ldots,n-1.\]Usage:
N_VScale(c, x, z);
-
void N_VAbs(N_Vector x, N_Vector z)
Sets the components of the
N_Vector
z to be the absolute values of the components of theN_Vector
x:\[z_i = |x_i|, \quad i=0,\ldots,n-1.\]Usage:
N_VAbs(x, z);
-
void N_VInv(N_Vector x, N_Vector z)
Sets the components of the
N_Vector
z to be the inverses of the components of theN_Vector
x:\[z_i = \frac{1}{x_i}, \quad i=0,\ldots,n-1.\]This routine may not check for division by 0. It should be called only with an x which is guaranteed to have all nonzero components.
Usage:
N_VInv(x, z);
-
void N_VAddConst(N_Vector x, sunrealtype b, N_Vector z)
Adds the
sunrealtype
scalar b to all components of x and returns the result in theN_Vector
z:\[z_i = x_i+b, \quad i=0,\ldots,n-1.\]Usage:
N_VAddConst(x, b, z);
-
sunrealtype N_VDotProd(N_Vector x, N_Vector z)
Returns the value of the dot-product of the vectors x and y:
\[d = \sum_{i=0}^{n-1} x_i y_i.\]Usage:
d = N_VDotProd(x, y);
-
sunrealtype N_VMaxNorm(N_Vector x)
Returns the value of the \(l_{\infty}\) norm of the
N_Vector
x:\[m = \max_{0\le i< n} |x_i|.\]Usage:
m = N_VMaxNorm(x);
-
sunrealtype N_VWrmsNorm(N_Vector x, N_Vector w)
Returns the weighted root-mean-square norm of the
N_Vector
x with (positive)sunrealtype
weight vector w:\[m = \sqrt{\left( \sum_{i=0}^{n-1} (x_i w_i)^2 \right) / n}\]Usage:
m = N_VWrmsNorm(x, w);
-
sunrealtype N_VWrmsNormMask(N_Vector x, N_Vector w, N_Vector id)
Returns the weighted root mean square norm of the
N_Vector
x withsunrealtype
weight vector w built using only the elements of x corresponding to positive elements of theN_Vector
id:\[m = \sqrt{\left( \sum_{i=0}^{n-1} (x_i w_i H(id_i))^2 \right) / n},\]where \(H(\alpha)=\begin{cases} 1 & \alpha>0\\ 0 & \alpha \leq 0\end{cases}\).
Usage:
m = N_VWrmsNormMask(x, w, id);
-
sunrealtype N_VMin(N_Vector x)
Returns the smallest element of the
N_Vector
x:\[m = \min_{0\le i< n} x_i.\]Usage:
m = N_VMin(x);
-
sunrealtype N_VWL2Norm(N_Vector x, N_Vector w)
Returns the weighted Euclidean \(l_2\) norm of the
N_Vector
x withsunrealtype
weight vector w:\[m = \sqrt{\sum_{i=0}^{n-1}\left(x_i w_i\right)^2}.\]Usage:
m = N_VWL2Norm(x, w);
-
sunrealtype N_VL1Norm(N_Vector x)
Returns the \(l_1\) norm of the
N_Vector
x:\[m = \sum_{i=0}^{n-1} |x_i|.\]Usage:
m = N_VL1Norm(x);
-
void N_VCompare(sunrealtype c, N_Vector x, N_Vector z)
Compares the components of the
N_Vector
x to thesunrealtype
scalar c and returns anN_Vector
z such that for all \(0\le i< n\),\[\begin{split}z_i = \begin{cases} 1.0 &\quad\text{if}\; |x_i| \ge c,\\ 0.0 &\quad\text{otherwise}\end{cases}.\end{split}\]Usage:
N_VCompare(c, x, z);
-
sunbooleantype N_VInvTest(N_Vector x, N_Vector z)
Sets the components of the
N_Vector
z to be the inverses of the components of theN_Vector
x, with prior testing for zero values:\[z_i = \frac{1}{x_i}, \quad i=0,\ldots,n-1.\]This routine returns a boolean assigned to
SUNTRUE
if all components of x are nonzero (successful inversion) and returnsSUNFALSE
otherwise.Usage:
t = N_VInvTest(x, z);
-
sunbooleantype N_VConstrMask(N_Vector c, N_Vector x, N_Vector m)
Performs the following constraint tests based on the values in \(c_i\):
\[\begin{split}\begin{array}{rllll} x_i &>& 0 \;&\text{if}\; &c_i = 2, \\ x_i &\ge& 0 \;&\text{if}\; &c_i = 1, \\ x_i &<& 0 \;&\text{if}\; &c_i = -2, \\ x_i &\le& 0 \;&\text{if}\; &c_i = -1. \end{array}\end{split}\]There is no constraint on \(x_i\) if \(c_i = 0\). This routine returns a boolean assigned to
SUNFALSE
if any element failed the constraint test and assigned toSUNTRUE
if all passed. It also sets a mask vector m, with elements equal to 1.0 where the constraint test failed, and 0.0 where the test passed. This routine is used only for constraint checking.Usage:
t = N_VConstrMask(c, x, m);
-
sunrealtype N_VMinQuotient(N_Vector num, N_Vector denom)
This routine returns the minimum of the quotients obtained by termwise dividing the elements of n by the elements in d:
\[\min_{0\le i< n} \frac{\text{num}_i}{\text{denom}_i}.\]A zero element in denom will be skipped. If no such quotients are found, then the large value
SUN_BIG_REAL
(defined in the header filesundials_types.h
) is returned.Usage:
minq = N_VMinQuotient(num, denom);
8.2.2. Fused operations
The following fused vector operations are optional. These
operations are intended to increase data reuse, reduce parallel
communication on distributed memory systems, and lower the number of
kernel launches on systems with accelerators. If a particular NVECTOR
implementation defines one of the fused vector operations as
NULL
, the NVECTOR interface will call one of the above standard
vector operations as necessary. As above, for each operation, we give
the name, usage of the function, and a description of its mathematical
operations below.
-
SUNErrCode N_VLinearCombination(int nv, sunrealtype *c, N_Vector *X, N_Vector z)
This routine computes the linear combination of nv vectors with \(n\) elements:
\[z_i = \sum_{j=0}^{nv-1} c_j x_{j,i}, \quad i=0,\ldots,n-1,\]where \(c\) is an array of \(nv\) scalars, \(x_j\) is a vector in the vector array X, and z is the output vector. If the output vector z is one of the vectors in X, then it must be the first vector in the vector array. The operation returns a
SUNErrCode
.Usage:
retval = N_VLinearCombination(nv, c, X, z);
-
SUNErrCode N_VScaleAddMulti(int nv, sunrealtype *c, N_Vector x, N_Vector *Y, N_Vector *Z)
This routine scales and adds one vector to nv vectors with \(n\) elements:
\[z_{j,i} = c_j x_i + y_{j,i}, \quad j=0,\ldots,nv-1 \quad i=0,\ldots,n-1,\]where c is an array of scalars, x is a vector, \(y_j\) is a vector in the vector array Y, and \(z_j\) is an output vector in the vector array Z. The operation returns a
SUNErrCode
.Usage:
retval = N_VScaleAddMulti(nv, c, x, Y, Z);
-
SUNErrCode N_VDotProdMulti(int nv, N_Vector x, N_Vector *Y, sunrealtype *d)
This routine computes the dot product of a vector with nv vectors having \(n\) elements:
\[d_j = \sum_{i=0}^{n-1} x_i y_{j,i}, \quad j=0,\ldots,nv-1,\]where d is an array of scalars containing the computed dot products, x is a vector, and \(y_j\) is a vector the vector array Y. The operation returns a
SUNErrCode
.Usage:
retval = N_VDotProdMulti(nv, x, Y, d);
8.2.3. Vector array operations
The following vector array operations are also optional. As with the
fused vector operations, these are intended to increase data reuse,
reduce parallel communication on distributed memory systems, and lower
the number of kernel launches on systems with accelerators. If a
particular NVECTOR implementation defines one of the fused or vector
array operations as NULL
, the NVECTOR interface will call one of
the above standard vector operations as necessary. As above, for each
operation, we give the name, usage of the function, and a description
of its mathematical operations below.
-
SUNErrCode N_VLinearSumVectorArray(int nv, sunrealtype a, N_Vector *X, sunrealtype b, N_Vector *Y, N_Vector *Z)
This routine computes the linear sum of two vector arrays of nv vectors with \(n\) elements:
\[z_{j,i} = a x_{j,i} + b y_{j,i}, \quad i=0,\ldots,n-1 \quad j=0,\ldots,nv-1,\]where a and b are scalars, \(x_j\) and \(y_j\) are vectors in the vector arrays X and Y respectively, and \(z_j\) is a vector in the output vector array Z. The operation returns a
SUNErrCode
.Usage:
retval = N_VLinearSumVectorArray(nv, a, X, b, Y, Z);
-
SUNErrCode N_VScaleVectorArray(int nv, sunrealtype *c, N_Vector *X, N_Vector *Z)
This routine scales each element in a vector of \(n\) elements in a vector array of nv vectors by a potentially different constant:
\[z_{j,i} = c_j x_{j,i}, \quad i=0,\ldots,n-1 \quad j=0,\ldots,nv-1,\]where c is an array of scalars, \(x_j\) is a vector in the vector array X, and \(z_j\) is a vector in the output vector array Z. The operation returns a
SUNErrCode
.Usage:
retval = N_VScaleVectorArray(nv, c, X, Z);
-
SUNErrCode N_VConstVectorArray(int nv, sunrealtype c, N_Vector *Z)
This routine sets each element in a vector of \(n\) elements in a vector array of nv vectors to the same value:
\[z_{j,i} = c, \quad i=0,\ldots,n-1 \quad j=0,\ldots,nv-1,\]where c is a scalar and \(z_j\) is a vector in the vector array Z. The operation returns a
SUNErrCode
.Usage:
retval = N_VConstVectorArray(nv, c, Z);
-
SUNErrCode N_VWrmsNormVectorArray(int nv, N_Vector *X, N_Vector *W, sunrealtype *m)
This routine computes the weighted root mean square norm of each vector in a vector array:
\[m_j = \left( \frac1n \sum_{i=0}^{n-1} \left(x_{j,i} w_{j,i}\right)^2\right)^{1/2}, \quad j=0,\ldots,nv-1,\]where \(x_j\) is a vector in the vector array X, \(w_j\) is a weight vector in the vector array W, and m is the output array of scalars containing the computed norms. The operation returns a
SUNErrCode
.Usage:
retval = N_VWrmsNormVectorArray(nv, X, W, m);
-
SUNErrCode N_VWrmsNormMaskVectorArray(int nv, N_Vector *X, N_Vector *W, N_Vector id, sunrealtype *m)
This routine computes the masked weighted root mean square norm of each vector in a vector array:
\[m_j = \left( \frac1n \sum_{i=0}^{n-1} \left(x_{j,i} w_{j,i} H(id_i)\right)^2 \right)^{1/2}, \quad j=0,\ldots,nv-1,\]where \(H(id_i)=1\) if \(id_i > 0\) and is zero otherwise, \(x_j\) is a vector in the vector array X, \(w_j\) is a weight vector in the vector array W, id is the mask vector, and m is the output array of scalars containing the computed norms. The operation returns a
SUNErrCode
.Usage:
retval = N_VWrmsNormMaskVectorArray(nv, X, W, id, m);
-
SUNErrCode N_VScaleAddMultiVectorArray(int nv, int nsum, sunrealtype *c, N_Vector *X, N_Vector **YY, N_Vector **ZZ)
This routine scales and adds a vector array of nv vectors to nsum other vector arrays:
\[z_{k,j,i} = c_k x_{j,i} + y_{k,j,i}, \quad i=0,\ldots,n-1 \quad j=0,\ldots,nv-1, \quad k=0,\ldots,nsum-1\]where c is an array of scalars, \(x_j\) is a vector in the vector array X, \(y_{k,j}\) is a vector in the array of vector arrays YY, and \(z_{k,j}\) is an output vector in the array of vector arrays ZZ. The operation returns a
SUNErrCode
.Usage:
retval = N_VScaleAddMultiVectorArray(nv, nsum, c, x, YY, ZZ);
-
SUNErrCode N_VLinearCombinationVectorArray(int nv, int nsum, sunrealtype *c, N_Vector **XX, N_Vector *Z)
This routine computes the linear combination of nsum vector arrays containing nv vectors:
\[z_{j,i} = \sum_{k=0}^{nsum-1} c_k x_{k,j,i}, \quad i=0,\ldots,n-1 \quad j=0,\ldots,nv-1,\]where c is an array of scalars, \(x_{k,j}\) is a vector in array of vector arrays XX, and \(z_{j,i}\) is an output vector in the vector array Z. If the output vector array is one of the vector arrays in XX, it must be the first vector array in XX. The operation returns a
SUNErrCode
.Usage:
retval = N_VLinearCombinationVectorArray(nv, nsum, c, XX, Z);
8.2.4. Local reduction operations
The following local reduction operations are also optional. As with
the fused and vector array operations, these are intended to reduce
parallel communication on distributed memory systems. If a particular
NVECTOR implementation defines one of the local reduction operations
as NULL
, the NVECTOR interface will call one of the above standard
vector operations as necessary. As above, for each operation, we give
the name, usage of the function, and a description of its mathematical
operations below.
-
sunrealtype N_VDotProdLocal(N_Vector x, N_Vector y)
This routine computes the MPI task-local portion of the ordinary dot product of x and y:
\[d=\sum_{i=0}^{n_{local}-1} x_i y_i,\]where \(n_{local}\) corresponds to the number of components in the vector on this MPI task (or \(n_{local}=n\) for MPI-unaware applications).
Usage:
d = N_VDotProdLocal(x, y);
-
sunrealtype N_VMaxNormLocal(N_Vector x)
This routine computes the MPI task-local portion of the maximum norm of the NVECTOR x:
\[m = \max_{0\le i< n_{local}} | x_i |,\]where \(n_{local}\) corresponds to the number of components in the vector on this MPI task (or \(n_{local}=n\) for MPI-unaware applications).
Usage:
m = N_VMaxNormLocal(x);
-
sunrealtype N_VMinLocal(N_Vector x)
This routine computes the smallest element of the MPI task-local portion of the NVECTOR x:
\[m = \min_{0\le i< n_{local}} x_i,\]where \(n_{local}\) corresponds to the number of components in the vector on this MPI task (or \(n_{local}=n\) for MPI-unaware applications).
Usage:
m = N_VMinLocal(x);
-
sunrealtype N_VL1NormLocal(N_Vector x)
This routine computes the MPI task-local portion of the \(l_1\) norm of the
N_Vector
x:\[n = \sum_{i=0}^{n_{local}-1} | x_i |,\]where \(n_{local}\) corresponds to the number of components in the vector on this MPI task (or \(n_{local}=n\) for MPI-unaware applications).
Usage:
n = N_VL1NormLocal(x);
-
sunrealtype N_VWSqrSumLocal(N_Vector x, N_Vector w)
This routine computes the MPI task-local portion of the weighted squared sum of the NVECTOR x with weight vector w:
\[s = \sum_{i=0}^{n_{local}-1} (x_i w_i)^2,\]where \(n_{local}\) corresponds to the number of components in the vector on this MPI task (or \(n_{local}=n\) for MPI-unaware applications).
Usage:
s = N_VWSqrSumLocal(x, w);
-
sunrealtype N_VWSqrSumMaskLocal(N_Vector x, N_Vector w, N_Vector id)
This routine computes the MPI task-local portion of the weighted squared sum of the NVECTOR x with weight vector w built using only the elements of x corresponding to positive elements of the NVECTOR id:
\[m = \sum_{i=0}^{n_{local}-1} (x_i w_i H(id_i))^2,\]where
\[\begin{split}H(\alpha) = \begin{cases} 1 & \alpha > 0 \\ 0 & \alpha \leq 0 \end{cases}\end{split}\]and \(n_{local}\) corresponds to the number of components in the vector on this MPI task (or \(n_{local}=n\) for MPI-unaware applications).
Usage:
s = N_VWSqrSumMaskLocal(x, w, id);
-
sunbooleantype N_VInvTestLocal(N_Vector x)
This routine sets the MPI task-local components of the NVECTOR z to be the inverses of the components of the NVECTOR x, with prior testing for zero values:
\[z_i = \frac{1}{x_i}, \: i=0,\ldots,n_{local}-1\]where \(n_{local}\) corresponds to the number of components in the vector on this MPI task (or \(n_{local}=n\) for MPI-unaware applications). This routine returns a boolean assigned to
SUNTRUE
if all task-local components of x are nonzero (successful inversion) and returnsSUNFALSE
otherwise.Usage:
t = N_VInvTestLocal(x);
-
sunbooleantype N_VConstrMaskLocal(N_Vector c, N_Vector x, N_Vector m)
Performs the following constraint tests based on the values in \(c_i\):
\[\begin{split}\begin{array}{rllll} x_i &>& 0 \;&\text{if}\; &c_i = 2, \\ x_i &\ge& 0 \;&\text{if}\; &c_i = 1, \\ x_i &<& 0 \;&\text{if}\; &c_i = -2, \\ x_i &\le& 0 \;&\text{if}\; &c_i = -1. \end{array}\end{split}\]for all MPI task-local components of the vectors. This routine returns a boolean assigned to
SUNFALSE
if any task-local element failed the constraint test and assigned toSUNTRUE
if all passed. It also sets a mask vector m, with elements equal to 1.0 where the constraint test failed, and 0.0 where the test passed. This routine is used only for constraint checking.Usage:
t = N_VConstrMaskLocal(c, x, m);
-
sunrealtype N_VMinQuotientLocal(N_Vector num, N_Vector denom)
This routine returns the minimum of the quotients obtained by term-wise dividing \(num_i\) by \(denom_i\), for all MPI task-local components of the vectors. A zero element in denom will be skipped. If no such quotients are found, then the large value
SUN_BIG_REAL
(defined in the header filesundials_types.h
) is returned.Usage:
minq = N_VMinQuotientLocal(num, denom);
8.2.5. Single Buffer Reduction Operations
The following optional operations are used to combine separate reductions into
a single MPI call by splitting the local computation and communication into
separate functions. These operations are used in low-synchronization
orthogonalization methods to reduce the number of MPI Allreduce
calls. If a
particular NVECTOR implementation does not define these operations additional
communication will be required.
-
SUNErrCode N_VDotProdMultiLocal(int nv, N_Vector x, N_Vector *Y, sunrealtype *d)
This routine computes the MPI task-local portion of the dot product of a vector \(x\) with nv vectors \(y_j\):
\[d_j = \sum_{i=0}^{n_{local}-1} x_i y_{j,i}, \quad j=0,\ldots,nv-1,\]where \(d\) is an array of scalars containing the computed dot products, \(x\) is a vector, \(y_j\) is a vector in the vector array Y, and \(n_{local}\) corresponds to the number of components in the vector on this MPI task. The operation returns a
SUNErrCode
.Usage:
retval = N_VDotProdMultiLocal(nv, x, Y, d);
-
SUNErrCode N_VDotProdMultiAllReduce(int nv, N_Vector x, sunrealtype *d)
This routine combines the MPI task-local portions of the dot product of a vector \(x\) with nv vectors:
retval = MPI_Allreduce(MPI_IN_PLACE, d, nv, MPI_SUNREALTYPE, MPI_SUM, comm)
where d is an array of nv scalars containing the local contributions to the dot product and comm is the MPI communicator associated with the vector x. The operation returns a
SUNErrCode
.Usage:
retval = N_VDotProdMultiAllReduce(nv, x, d);
8.2.6. Exchange operations
The following vector exchange operations are also optional and are
intended only for use when interfacing with the XBraid library for
parallel-in-time integration. In that setting these operations are
required but are otherwise unused by SUNDIALS packages and may be set
to NULL
. For each operation, we give the function signature, a
description of the expected behavior, and an example of the function
usage.
-
SUNErrCode N_VBufSize(N_Vector x, sunindextype *size)
This routine returns the buffer size need to exchange in the data in the vector x between computational nodes.
Usage:
flag = N_VBufSize(x, &buf_size)
-
SUNErrCode N_VBufPack(N_Vector x, void *buf)
This routine fills the exchange buffer buf with the vector data in x.
Usage:
flag = N_VBufPack(x, &buf)
-
SUNErrCode N_VBufUnpack(N_Vector x, void *buf)
This routine unpacks the data in the exchange buffer buf into the vector x.
Usage:
flag = N_VBufUnpack(x, buf)
8.2.7. Output operations
The following optional vector operations are for writing vector data either to
stdout
or to a given file.