Preconditioner Routines for Banded Problems

This module provides a general-purpose banded preconditioner for use with the daskr solver, with the Krylov linear system method.

When using daskr to solve a problem $G(t,y,\dot{y}) = 0$, whose Jacobian $J = \partial G/ \partial y + c_J \partial G/ \partial \dot{y}$, where $c_J$ is a scalar, is either banded or approximately equal to a banded matrix, the routines jac_banpre and psol_banpre can be used to generate a banded approximation to $J$ as the preconditioner and to solve the resulting banded linear system, in conjunction with the Krylov method option.

Other than the user-supplied residual routine res defining $G(t,y,\dot{y})$, the only other inputs required by these routines are the half-bandwidth parameters $\mathrm{ml}$ and $\mathrm{mu}$ of the approximate banded Jacobian. If the system size is $\mathrm{neq}$, the half-bandwidths are defined as integers between 0 and $\mathrm{neq} -1$ such that only elements with indices $(i,j)$ satisfying $-\mathrm{ml} \le j, i \le \mathrm{mu}$ are to be retained in the preconditioner. For example, if $\mathrm{ml} = \mathrm{mu} = 0$, a diagonal matrix will be generated as the preconditioner. The banded preconditioner is obtained by difference quotient approximations. If the true problem Jacobian is not banded but is approximately equal to a matrix that is banded, the procedure used here will have the effect of lumping the elements outside of the band onto the elements within the band.

Usage

To use these routines in conjunction with daskr, the user's calling program should include the following, in addition to setting the other daskr input parameters:

• Dimension the array ipar to have length at least 2, and load the half-bandwidths into ipar as ipar(1) = ml and ipar(2) = mu. Array ipar is used to communicate these parameters to jac_banpre and psol_banpre. If the user program also uses ipar for communication with res, that data should be located beyond the first 2 positions.

• Set info(12) = 1 to select the Krylov iterative method and info(15) = 1 to indicate that a jac routine exists. Then in the call to daskr, pass the procedure names jac_banpre and psol_banpre as the arguments jac and psol, respectively.

• The daskr work arrays rwork and iwork must include segments rwp and iwp for use by jac_banpre and psol_banpre. The lengths of these arrays depend on the problem size and half-bandwidths, as as shown in the table below. Note the integer divide in lrwp. Load these lengths in iwork as iwork(27) = lrwp and iwork(28) = liwp, and include these values in the declared size of rwork and iwork, respectively.

Example

The program example_heat demonstrates the use of this preconditioner.