version 1.1, 2010/08/07 13:21:04
|
version 1.14, 2018/05/29 06:55:17
|
Line 1
|
Line 1
|
|
*> \brief \b DLA_GBRFSX_EXTENDED improves the computed solution to a system of linear equations for general banded matrices by performing extra-precise iterative refinement and provides error bounds and backward error estimates for the solution. |
|
* |
|
* =========== DOCUMENTATION =========== |
|
* |
|
* Online html documentation available at |
|
* http://www.netlib.org/lapack/explore-html/ |
|
* |
|
*> \htmlonly |
|
*> Download DLA_GBRFSX_EXTENDED + dependencies |
|
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dla_gbrfsx_extended.f"> |
|
*> [TGZ]</a> |
|
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dla_gbrfsx_extended.f"> |
|
*> [ZIP]</a> |
|
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dla_gbrfsx_extended.f"> |
|
*> [TXT]</a> |
|
*> \endhtmlonly |
|
* |
|
* Definition: |
|
* =========== |
|
* |
|
* SUBROUTINE DLA_GBRFSX_EXTENDED( PREC_TYPE, TRANS_TYPE, N, KL, KU, |
|
* NRHS, AB, LDAB, AFB, LDAFB, IPIV, |
|
* COLEQU, C, B, LDB, Y, LDY, |
|
* BERR_OUT, N_NORMS, ERR_BNDS_NORM, |
|
* ERR_BNDS_COMP, RES, AYB, DY, |
|
* Y_TAIL, RCOND, ITHRESH, RTHRESH, |
|
* DZ_UB, IGNORE_CWISE, INFO ) |
|
* |
|
* .. Scalar Arguments .. |
|
* INTEGER INFO, LDAB, LDAFB, LDB, LDY, N, KL, KU, NRHS, |
|
* $ PREC_TYPE, TRANS_TYPE, N_NORMS, ITHRESH |
|
* LOGICAL COLEQU, IGNORE_CWISE |
|
* DOUBLE PRECISION RTHRESH, DZ_UB |
|
* .. |
|
* .. Array Arguments .. |
|
* INTEGER IPIV( * ) |
|
* DOUBLE PRECISION AB( LDAB, * ), AFB( LDAFB, * ), B( LDB, * ), |
|
* $ Y( LDY, * ), RES(*), DY(*), Y_TAIL(*) |
|
* DOUBLE PRECISION C( * ), AYB(*), RCOND, BERR_OUT(*), |
|
* $ ERR_BNDS_NORM( NRHS, * ), |
|
* $ ERR_BNDS_COMP( NRHS, * ) |
|
* .. |
|
* |
|
* |
|
*> \par Purpose: |
|
* ============= |
|
*> |
|
*> \verbatim |
|
*> |
|
*> |
|
*> DLA_GBRFSX_EXTENDED improves the computed solution to a system of |
|
*> linear equations by performing extra-precise iterative refinement |
|
*> and provides error bounds and backward error estimates for the solution. |
|
*> This subroutine is called by DGBRFSX to perform iterative refinement. |
|
*> In addition to normwise error bound, the code provides maximum |
|
*> componentwise error bound if possible. See comments for ERR_BNDS_NORM |
|
*> and ERR_BNDS_COMP for details of the error bounds. Note that this |
|
*> subroutine is only resonsible for setting the second fields of |
|
*> ERR_BNDS_NORM and ERR_BNDS_COMP. |
|
*> \endverbatim |
|
* |
|
* Arguments: |
|
* ========== |
|
* |
|
*> \param[in] PREC_TYPE |
|
*> \verbatim |
|
*> PREC_TYPE is INTEGER |
|
*> Specifies the intermediate precision to be used in refinement. |
|
*> The value is defined by ILAPREC(P) where P is a CHARACTER and |
|
*> P = 'S': Single |
|
*> = 'D': Double |
|
*> = 'I': Indigenous |
|
*> = 'X', 'E': Extra |
|
*> \endverbatim |
|
*> |
|
*> \param[in] TRANS_TYPE |
|
*> \verbatim |
|
*> TRANS_TYPE is INTEGER |
|
*> Specifies the transposition operation on A. |
|
*> The value is defined by ILATRANS(T) where T is a CHARACTER and |
|
*> T = 'N': No transpose |
|
*> = 'T': Transpose |
|
*> = 'C': Conjugate transpose |
|
*> \endverbatim |
|
*> |
|
*> \param[in] N |
|
*> \verbatim |
|
*> N is INTEGER |
|
*> The number of linear equations, i.e., the order of the |
|
*> matrix A. N >= 0. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] KL |
|
*> \verbatim |
|
*> KL is INTEGER |
|
*> The number of subdiagonals within the band of A. KL >= 0. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] KU |
|
*> \verbatim |
|
*> KU is INTEGER |
|
*> The number of superdiagonals within the band of A. KU >= 0 |
|
*> \endverbatim |
|
*> |
|
*> \param[in] NRHS |
|
*> \verbatim |
|
*> NRHS is INTEGER |
|
*> The number of right-hand-sides, i.e., the number of columns of the |
|
*> matrix B. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] AB |
|
*> \verbatim |
|
*> AB is DOUBLE PRECISION array, dimension (LDAB,N) |
|
*> On entry, the N-by-N matrix AB. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] LDAB |
|
*> \verbatim |
|
*> LDAB is INTEGER |
|
*> The leading dimension of the array AB. LDBA >= max(1,N). |
|
*> \endverbatim |
|
*> |
|
*> \param[in] AFB |
|
*> \verbatim |
|
*> AFB is DOUBLE PRECISION array, dimension (LDAFB,N) |
|
*> The factors L and U from the factorization |
|
*> A = P*L*U as computed by DGBTRF. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] LDAFB |
|
*> \verbatim |
|
*> LDAFB is INTEGER |
|
*> The leading dimension of the array AF. LDAFB >= max(1,N). |
|
*> \endverbatim |
|
*> |
|
*> \param[in] IPIV |
|
*> \verbatim |
|
*> IPIV is INTEGER array, dimension (N) |
|
*> The pivot indices from the factorization A = P*L*U |
|
*> as computed by DGBTRF; row i of the matrix was interchanged |
|
*> with row IPIV(i). |
|
*> \endverbatim |
|
*> |
|
*> \param[in] COLEQU |
|
*> \verbatim |
|
*> COLEQU is LOGICAL |
|
*> If .TRUE. then column equilibration was done to A before calling |
|
*> this routine. This is needed to compute the solution and error |
|
*> bounds correctly. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] C |
|
*> \verbatim |
|
*> C is DOUBLE PRECISION array, dimension (N) |
|
*> The column scale factors for A. If COLEQU = .FALSE., C |
|
*> is not accessed. If C is input, each element of 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. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] B |
|
*> \verbatim |
|
*> B is DOUBLE PRECISION array, dimension (LDB,NRHS) |
|
*> The right-hand-side matrix B. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] LDB |
|
*> \verbatim |
|
*> LDB is INTEGER |
|
*> The leading dimension of the array B. LDB >= max(1,N). |
|
*> \endverbatim |
|
*> |
|
*> \param[in,out] Y |
|
*> \verbatim |
|
*> Y is DOUBLE PRECISION array, dimension (LDY,NRHS) |
|
*> On entry, the solution matrix X, as computed by DGBTRS. |
|
*> On exit, the improved solution matrix Y. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] LDY |
|
*> \verbatim |
|
*> LDY is INTEGER |
|
*> The leading dimension of the array Y. LDY >= max(1,N). |
|
*> \endverbatim |
|
*> |
|
*> \param[out] BERR_OUT |
|
*> \verbatim |
|
*> BERR_OUT is DOUBLE PRECISION array, dimension (NRHS) |
|
*> On exit, BERR_OUT(j) contains the componentwise relative backward |
|
*> error for right-hand-side j from the formula |
|
*> max(i) ( abs(RES(i)) / ( abs(op(A_s))*abs(Y) + abs(B_s) )(i) ) |
|
*> where abs(Z) is the componentwise absolute value of the matrix |
|
*> or vector Z. This is computed by DLA_LIN_BERR. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] N_NORMS |
|
*> \verbatim |
|
*> N_NORMS is INTEGER |
|
*> Determines which error bounds to return (see ERR_BNDS_NORM |
|
*> and ERR_BNDS_COMP). |
|
*> If N_NORMS >= 1 return normwise error bounds. |
|
*> If N_NORMS >= 2 return componentwise error bounds. |
|
*> \endverbatim |
|
*> |
|
*> \param[in,out] ERR_BNDS_NORM |
|
*> \verbatim |
|
*> ERR_BNDS_NORM is DOUBLE PRECISION array, dimension (NRHS, N_ERR_BNDS) |
|
*> For each right-hand side, this array 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 ith solution vector: |
|
*> max_j (abs(XTRUE(j,i) - X(j,i))) |
|
*> ------------------------------ |
|
*> max_j abs(X(j,i)) |
|
*> |
|
*> The array is indexed by the type of error information as described |
|
*> below. There currently are up to three pieces of information |
|
*> returned. |
|
*> |
|
*> The first index in ERR_BNDS_NORM(i,:) corresponds to the ith |
|
*> right-hand side. |
|
*> |
|
*> The second index in ERR_BNDS_NORM(:,err) contains the following |
|
*> three fields: |
|
*> err = 1 "Trust/don't trust" boolean. Trust the answer if the |
|
*> reciprocal condition number is less than the threshold |
|
*> sqrt(n) * slamch('Epsilon'). |
|
*> |
|
*> 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('Epsilon'). 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('Epsilon') to determine if the error |
|
*> estimate is "guaranteed". These reciprocal condition |
|
*> numbers are 1 / (norm(Z^{-1},inf) * norm(Z,inf)) for some |
|
*> appropriately scaled matrix Z. |
|
*> 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. |
|
*> |
|
*> This subroutine is only responsible for setting the second field |
|
*> above. |
|
*> See Lapack Working Note 165 for further details and extra |
|
*> cautions. |
|
*> \endverbatim |
|
*> |
|
*> \param[in,out] ERR_BNDS_COMP |
|
*> \verbatim |
|
*> ERR_BNDS_COMP is DOUBLE PRECISION array, dimension (NRHS, N_ERR_BNDS) |
|
*> For each right-hand side, this array 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 ith solution vector: |
|
*> abs(XTRUE(j,i) - X(j,i)) |
|
*> max_j ---------------------- |
|
*> abs(X(j,i)) |
|
*> |
|
*> The array is indexed by the right-hand side i (on which the |
|
*> componentwise relative error depends), and the type of error |
|
*> information as described below. There currently are up to three |
|
*> pieces of information returned for each right-hand side. If |
|
*> componentwise accuracy is not requested (PARAMS(3) = 0.0), then |
|
*> ERR_BNDS_COMP is not accessed. If N_ERR_BNDS .LT. 3, then at most |
|
*> the first (:,N_ERR_BNDS) entries are returned. |
|
*> |
|
*> The first index in ERR_BNDS_COMP(i,:) corresponds to the ith |
|
*> right-hand side. |
|
*> |
|
*> The second index in ERR_BNDS_COMP(:,err) contains the following |
|
*> three fields: |
|
*> err = 1 "Trust/don't trust" boolean. Trust the answer if the |
|
*> reciprocal condition number is less than the threshold |
|
*> sqrt(n) * slamch('Epsilon'). |
|
*> |
|
*> 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('Epsilon'). This error bound should only |
|
*> be trusted if the previous boolean is true. |
|
*> |
|
*> err = 3 Reciprocal condition number: Estimated componentwise |
|
*> reciprocal condition number. Compared with the threshold |
|
*> sqrt(n) * slamch('Epsilon') to determine if the error |
|
*> estimate is "guaranteed". These reciprocal condition |
|
*> numbers are 1 / (norm(Z^{-1},inf) * norm(Z,inf)) for some |
|
*> appropriately scaled matrix Z. |
|
*> Let Z = S*(A*diag(x)), where x is the solution for the |
|
*> current right-hand side and S scales each row of |
|
*> A*diag(x) by a power of the radix so all absolute row |
|
*> sums of Z are approximately 1. |
|
*> |
|
*> This subroutine is only responsible for setting the second field |
|
*> above. |
|
*> See Lapack Working Note 165 for further details and extra |
|
*> cautions. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] RES |
|
*> \verbatim |
|
*> RES is DOUBLE PRECISION array, dimension (N) |
|
*> Workspace to hold the intermediate residual. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] AYB |
|
*> \verbatim |
|
*> AYB is DOUBLE PRECISION array, dimension (N) |
|
*> Workspace. This can be the same workspace passed for Y_TAIL. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] DY |
|
*> \verbatim |
|
*> DY is DOUBLE PRECISION array, dimension (N) |
|
*> Workspace to hold the intermediate solution. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] Y_TAIL |
|
*> \verbatim |
|
*> Y_TAIL is DOUBLE PRECISION array, dimension (N) |
|
*> Workspace to hold the trailing bits of the intermediate solution. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] RCOND |
|
*> \verbatim |
|
*> RCOND is DOUBLE PRECISION |
|
*> Reciprocal scaled condition number. This is an estimate of the |
|
*> reciprocal Skeel condition number of the matrix A after |
|
*> equilibration (if done). If this is less than the machine |
|
*> precision (in particular, if it is zero), 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. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] ITHRESH |
|
*> \verbatim |
|
*> ITHRESH is INTEGER |
|
*> The maximum number of residual computations allowed for |
|
*> refinement. The default is 10. For 'aggressive' set to 100 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. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] RTHRESH |
|
*> \verbatim |
|
*> RTHRESH is DOUBLE PRECISION |
|
*> Determines when to stop refinement if the error estimate stops |
|
*> decreasing. Refinement will stop when the next solution no longer |
|
*> satisfies norm(dx_{i+1}) < RTHRESH * norm(dx_i) where norm(Z) is |
|
*> the infinity norm of Z. RTHRESH satisfies 0 < RTHRESH <= 1. The |
|
*> default value is 0.5. For 'aggressive' set to 0.9 to permit |
|
*> convergence on extremely ill-conditioned matrices. See LAWN 165 |
|
*> for more details. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] DZ_UB |
|
*> \verbatim |
|
*> DZ_UB is DOUBLE PRECISION |
|
*> Determines when to start considering componentwise convergence. |
|
*> Componentwise convergence is only considered after each component |
|
*> of the solution Y is stable, which we definte as the relative |
|
*> change in each component being less than DZ_UB. The default value |
|
*> is 0.25, requiring the first bit to be stable. See LAWN 165 for |
|
*> more details. |
|
*> \endverbatim |
|
*> |
|
*> \param[in] IGNORE_CWISE |
|
*> \verbatim |
|
*> IGNORE_CWISE is LOGICAL |
|
*> If .TRUE. then ignore componentwise convergence. Default value |
|
*> is .FALSE.. |
|
*> \endverbatim |
|
*> |
|
*> \param[out] INFO |
|
*> \verbatim |
|
*> INFO is INTEGER |
|
*> = 0: Successful exit. |
|
*> < 0: if INFO = -i, the ith argument to DGBTRS had an illegal |
|
*> value |
|
*> \endverbatim |
|
* |
|
* Authors: |
|
* ======== |
|
* |
|
*> \author Univ. of Tennessee |
|
*> \author Univ. of California Berkeley |
|
*> \author Univ. of Colorado Denver |
|
*> \author NAG Ltd. |
|
* |
|
*> \date June 2017 |
|
* |
|
*> \ingroup doubleGBcomputational |
|
* |
|
* ===================================================================== |
SUBROUTINE DLA_GBRFSX_EXTENDED( PREC_TYPE, TRANS_TYPE, N, KL, KU, |
SUBROUTINE DLA_GBRFSX_EXTENDED( PREC_TYPE, TRANS_TYPE, N, KL, KU, |
$ NRHS, AB, LDAB, AFB, LDAFB, IPIV, |
$ NRHS, AB, LDAB, AFB, LDAFB, IPIV, |
$ COLEQU, C, B, LDB, Y, LDY, |
$ COLEQU, C, B, LDB, Y, LDY, |
Line 6
|
Line 411
|
$ Y_TAIL, RCOND, ITHRESH, RTHRESH, |
$ Y_TAIL, RCOND, ITHRESH, RTHRESH, |
$ DZ_UB, IGNORE_CWISE, INFO ) |
$ DZ_UB, IGNORE_CWISE, INFO ) |
* |
* |
* -- LAPACK routine (version 3.2.1) -- |
* -- LAPACK computational routine (version 3.7.1) -- |
* -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and -- |
* -- LAPACK is a software package provided by Univ. of Tennessee, -- |
* -- Jason Riedy of Univ. of California Berkeley. -- |
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- |
* -- April 2009 -- |
* June 2017 |
* |
|
* -- LAPACK is a software package provided by Univ. of Tennessee, -- |
|
* -- Univ. of California Berkeley and NAG Ltd. -- |
|
* |
* |
IMPLICIT NONE |
|
* .. |
|
* .. Scalar Arguments .. |
* .. Scalar Arguments .. |
INTEGER INFO, LDAB, LDAFB, LDB, LDY, N, KL, KU, NRHS, |
INTEGER INFO, LDAB, LDAFB, LDB, LDY, N, KL, KU, NRHS, |
$ PREC_TYPE, TRANS_TYPE, N_NORMS, ITHRESH |
$ PREC_TYPE, TRANS_TYPE, N_NORMS, ITHRESH |
Line 31
|
Line 431
|
$ ERR_BNDS_COMP( NRHS, * ) |
$ ERR_BNDS_COMP( NRHS, * ) |
* .. |
* .. |
* |
* |
* Purpose |
|
* ======= |
|
* |
|
* DLA_GBRFSX_EXTENDED improves the computed solution to a system of |
|
* linear equations by performing extra-precise iterative refinement |
|
* and provides error bounds and backward error estimates for the solution. |
|
* This subroutine is called by DGBRFSX to perform iterative refinement. |
|
* In addition to normwise error bound, the code provides maximum |
|
* componentwise error bound if possible. See comments for ERR_BNDS_NORM |
|
* and ERR_BNDS_COMP for details of the error bounds. Note that this |
|
* subroutine is only resonsible for setting the second fields of |
|
* ERR_BNDS_NORM and ERR_BNDS_COMP. |
|
* |
|
* Arguments |
|
* ========= |
|
* |
|
* PREC_TYPE (input) INTEGER |
|
* Specifies the intermediate precision to be used in refinement. |
|
* The value is defined by ILAPREC(P) where P is a CHARACTER and |
|
* P = 'S': Single |
|
* = 'D': Double |
|
* = 'I': Indigenous |
|
* = 'X', 'E': Extra |
|
* |
|
* TRANS_TYPE (input) INTEGER |
|
* Specifies the transposition operation on A. |
|
* The value is defined by ILATRANS(T) where T is a CHARACTER and |
|
* T = 'N': No transpose |
|
* = 'T': Transpose |
|
* = 'C': Conjugate transpose |
|
* |
|
* N (input) INTEGER |
|
* The number of linear equations, i.e., the order of the |
|
* matrix A. N >= 0. |
|
* |
|
* KL (input) INTEGER |
|
* The number of subdiagonals within the band of A. KL >= 0. |
|
* |
|
* KU (input) INTEGER |
|
* The number of superdiagonals within the band of A. KU >= 0 |
|
* |
|
* NRHS (input) INTEGER |
|
* The number of right-hand-sides, i.e., the number of columns of the |
|
* matrix B. |
|
* |
|
* A (input) DOUBLE PRECISION array, dimension (LDA,N) |
|
* On entry, the N-by-N matrix A. |
|
* |
|
* LDA (input) INTEGER |
|
* The leading dimension of the array A. LDA >= max(1,N). |
|
* |
|
* AF (input) DOUBLE PRECISION array, dimension (LDAF,N) |
|
* The factors L and U from the factorization |
|
* A = P*L*U as computed by DGBTRF. |
|
* |
|
* LDAF (input) INTEGER |
|
* The leading dimension of the array AF. LDAF >= max(1,N). |
|
* |
|
* IPIV (input) INTEGER array, dimension (N) |
|
* The pivot indices from the factorization A = P*L*U |
|
* as computed by DGBTRF; row i of the matrix was interchanged |
|
* with row IPIV(i). |
|
* |
|
* COLEQU (input) LOGICAL |
|
* If .TRUE. then column equilibration was done to A before calling |
|
* this routine. This is needed to compute the solution and error |
|
* bounds correctly. |
|
* |
|
* C (input) DOUBLE PRECISION array, dimension (N) |
|
* The column scale factors for A. If COLEQU = .FALSE., C |
|
* is not accessed. If C is input, each element of 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. |
|
* |
|
* B (input) DOUBLE PRECISION array, dimension (LDB,NRHS) |
|
* The right-hand-side matrix B. |
|
* |
|
* LDB (input) INTEGER |
|
* The leading dimension of the array B. LDB >= max(1,N). |
|
* |
|
* Y (input/output) DOUBLE PRECISION array, dimension |
|
* (LDY,NRHS) |
|
* On entry, the solution matrix X, as computed by DGBTRS. |
|
* On exit, the improved solution matrix Y. |
|
* |
|
* LDY (input) INTEGER |
|
* The leading dimension of the array Y. LDY >= max(1,N). |
|
* |
|
* BERR_OUT (output) DOUBLE PRECISION array, dimension (NRHS) |
|
* On exit, BERR_OUT(j) contains the componentwise relative backward |
|
* error for right-hand-side j from the formula |
|
* max(i) ( abs(RES(i)) / ( abs(op(A_s))*abs(Y) + abs(B_s) )(i) ) |
|
* where abs(Z) is the componentwise absolute value of the matrix |
|
* or vector Z. This is computed by DLA_LIN_BERR. |
|
* |
|
* N_NORMS (input) INTEGER |
|
* Determines which error bounds to return (see ERR_BNDS_NORM |
|
* and ERR_BNDS_COMP). |
|
* If N_NORMS >= 1 return normwise error bounds. |
|
* If N_NORMS >= 2 return componentwise error bounds. |
|
* |
|
* ERR_BNDS_NORM (input/output) DOUBLE PRECISION array, dimension |
|
* (NRHS, N_ERR_BNDS) |
|
* For each right-hand side, this array 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 ith solution vector: |
|
* max_j (abs(XTRUE(j,i) - X(j,i))) |
|
* ------------------------------ |
|
* max_j abs(X(j,i)) |
|
* |
|
* The array is indexed by the type of error information as described |
|
* below. There currently are up to three pieces of information |
|
* returned. |
|
* |
|
* The first index in ERR_BNDS_NORM(i,:) corresponds to the ith |
|
* right-hand side. |
|
* |
|
* The second index in ERR_BNDS_NORM(:,err) contains the following |
|
* three fields: |
|
* err = 1 "Trust/don't trust" boolean. Trust the answer if the |
|
* reciprocal condition number is less than the threshold |
|
* sqrt(n) * slamch('Epsilon'). |
|
* |
|
* 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('Epsilon'). 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('Epsilon') to determine if the error |
|
* estimate is "guaranteed". These reciprocal condition |
|
* numbers are 1 / (norm(Z^{-1},inf) * norm(Z,inf)) for some |
|
* appropriately scaled matrix Z. |
|
* 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. |
|
* |
|
* This subroutine is only responsible for setting the second field |
|
* above. |
|
* See Lapack Working Note 165 for further details and extra |
|
* cautions. |
|
* |
|
* ERR_BNDS_COMP (input/output) DOUBLE PRECISION array, dimension |
|
* (NRHS, N_ERR_BNDS) |
|
* For each right-hand side, this array 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 ith solution vector: |
|
* abs(XTRUE(j,i) - X(j,i)) |
|
* max_j ---------------------- |
|
* abs(X(j,i)) |
|
* |
|
* The array is indexed by the right-hand side i (on which the |
|
* componentwise relative error depends), and the type of error |
|
* information as described below. There currently are up to three |
|
* pieces of information returned for each right-hand side. If |
|
* componentwise accuracy is not requested (PARAMS(3) = 0.0), then |
|
* ERR_BNDS_COMP is not accessed. If N_ERR_BNDS .LT. 3, then at most |
|
* the first (:,N_ERR_BNDS) entries are returned. |
|
* |
|
* The first index in ERR_BNDS_COMP(i,:) corresponds to the ith |
|
* right-hand side. |
|
* |
|
* The second index in ERR_BNDS_COMP(:,err) contains the following |
|
* three fields: |
|
* err = 1 "Trust/don't trust" boolean. Trust the answer if the |
|
* reciprocal condition number is less than the threshold |
|
* sqrt(n) * slamch('Epsilon'). |
|
* |
|
* 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('Epsilon'). This error bound should only |
|
* be trusted if the previous boolean is true. |
|
* |
|
* err = 3 Reciprocal condition number: Estimated componentwise |
|
* reciprocal condition number. Compared with the threshold |
|
* sqrt(n) * slamch('Epsilon') to determine if the error |
|
* estimate is "guaranteed". These reciprocal condition |
|
* numbers are 1 / (norm(Z^{-1},inf) * norm(Z,inf)) for some |
|
* appropriately scaled matrix Z. |
|
* Let Z = S*(A*diag(x)), where x is the solution for the |
|
* current right-hand side and S scales each row of |
|
* A*diag(x) by a power of the radix so all absolute row |
|
* sums of Z are approximately 1. |
|
* |
|
* This subroutine is only responsible for setting the second field |
|
* above. |
|
* See Lapack Working Note 165 for further details and extra |
|
* cautions. |
|
* |
|
* RES (input) DOUBLE PRECISION array, dimension (N) |
|
* Workspace to hold the intermediate residual. |
|
* |
|
* AYB (input) DOUBLE PRECISION array, dimension (N) |
|
* Workspace. This can be the same workspace passed for Y_TAIL. |
|
* |
|
* DY (input) DOUBLE PRECISION array, dimension (N) |
|
* Workspace to hold the intermediate solution. |
|
* |
|
* Y_TAIL (input) DOUBLE PRECISION array, dimension (N) |
|
* Workspace to hold the trailing bits of the intermediate solution. |
|
* |
|
* RCOND (input) DOUBLE PRECISION |
|
* Reciprocal scaled condition number. This is an estimate of the |
|
* reciprocal Skeel condition number of the matrix A after |
|
* equilibration (if done). If this is less than the machine |
|
* precision (in particular, if it is zero), 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. |
|
* |
|
* ITHRESH (input) INTEGER |
|
* The maximum number of residual computations allowed for |
|
* refinement. The default is 10. For 'aggressive' set to 100 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. |
|
* |
|
* RTHRESH (input) DOUBLE PRECISION |
|
* Determines when to stop refinement if the error estimate stops |
|
* decreasing. Refinement will stop when the next solution no longer |
|
* satisfies norm(dx_{i+1}) < RTHRESH * norm(dx_i) where norm(Z) is |
|
* the infinity norm of Z. RTHRESH satisfies 0 < RTHRESH <= 1. The |
|
* default value is 0.5. For 'aggressive' set to 0.9 to permit |
|
* convergence on extremely ill-conditioned matrices. See LAWN 165 |
|
* for more details. |
|
* |
|
* DZ_UB (input) DOUBLE PRECISION |
|
* Determines when to start considering componentwise convergence. |
|
* Componentwise convergence is only considered after each component |
|
* of the solution Y is stable, which we definte as the relative |
|
* change in each component being less than DZ_UB. The default value |
|
* is 0.25, requiring the first bit to be stable. See LAWN 165 for |
|
* more details. |
|
* |
|
* IGNORE_CWISE (input) LOGICAL |
|
* If .TRUE. then ignore componentwise convergence. Default value |
|
* is .FALSE.. |
|
* |
|
* INFO (output) INTEGER |
|
* = 0: Successful exit. |
|
* < 0: if INFO = -i, the ith argument to DGBTRS had an illegal |
|
* value |
|
* |
|
* ===================================================================== |
* ===================================================================== |
* |
* |
* .. Local Scalars .. |
* .. Local Scalars .. |