naginterfaces.library.ode.ivp_​stiff_​sparjac_​setup

naginterfaces.library.ode.ivp_stiff_sparjac_setup(neq, neqmax, jceval, nwkjac, ia, ja, comm, njcpvt, sens, u, eta, lblock, isplit=73)

ivp_stiff_sparjac_setup is a setup function which must be called prior to an integrator in submodule ode, if sparse matrix linear algebra is required.

For full information please refer to the NAG Library document for d02nu

https://www.nag.com/numeric/nl/nagdoc_29.3/flhtml/d02/d02nuf.html

Parameters

neqint

The number of differential equations.

neqmaxint

A bound on the maximum number of differential equations to be solved during the integration.

jcevalstr, length 1

Specifies the technique to be used to compute the Jacobian.

$jceval=\text{'N'}$

The sparsity structure and the value of the Jacobian are to be determined numerically by the integrator.

$jceval=\text{'S'}$

The sparsity structure of the Jacobian is supplied in the arrays $ia$ and $ja$ but its value is to be determined numerically. This is the recommended mode of operation unless it is a simple matter to supply the Jacobian.

$jceval=\text{'A'}$

The Jacobian will be evaluated by calls to $\text{jac}$. The sparsity structure will be estimated by calls to $\text{jac}$; that is, no explicit sparsity structure need be supplied in the arrays $ia$ and $ja$.

$jceval=\text{'F'}$

The sparsity structure of the Jacobian is supplied in $ia$ and $ja$, and its value will be determined by calls to $\text{jac}$. This is the recommended mode of operation if the $\text{jac}$ is simple to form.

$jceval=\text{'D'}$

The default choice is to be made. In this case ‘D’ is interpreted as ‘S’.

If the sparsity structure is supplied in arrays $ia$ and $ja$, any evidence from the numerical or analytical formation of the Jacobian that this structure is not correct, is ignored.

Only the first character of the actual argument $jceval$ is passed to ivp_stiff_sparjac_setup; hence it is permissible for the actual argument to be more descriptive, e.g., ‘Numerical’, ‘Structural’, ‘Analytical’, ‘Full information’ or ‘Default’ in a call to ivp_stiff_sparjac_setup.

If the option $jceval=\text{'N'}$, $\text{'S'}$ or $\text{'D'}$ is used then the actual argument corresponding to $\text{jac}$ in the call to ivp_stiff_exp_sparjac() or ivp_stiff_imp_sparjac() must be specified as None.

If integration is to be performed by reverse communication (ivp_stiff_exp_revcom() or ivp_stiff_imp_revcom()) then $jceval$ should be set to either ‘N’ or ‘A’.

In this case $ia$ and $ja$ are not used and their lengths may be set to $1$.

nwkjacint

The size of the array $\text{wkjac}$, which you are supplying to the integrator, as declared in the (sub)program from which ivp_stiff_sparjac_setup is called.

Suggested value: $nwkjac=4\times neqmax$ if $jceval=\text{'N'}$ or $\text{'A'}$. If $nwkjac$ is less than this estimate, a message is printed on the file object associated with the advisory I/O unit (see FileObjManager), and execution continues.

iaint, array-like, shape $\left(\text{nia}\right)$

If $jceval=\text{'S'}$, $\text{'F'}$ or $\text{'D'}$, $ia$ must contain details of the sparsity pattern to be used for the Jacobian. See $ja$.

$ia$ is not used if $jceval=\text{'N'}$ or $\text{'A'}$.

jaint, array-like, shape $\left(\text{nja}\right)$

If $jceval=\text{'S'}$, $\text{'F'}$ or $\text{'D'}$, $ja$ must contain details of the sparsity pattern to be used for the Jacobian. $ja$ contains the row indices where nonzero elements occur, reading in column-wise order, and $ia$ contains the starting locations in $ja$ of the descriptions of columns $1,2,\dots ,neq$ in that order, with $ia\left[0\right]=1$. Thus for each column index $j=1,2,\dots ,neq$, the values of the row index $i$ in column $j$ where a nonzero element may occur are given by

$$i=ja\left[k\right]$$

where $ia[j-1]\le k<ia\left[j\right]$.

Thus the total number of nonzeros, $\text{nelement}$, must be $ia\left[neq\right]-1$.

For example, for the following matrix

$$\begin{array}{c}\left(\begin{array}{ccccc}x& 0& x& 0& 0\\ 0& x& x& x& 0\\ x& x& x& 0& 0\\ x& 0& 0& x& x\\ 0& 0& 0& x& x\end{array}\right)\end{array}$$

where $x$ represents nonzero elements (13 in all) the arrays $ia$ and $ja$ should be

$$\begin{array}{c}\begin{array}{cccccccccccccc}ia[k-1]& 1& 4& 6& 9& 12& 14& & & & & & & \\ ja[k-1]& 1& 3& 4& 2& 3& 1& 2& 3& 2& 4& 5& 4& 5\end{array}\end{array}$$

$ja$ is not used if $jceval=\text{'N'}$ or $\text{'A'}$.

commdict, communication object, modified in place

Communication structure.

This argument must have been initialized by a prior call to ivp_stiff_bdf() or ivp_stiff_blend().

njcpvtint

