Contents

# ?gerfsx

Uses extra precise iterative refinement to improve the solution to the system of linear equations with a general coefficient matrix A and provides error bounds and backward error estimates.

## Syntax

Include Files
• mkl.h
Description
The routine improves the computed solution to a system of linear equations and provides error bounds and backward error estimates for the solution. In addition to a normwise error bound, the code provides a maximum componentwise error bound, if possible. See comments for
err_bnds_norm
and
err_bnds_comp
for details of the error bounds.
The original system of linear equations may have been equilibrated before calling this routine, as described by the parameters
equed
,
r
, and
c
below. In this case, the solution and error bounds returned are for the original unequilibrated system.
Input Parameters
matrix_layout
Specifies whether two-dimensional array storage is row major (
LAPACK_ROW_MAJOR
) or column major (
LAPACK_COL_MAJOR
).
trans
Must be
'N'
,
'T'
, or
'C'
.
Specifies the form of the system of equations:
If
trans
=
'N'
, the system has the form
A
*
X
=
B
(No transpose).
If
trans
=
'T'
, the system has the form
A
T
*
X
=
B
(Transpose).
If
trans
=
'C'
, the system has the form
A
H
*
X
=
B
(Conjugate transpose for complex flavors, Transpose for real flavors).
equed
Must be
'N'
,
'R'
,
'C'
, or
'B'
.
Specifies the form of equilibration that was done to
A
before calling this routine.
If
equed
=
'N'
, no equilibration was done.
If
equed
=
'R'
, row equilibration was done, that is,
A
has been premultiplied by
diag
(
r
).
If
equed
=
'C'
, column equilibration was done, that is,
A
has been postmultiplied by
diag
(
c
).
If
equed
=
'B'
, both row and column equilibration was done, that is,
A
has been replaced by
diag
(
r
)*
A
*
diag
(
c
)
. The right-hand side
B
has been changed accordingly.
n
The number of linear equations; the order of the matrix
A
;
n
0.
nrhs
The number of right-hand sides; the number of columns of the matrices
B
and
X
;
nrhs
0.
a
,
af
,
b
Arrays:
a
(size max(1,
lda
*
n
))
,
af
(size max(1,
ldaf
*
n
))
,
b
(size max(1,
ldb
*
nrhs
) for column major layout and max(1,
ldb
*
n
) for row major layout)
.
The array
a
contains the original
n
-by-
n
matrix
A
.
The array
af
contains the factored form of the matrix
A
, that is, the factors
L
and
U
from the factorization
A
=
P*L*U
as computed by
?getrf
.
The array
b
contains the matrix
B
whose columns are the right-hand sides for the systems of equations.
lda
a
;
lda
max(1,
n
)
.
ldaf
af
;
ldaf
max(1,
n
)
.
ipiv
Array, size at least
max(1,
n
)
. Contains the pivot indices as computed by
?getrf
; for row
1
i
n
, row
i
of the matrix was interchanged with row
ipiv
(
i
)
.
r
,
c
Arrays:
r
(size
n
),
c
(size
n
). The array
r
contains the row scale factors for
A
, and the array
c
contains the column scale factors for
A
.
equed
=
'R'
or
'B'
,
A
is multiplied on the left by
diag
(
r
); if
equed
=
'N'
or
'C'
,
r
is not accessed.
If
equed
=
'R'
or
'B'
, each element of
r
must be positive.
If
equed
=
'C'
or
'B'
,
A
is multiplied on the right by
diag
(
c
); if
equed
=
'N'
or
'R'
,
c
is not accessed.
If
equed
=
'C'
or
'B'
, each element of
c
must be positive.
Each element of
r
or
c
should be a power of the radix to ensure a reliable solution and error estimates. Scaling by powers of the radix does not cause rounding errors unless the result underflows or overflows. Rounding errors during scaling lead to refining with a matrix that is not equivalent to the input matrix, producing error estimates that may not be reliable.
ldb
The leading dimension of the array
b
;
ldb
max(1,
n
) for column major layout and
ldb
nrhs
for row major layout
.
x
Array,
of size max(1,
ldx
*
nrhs
) for column major layout and max(1,
ldx
*
n
) for row major layout
.
The solution matrix
X
as computed by
?getrs
ldx
The leading dimension of the output array
x
;
ldx
max(1,
n
) for column major layout and
ldx
nrhs
for row major layout
.
n_err_bnds
Number of error bounds to return for each right hand side and each type (normwise or componentwise). See
err_bnds_norm
and
err_bnds_comp
descriptions in
Output Arguments
section below.
nparams
Specifies the number of parameters set in
params
. If
0, the
params
array is never referenced and default values are used.
params
Array, size
nparams
. Specifies algorithm parameters. If an entry is less than 0.0, that entry is filled with the default value used for that parameter. Only positions up to
nparams
are accessed; defaults are used for higher-numbered parameters. If defaults are acceptable, you can pass
nparams
= 0, which prevents the source code from accessing the
params
argument.
params
[0]
: Whether to perform iterative refinement or not. Default: 1.0
=0.0
No refinement is performed and no error bounds are computed.
=1.0
Use the double-precision refinement algorithm, possibly with doubled-single computations if the compilation environment does not support double precision.
(Other values are reserved for future use.)
params
[1]
: Maximum number of residual computations allowed for refinement.
Default
10.0
Aggressive
Set to 100.0 to permit convergence using approximate factorizations or factorizations other than
LU
. If the factorization uses a technique other than Gaussian elimination, the guarantees in
err_bnds_norm
and
err_bnds_comp
may no longer be trustworthy.
params
[2]
: Flag determining if the code will attempt to find a solution with a small componentwise relative error in the double-precision algorithm. Positive is true, 0.0 is false. Default: 1.0 (attempt componentwise convergence).
Output Parameters
x
The improved solution matrix
X
.
rcond
Reciprocal scaled condition number. An estimate of the reciprocal Skeel condition number of the matrix
A
after equilibration (if done). If
rcond
is less than the machine precision, in particular, if
rcond
= 0, the matrix is singular to working precision. Note that the error may still be small even if this number is very small and the matrix appears ill-conditioned.
berr
Array, size at least
max(1,
nrhs
)
. Contains the componentwise relative backward error for each solution vector
x
j
, that is, the smallest relative change in any element of
A
or
B
that makes
x
j
an exact solution.
err_bnds_norm
Array of size
nrhs
*
n_err_bnds
. For each right-hand side, contains information about various error bounds and condition numbers corresponding to the normwise relative error
, which is defined as follows:
Normwise relative error in the
i
-th solution vector
The array is indexed by the type of error information as described below. There are currently up to three pieces of information returned.
err
=1
"Trust/don't trust" boolean. Trust the answer if the reciprocal condition number is less than the threshold
sqrt
(
n
)*
slamch
(
ε
)
for single precision flavors and
sqrt
(
n
)*
dlamch
(
ε
)
for double precision flavors.
err
=2
"Guaranteed" error bound. The estimated forward error, almost certainly within a factor of 10 of the true error so long as the next entry is greater than the threshold
sqrt
(
n
)*
slamch
(
ε
)
for single precision flavors and
sqrt
(
n
)*
dlamch
(
ε
)
for double precision flavors. This error bound should only be trusted if the previous boolean is true.
err
=3
Reciprocal condition number. Estimated normwise reciprocal condition number. Compared with the threshold
sqrt
(
n
)*
slamch
(
ε
)
for single precision flavors and
sqrt
(
n
)*
dlamch
(
ε
)
for double precision flavors to determine if the error estimate is "guaranteed". These reciprocal condition numbers for some appropriately scaled matrix
Z
are:
Let
z
=
s
*
a
, where
s
scales each row by a power of the radix so all absolute row sums of
z
are approximately 1.
The information for right-hand side
i
, where 1
i
nrhs
, and type of error
err
is stored in:
• Column major layout:
err_bnds_norm
[(
err
- 1)*
nrhs
+
i
- 1]
.
• Row major layout:
err_bnds_norm
[
err
- 1 + (
i
- 1)*
n_err_bnds
]
err_bnds_comp
Array of size
nrhs
*
n_err_bnds
. For each right-hand side, contains information about various error bounds and condition numbers corresponding to the componentwise relative error
, which is defined as follows:
Componentwise relative error in the
i
-th solution vector:
The array is indexed by the type of error information as described below. There are currently up to three pieces of information returned for each right-hand side. If componentwise accuracy is not requested (
params[2]
= 0.0), then
err_bnds_comp
is not accessed. If
n_err_bnds
< 3, then at most the first
n_err_bnds
columns of the
err_bnds_comp
array are returned.
err
=1
"Trust/don't trust" boolean. Trust the answer if the reciprocal condition number is less than the threshold
sqrt
(
n
)*
slamch
(