File:  [local] / rpl / lapack / lapack / dtrcon.f
Revision 1.18: download - view: text, annotated - select for diffs - revision graph
Mon Aug 7 08:39:13 2023 UTC (9 months, 1 week ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_35, rpl-4_1_34, HEAD
Première mise à jour de lapack et blas.

    1: *> \brief \b DTRCON
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at
    6: *            http://www.netlib.org/lapack/explore-html/
    7: *
    8: *> \htmlonly
    9: *> Download DTRCON + dependencies
   10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dtrcon.f">
   11: *> [TGZ]</a>
   12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dtrcon.f">
   13: *> [ZIP]</a>
   14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dtrcon.f">
   15: *> [TXT]</a>
   16: *> \endhtmlonly
   17: *
   18: *  Definition:
   19: *  ===========
   20: *
   21: *       SUBROUTINE DTRCON( NORM, UPLO, DIAG, N, A, LDA, RCOND, WORK,
   22: *                          IWORK, INFO )
   23: *
   24: *       .. Scalar Arguments ..
   25: *       CHARACTER          DIAG, NORM, UPLO
   26: *       INTEGER            INFO, LDA, N
   27: *       DOUBLE PRECISION   RCOND
   28: *       ..
   29: *       .. Array Arguments ..
   30: *       INTEGER            IWORK( * )
   31: *       DOUBLE PRECISION   A( LDA, * ), WORK( * )
   32: *       ..
   33: *
   34: *
   35: *> \par Purpose:
   36: *  =============
   37: *>
   38: *> \verbatim
   39: *>
   40: *> DTRCON estimates the reciprocal of the condition number of a
   41: *> triangular matrix A, in either the 1-norm or the infinity-norm.
   42: *>
   43: *> The norm of A is computed and an estimate is obtained for
   44: *> norm(inv(A)), then the reciprocal of the condition number is
   45: *> computed as
   46: *>    RCOND = 1 / ( norm(A) * norm(inv(A)) ).
   47: *> \endverbatim
   48: *
   49: *  Arguments:
   50: *  ==========
   51: *
   52: *> \param[in] NORM
   53: *> \verbatim
   54: *>          NORM is CHARACTER*1
   55: *>          Specifies whether the 1-norm condition number or the
   56: *>          infinity-norm condition number is required:
   57: *>          = '1' or 'O':  1-norm;
   58: *>          = 'I':         Infinity-norm.
   59: *> \endverbatim
   60: *>
   61: *> \param[in] UPLO
   62: *> \verbatim
   63: *>          UPLO is CHARACTER*1
   64: *>          = 'U':  A is upper triangular;
   65: *>          = 'L':  A is lower triangular.
   66: *> \endverbatim
   67: *>
   68: *> \param[in] DIAG
   69: *> \verbatim
   70: *>          DIAG is CHARACTER*1
   71: *>          = 'N':  A is non-unit triangular;
   72: *>          = 'U':  A is unit triangular.
   73: *> \endverbatim
   74: *>
   75: *> \param[in] N
   76: *> \verbatim
   77: *>          N is INTEGER
   78: *>          The order of the matrix A.  N >= 0.
   79: *> \endverbatim
   80: *>
   81: *> \param[in] A
   82: *> \verbatim
   83: *>          A is DOUBLE PRECISION array, dimension (LDA,N)
   84: *>          The triangular matrix A.  If UPLO = 'U', the leading N-by-N
   85: *>          upper triangular part of the array A contains the upper
   86: *>          triangular matrix, and the strictly lower triangular part of
   87: *>          A is not referenced.  If UPLO = 'L', the leading N-by-N lower
   88: *>          triangular part of the array A contains the lower triangular
   89: *>          matrix, and the strictly upper triangular part of A is not
   90: *>          referenced.  If DIAG = 'U', the diagonal elements of A are
   91: *>          also not referenced and are assumed to be 1.
   92: *> \endverbatim
   93: *>
   94: *> \param[in] LDA
   95: *> \verbatim
   96: *>          LDA is INTEGER
   97: *>          The leading dimension of the array A.  LDA >= max(1,N).
   98: *> \endverbatim
   99: *>
  100: *> \param[out] RCOND
  101: *> \verbatim
  102: *>          RCOND is DOUBLE PRECISION
  103: *>          The reciprocal of the condition number of the matrix A,
  104: *>          computed as RCOND = 1/(norm(A) * norm(inv(A))).
  105: *> \endverbatim
  106: *>
  107: *> \param[out] WORK
  108: *> \verbatim
  109: *>          WORK is DOUBLE PRECISION array, dimension (3*N)
  110: *> \endverbatim
  111: *>
  112: *> \param[out] IWORK
  113: *> \verbatim
  114: *>          IWORK is INTEGER array, dimension (N)
  115: *> \endverbatim
  116: *>
  117: *> \param[out] INFO
  118: *> \verbatim
  119: *>          INFO is INTEGER
  120: *>          = 0:  successful exit
  121: *>          < 0:  if INFO = -i, the i-th argument had an illegal value
  122: *> \endverbatim
  123: *
  124: *  Authors:
  125: *  ========
  126: *
  127: *> \author Univ. of Tennessee
  128: *> \author Univ. of California Berkeley
  129: *> \author Univ. of Colorado Denver
  130: *> \author NAG Ltd.
  131: *
  132: *> \ingroup doubleOTHERcomputational
  133: *
  134: *  =====================================================================
  135:       SUBROUTINE DTRCON( NORM, UPLO, DIAG, N, A, LDA, RCOND, WORK,
  136:      $                   IWORK, INFO )
  137: *
  138: *  -- LAPACK computational routine --
  139: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  140: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  141: *
  142: *     .. Scalar Arguments ..
  143:       CHARACTER          DIAG, NORM, UPLO
  144:       INTEGER            INFO, LDA, N
  145:       DOUBLE PRECISION   RCOND
  146: *     ..
  147: *     .. Array Arguments ..
  148:       INTEGER            IWORK( * )
  149:       DOUBLE PRECISION   A( LDA, * ), WORK( * )
  150: *     ..
  151: *
  152: *  =====================================================================
  153: *
  154: *     .. Parameters ..
  155:       DOUBLE PRECISION   ONE, ZERO
  156:       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
  157: *     ..
  158: *     .. Local Scalars ..
  159:       LOGICAL            NOUNIT, ONENRM, UPPER
  160:       CHARACTER          NORMIN
  161:       INTEGER            IX, KASE, KASE1
  162:       DOUBLE PRECISION   AINVNM, ANORM, SCALE, SMLNUM, XNORM
  163: *     ..
  164: *     .. Local Arrays ..
  165:       INTEGER            ISAVE( 3 )
  166: *     ..
  167: *     .. External Functions ..
  168:       LOGICAL            LSAME
  169:       INTEGER            IDAMAX
  170:       DOUBLE PRECISION   DLAMCH, DLANTR
  171:       EXTERNAL           LSAME, IDAMAX, DLAMCH, DLANTR
  172: *     ..
  173: *     .. External Subroutines ..
  174:       EXTERNAL           DLACN2, DLATRS, DRSCL, XERBLA
  175: *     ..
  176: *     .. Intrinsic Functions ..
  177:       INTRINSIC          ABS, DBLE, MAX
  178: *     ..
  179: *     .. Executable Statements ..
  180: *
  181: *     Test the input parameters.
  182: *
  183:       INFO = 0
  184:       UPPER = LSAME( UPLO, 'U' )
  185:       ONENRM = NORM.EQ.'1' .OR. LSAME( NORM, 'O' )
  186:       NOUNIT = LSAME( DIAG, 'N' )
  187: *
  188:       IF( .NOT.ONENRM .AND. .NOT.LSAME( NORM, 'I' ) ) THEN
  189:          INFO = -1
  190:       ELSE IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
  191:          INFO = -2
  192:       ELSE IF( .NOT.NOUNIT .AND. .NOT.LSAME( DIAG, 'U' ) ) THEN
  193:          INFO = -3
  194:       ELSE IF( N.LT.0 ) THEN
  195:          INFO = -4
  196:       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
  197:          INFO = -6
  198:       END IF
  199:       IF( INFO.NE.0 ) THEN
  200:          CALL XERBLA( 'DTRCON', -INFO )
  201:          RETURN
  202:       END IF
  203: *
  204: *     Quick return if possible
  205: *
  206:       IF( N.EQ.0 ) THEN
  207:          RCOND = ONE
  208:          RETURN
  209:       END IF
  210: *
  211:       RCOND = ZERO
  212:       SMLNUM = DLAMCH( 'Safe minimum' )*DBLE( MAX( 1, N ) )
  213: *
  214: *     Compute the norm of the triangular matrix A.
  215: *
  216:       ANORM = DLANTR( NORM, UPLO, DIAG, N, N, A, LDA, WORK )
  217: *
  218: *     Continue only if ANORM > 0.
  219: *
  220:       IF( ANORM.GT.ZERO ) THEN
  221: *
  222: *        Estimate the norm of the inverse of A.
  223: *
  224:          AINVNM = ZERO
  225:          NORMIN = 'N'
  226:          IF( ONENRM ) THEN
  227:             KASE1 = 1
  228:          ELSE
  229:             KASE1 = 2
  230:          END IF
  231:          KASE = 0
  232:    10    CONTINUE
  233:          CALL DLACN2( N, WORK( N+1 ), WORK, IWORK, AINVNM, KASE, ISAVE )
  234:          IF( KASE.NE.0 ) THEN
  235:             IF( KASE.EQ.KASE1 ) THEN
  236: *
  237: *              Multiply by inv(A).
  238: *
  239:                CALL DLATRS( UPLO, 'No transpose', DIAG, NORMIN, N, A,
  240:      $                      LDA, WORK, SCALE, WORK( 2*N+1 ), INFO )
  241:             ELSE
  242: *
  243: *              Multiply by inv(A**T).
  244: *
  245:                CALL DLATRS( UPLO, 'Transpose', DIAG, NORMIN, N, A, LDA,
  246:      $                      WORK, SCALE, WORK( 2*N+1 ), INFO )
  247:             END IF
  248:             NORMIN = 'Y'
  249: *
  250: *           Multiply by 1/SCALE if doing so will not cause overflow.
  251: *
  252:             IF( SCALE.NE.ONE ) THEN
  253:                IX = IDAMAX( N, WORK, 1 )
  254:                XNORM = ABS( WORK( IX ) )
  255:                IF( SCALE.LT.XNORM*SMLNUM .OR. SCALE.EQ.ZERO )
  256:      $            GO TO 20
  257:                CALL DRSCL( N, SCALE, WORK, 1 )
  258:             END IF
  259:             GO TO 10
  260:          END IF
  261: *
  262: *        Compute the estimate of the reciprocal condition number.
  263: *
  264:          IF( AINVNM.NE.ZERO )
  265:      $      RCOND = ( ONE / ANORM ) / AINVNM
  266:       END IF
  267: *
  268:    20 CONTINUE
  269:       RETURN
  270: *
  271: *     End of DTRCON
  272: *
  273:       END

CVSweb interface <joel.bertrand@systella.fr>