The length of the array $comm$[‘jacpvt’], which you are supplying to the integrator, as dimensioned in the sub(program) from which ivp_stiff_sparjac_setup is called.

Suggested value: $njcpvt=20\times neqmax$ if $jceval=\text{'N'}$ or $\text{'A'}$. If $njcpvt$ is less than this estimate, a message is printed on the file object associated with the advisory I/O unit (see FileObjManager), and execution continues.

sensfloat

A threshold argument used to determine whether or not a matrix element is zero; when $sens$ is set to $0.0$ on entry, the function will use $sens=100.0\times \text{machine precision}$. Otherwise the absolute value of $sens$ is used.

ufloat

Should have a value between $0.0$ and $0.9999$. Otherwise a default value of $0.1$ is used. When the sparsity pattern has been evaluated, the first Jacobian computed is decomposed with $u$ governing the choice of pivots; subsequent Jacobian decompositions use the same pattern of decomposition until the sparsity pattern is re-evaluated. When searching a row for a pivot, any element is excluded from the search which is less than $u$ times the largest of those elements in the row available as pivots. Thus decreasing $u$ biases the algorithm towards maintaining sparsity at the expense of numerical stability.

etafloat

A relative pivot threshold, below which on subsequent decompositions (as described under $u$), an internal error is provoked.

$eta>1.0$

No check on pivot size is made.

$eta\le 0.0$

The default value $eta=1.0e-4$ is used.

lblockbool

Indicates if preordering is used before decomposition.

If $lblock=True$, on entry, the Jacobian matrix is preordered to block lower triangular form before a decomposition is performed (this is the recommended mode).

If you know the structure of the Jacobian to be irreducible, that is not permutable to block lower triangular form, you should set $lblock=False$.

For example, a Jacobian arising from using the method of lines for parabolic partial differential equations would normally be irreducible. (See the specification of ivp_stiff_sparjac_diag() for optional output concerning $lblock$.)

isplitint, optional

This argument is used for splitting the integer workspace $comm$[‘jacpvt’] to effect an efficient decomposition. It must satisfy $1\le isplit\le 99$. If $isplit$ lies outside this range on entry, a default value of $73$ is used. An appropriate value for $isplit$ for subsequent runs on similar problems is available via the optional output ivp_stiff_sparjac_diag().

Raises

NagValueError

(errno $1$)

On entry, $nwkjac=⟨value⟩$ and $ia\left[neq\right]-1+2\times neq=⟨value⟩$.

Constraint: $nwkjac\ge ia\left[neq\right]-1+2\times neq$.

(errno $1$)

On entry, $ja$ defines duplicate elements in row $⟨value⟩$ and column $⟨value⟩$.

(errno $1$)

On entry, $i=⟨value⟩$, $ja\left[i\right]=⟨value⟩$ and $neq=⟨value⟩$.

Constraint: $1\le ja\left[i\right]\le neq$ for all $i$.

(errno $1$)

On entry, $i=⟨value⟩$, $ia\left[i\right]=⟨value⟩$, $i-1=⟨value⟩$, $ia[i-1]=⟨value⟩$, $neq=⟨value⟩$.

Constraint: $ia[⟨value⟩]-ia[⟨value⟩]\le neq$.

(errno $1$)

On entry, $ia[⟨value⟩]=⟨value⟩$ and $ia[⟨value⟩]=⟨value⟩$.

Constraint: $ia[⟨value⟩]\ge ia[⟨value⟩]$.

(errno $1$)

On entry, $ia\left[0\right]=⟨value⟩$.

Constraint: $ia\left[0\right]=1$.

(errno $1$)

On entry, $njcpvt=⟨value⟩$ and $3\times \left(ia\right[neq]-1)+14\times neq+1=⟨value⟩$.

Constraint: $njcpvt\ge 3\times \left(ia\right[neq]-1)+14\times neq+1$.

(errno $1$)

On entry, $\text{nja}=⟨value⟩$ and $ia\left[neq\right]-1=⟨value⟩$.

Constraint: $\text{nja}\ge ia\left[neq\right]-1$.

(errno $1$)

On entry, $ia\left[neq\right]-1=⟨value⟩$ and $neq=⟨value⟩$.

Constraint: $ia\left[neq\right]-1\ge neq$.

(errno $1$)

On entry, $\text{nia}=⟨value⟩$ and $neq+1=⟨value⟩$.

Constraint: $\text{nia}\ge neq+1$.

(errno $1$)

On entry, $jceval=⟨value⟩$.

Constraint: $jceval=\text{'A'}$, $\text{'N'}$, $\text{'S'}$, $\text{'F'}$ or $\text{'D'}$.

(errno $1$)

On entry, $neq=⟨value⟩$ and $neqmax=⟨value⟩$.

Constraint: $neq\le neqmax$.

(errno $1$)

On entry, $neqmax=⟨value⟩$.

Constraint: $neqmax\ge 1$.

(errno $1$)

On entry, $neq=⟨value⟩$.

Constraint: $neq\ge 1$.

Notes

No equivalent traditional C interface for this routine exists in the NAG Library.

ivp_stiff_sparjac_setup defines the linear algebra to be used as sparse matrix linear algebra, permits you to specify the method for calculating the Jacobian and its structure, and checks the validity of certain input values.