CXML
dsskyr
Symmetric sparse iterative refinement using skyline storage scheme
FORMAT
DSSKYR
(n, au, auf, iaudiag, nau, b, ldb, x, ldx, ferr, berr, nbx, iparam,
rparam, iwrk, rwrk, ierror)
Arguments
n integer*4
On entry, the order of the matrix A.
On exit, n is unchanged.
au real*8
On entry, an array of length at least nau, containing
the matrix A stored in the skyline storage scheme,
using either the profile-in or the diagonal-out storage
mode.
On exit, au is unchanged.
auf real*8
On entry, an array of length at least nau, containing
transp(U)*D*U factorization of the the matrix A stored
in the skyline storage scheme, using either the
profile-in or the diagonal-out storage mode. The
factorization has been obtained by a prior call to the
routine DSSKYF. auf must remain unchanged between
calls to the routines DSSKYF and DSSKYR.
On exit, auf is unchanged.
iaudiag integer*4
On entry, an array of length at least n for the
profile-in storage mode and (n+1) for the diagonal-out
storage mode, containing the pointers to the locations
of the diagonal elements in arrays AU and AUF.
On exit, iaudiag is unchanged.
nau integer*4
On entry, the number of elements in array AU. nau is
also the envelope size of the symmetric part of the
matrix A. For the profile-in storage mode, nau =
IAUDIAG(n). For the diagonal-out storage mode, nau =
IAUDIAG(n+1) - 1.
On exit, nau is unchanged.
b real*8
On entry, a two dimensional array B of order ldb by at
least nbx, containing the nbx right sides.
On exit, b is unchanged.
ldb integer*4
On entry, the leading dimension of array B. ldb >=n.
On exit, ldb is unchanged.
x real*8
On entry, a two dimensional array X of order ldx by at
least nbx, containing the nbx solution vectors obtained
after a call to the routine DSSKYS.
On exit, x contains the improved solutions obtained
after iterative refinement.
ldx integer*4
On entry, the leading dimension of array X. ldx >=n.
On exit, ldx is unchanged.
ferr real*8
On entry, an array FERR of length at least nbx, whose
elements are unspecified variables.
On exit, ferr contains the estimated error bounds for
each of the nbx solution vectors.
berr
On entry, an array BERR of length at least nbz, whose
elements are unspecified variables.
On exit, berr contains the component-wise relative
backward error for each of the nbz solution vectors.
nbx
On entry, the number of right sides.
On exit, nbz is unchanged.
iparam integer*4
An array of length at least 100, containing the integer
parameters for the iterative refinement and error
bounds calculation.
iparam(1): niparam
On entry, defines the length of the array IPARAM.
niparam >= 100.
On exit, iparam(1) is unchanged.
iparam(2): nrparam
On entry, defines the length of the array RPARAM. As
the real parameter array is not used at present,
nrparam can be unspecified.
On exit, iparam(2) is unchanged.
iparam(3): niwrk
On entry, defines the size of the integer work array,
IWRK. niwrk >= 3n.
On exit, iparam(3) is unchanged.
iparam(4): nrwrk
On entry, defines the size of the real work array,
RWRK. nrwrk >= 3n.
On exit, iparam(4) is unchanged.
iparam(5): iounit
On entry, defines the I/O unit number for printing
error messages and information from the routine DSSKYR.
The I/O unit must be opened in the calling subprogram.
If iounit <= 0, no output is generated.
On exit, iparam(5) is unchanged.
iparam(6): iolevel
On entry, defines the message level that determines the
amount of information printed out to iounit, when
iounit > 0.
iolevel = 0 : fatal error messages only
iolevel = 1 : error messages and minimal information
iolevel = 2 : error messages and detailed information
On exit, iparam(6) is unchanged.
iparam(7): idefault
On entry, defines if the default values should be used
in arrays IPARAM and RPARAM. If idefault = 0, then the
following default values are assigned:
IPARAM(1) = niparam = 100
IPARAM(6) = iolevel = 0
IPARAM(8) = istore = 1
IPARAM(9) = itmax = 5
If idefault = 1, then you must assign values to the
above variables before the call to the DSSKYR routine.
On exit, iparam(7) is unchanged.
iparam(8): istore
On entry, defines the type of storage scheme used for
the skyline matrix. If istore = 1, the matrix A is
stored using the profile-in storage mode; if istore =
2, the matrix A is stored using the diagonal-out
storage mode. Default: istore = 1.
On exit, iparam(8) is unchanged.
iparam(9): itmax
On entry, defines the maximum number of iterations for
the iterative refinement process. Default: itmax = 5.
On exit, iparam(9) is unchanged.
rparam real*8
An array of length at least 100, containing the real
parameters for the iterative refinement and error
bounds calculation.
On exit, rparam is unchanged. rparam is not used by
the routine DSSKYR at present, but is reserved for
future use. It can be a dummy variable.
iwrk integer*4
On entry, an array of length at least 3n used for
integer workspace. The first 2n elements of the array
IWRK, generated by the routine DSSKYF, should be passed
unchanged to the routine DSSKYR.
On exit, the first 2n elements of iwrk are unchanged.
rwrk real*8
On entry, an array of length at least 3n used for real
workspace.
On exit, the first 3n elements of rwrk are
overwritten.
ierror integer*4
On entry, an unspecified variable.
On exit, ierror contains the error flag. A value of
zero indicates a normal exit from the routine DSSKYR.
Description
DSSKYR obtains an improved solution to the system
A X = B
via iterative refinement. This is done by calculating the matrix of
residuals R using the matrix of solutions X_hat obtained from DSSKYS, and
obtaining a new matrix of solutions X(new) as follows:
R = B - A * X_hat
delta_X = inverse(A) * R
and
X(new) = X_hat + delta_X
The process of iterative refinement therefore requires both the original
matrix A as well as the transp(U)*D*U factorization obtained via the
routine DSSKYF. Since this routine overwrites the matrix A by the
factorization, a copy of the matrix must be made prior to the call to
DSSKYF. Further, both the right sides B and the solution vectors X are
required during iterative refinement. Since the solution process in the
routine DSSKYS overwrites the right sides with the solution vectors, a copy
of the right sides must be made prior to the call to the routine DSSKYS.
In addition to the iterative refinement of the solution vectors, the
routine DSSKYR also provides the component-wise relative backward error,
berr and the estimated forward error bound, ferr, for each solution vector
[Arioli, Demmel, Duff 1989, Anderson et. al. 1992]. berr is the smallest
relative change in any entry of A or B that makes x_hat an exact solution.
ferr bounds the magnitude of the largest entry in x_hat - x(true) divided
by the magnitude of the largest entry in x_hat.
The process of iterative refinement is continued as long as all of the
following conditions are satisfied [Arioli, Demmel, Duff 1989]:
• The number of iterations of the iterative refinement process is less
than IPARAM(9) = itmax.
• berr reduces by at least a factor of 2 during the previous iteration.
• berr is larger than the machine precision.
The routine DSSKYR is called after a call to the routine DSSKYF to obtain
the transp(U)*D*U factorization and a call to the routine DSSKYS to obtain
the solution X. The first 2n elements of the integer workspace array IWRK,
generated by DSSKYF, contain information for use by DSSKYR and therefore
must remain unchanged between the calls to the routines DSSKYF and DSSKYR.
The real work array, RWRK, is not used at present. The storage scheme used
in the routines DSSKYF, DSSKYS, and DSSKYR must be identical.
CXML Home Page Index of CXML Routines