File:  [local] / rpl / lapack / lapack / zla_gercond_c.f
Revision 1.12: download - view: text, annotated - select for diffs - revision graph
Mon Jan 27 09:28:36 2014 UTC (10 years, 3 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_24, rpl-4_1_23, rpl-4_1_22, rpl-4_1_21, rpl-4_1_20, rpl-4_1_19, rpl-4_1_18, rpl-4_1_17, HEAD
Cohérence.

    1: *> \brief \b ZLA_GERCOND_C computes the infinity norm condition number of op(A)*inv(diag(c)) for general matrices.
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at 
    6: *            http://www.netlib.org/lapack/explore-html/ 
    7: *
    8: *> \htmlonly
    9: *> Download ZLA_GERCOND_C + dependencies 
   10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zla_gercond_c.f"> 
   11: *> [TGZ]</a> 
   12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zla_gercond_c.f"> 
   13: *> [ZIP]</a> 
   14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zla_gercond_c.f"> 
   15: *> [TXT]</a>
   16: *> \endhtmlonly 
   17: *
   18: *  Definition:
   19: *  ===========
   20: *
   21: *       DOUBLE PRECISION FUNCTION ZLA_GERCOND_C( TRANS, N, A, LDA, AF, 
   22: *                                                LDAF, IPIV, C, CAPPLY,
   23: *                                                INFO, WORK, RWORK )
   24:    25: *       .. Scalar Aguments ..
   26: *       CHARACTER          TRANS
   27: *       LOGICAL            CAPPLY
   28: *       INTEGER            N, LDA, LDAF, INFO
   29: *       ..
   30: *       .. Array Arguments ..
   31: *       INTEGER            IPIV( * )
   32: *       COMPLEX*16         A( LDA, * ), AF( LDAF, * ), WORK( * )
   33: *       DOUBLE PRECISION   C( * ), RWORK( * )
   34: *       ..
   35: *  
   36: *
   37: *> \par Purpose:
   38: *  =============
   39: *>
   40: *> \verbatim
   41: *>
   42: *>    ZLA_GERCOND_C computes the infinity norm condition number of
   43: *>    op(A) * inv(diag(C)) where C is a DOUBLE PRECISION vector.
   44: *> \endverbatim
   45: *
   46: *  Arguments:
   47: *  ==========
   48: *
   49: *> \param[in] TRANS
   50: *> \verbatim
   51: *>          TRANS is CHARACTER*1
   52: *>     Specifies the form of the system of equations:
   53: *>       = 'N':  A * X = B     (No transpose)
   54: *>       = 'T':  A**T * X = B  (Transpose)
   55: *>       = 'C':  A**H * X = B  (Conjugate Transpose = Transpose)
   56: *> \endverbatim
   57: *>
   58: *> \param[in] N
   59: *> \verbatim
   60: *>          N is INTEGER
   61: *>     The number of linear equations, i.e., the order of the
   62: *>     matrix A.  N >= 0.
   63: *> \endverbatim
   64: *>
   65: *> \param[in] A
   66: *> \verbatim
   67: *>          A is COMPLEX*16 array, dimension (LDA,N)
   68: *>     On entry, the N-by-N matrix A
   69: *> \endverbatim
   70: *>
   71: *> \param[in] LDA
   72: *> \verbatim
   73: *>          LDA is INTEGER
   74: *>     The leading dimension of the array A.  LDA >= max(1,N).
   75: *> \endverbatim
   76: *>
   77: *> \param[in] AF
   78: *> \verbatim
   79: *>          AF is COMPLEX*16 array, dimension (LDAF,N)
   80: *>     The factors L and U from the factorization
   81: *>     A = P*L*U as computed by ZGETRF.
   82: *> \endverbatim
   83: *>
   84: *> \param[in] LDAF
   85: *> \verbatim
   86: *>          LDAF is INTEGER
   87: *>     The leading dimension of the array AF.  LDAF >= max(1,N).
   88: *> \endverbatim
   89: *>
   90: *> \param[in] IPIV
   91: *> \verbatim
   92: *>          IPIV is INTEGER array, dimension (N)
   93: *>     The pivot indices from the factorization A = P*L*U
   94: *>     as computed by ZGETRF; row i of the matrix was interchanged
   95: *>     with row IPIV(i).
   96: *> \endverbatim
   97: *>
   98: *> \param[in] C
   99: *> \verbatim
  100: *>          C is DOUBLE PRECISION array, dimension (N)
  101: *>     The vector C in the formula op(A) * inv(diag(C)).
  102: *> \endverbatim
  103: *>
  104: *> \param[in] CAPPLY
  105: *> \verbatim
  106: *>          CAPPLY is LOGICAL
  107: *>     If .TRUE. then access the vector C in the formula above.
  108: *> \endverbatim
  109: *>
  110: *> \param[out] INFO
  111: *> \verbatim
  112: *>          INFO is INTEGER
  113: *>       = 0:  Successful exit.
  114: *>     i > 0:  The ith argument is invalid.
  115: *> \endverbatim
  116: *>
  117: *> \param[in] WORK
  118: *> \verbatim
  119: *>          WORK is COMPLEX*16 array, dimension (2*N).
  120: *>     Workspace.
  121: *> \endverbatim
  122: *>
  123: *> \param[in] RWORK
  124: *> \verbatim
  125: *>          RWORK is DOUBLE PRECISION array, dimension (N).
  126: *>     Workspace.
  127: *> \endverbatim
  128: *
  129: *  Authors:
  130: *  ========
  131: *
  132: *> \author Univ. of Tennessee 
  133: *> \author Univ. of California Berkeley 
  134: *> \author Univ. of Colorado Denver 
  135: *> \author NAG Ltd. 
  136: *
  137: *> \date September 2012
  138: *
  139: *> \ingroup complex16GEcomputational
  140: *
  141: *  =====================================================================
  142:       DOUBLE PRECISION FUNCTION ZLA_GERCOND_C( TRANS, N, A, LDA, AF, 
  143:      $                                         LDAF, IPIV, C, CAPPLY,
  144:      $                                         INFO, WORK, RWORK )
  145: *
  146: *  -- LAPACK computational routine (version 3.4.2) --
  147: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  148: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  149: *     September 2012
  150: *
  151: *     .. Scalar Aguments ..
  152:       CHARACTER          TRANS
  153:       LOGICAL            CAPPLY
  154:       INTEGER            N, LDA, LDAF, INFO
  155: *     ..
  156: *     .. Array Arguments ..
  157:       INTEGER            IPIV( * )
  158:       COMPLEX*16         A( LDA, * ), AF( LDAF, * ), WORK( * )
  159:       DOUBLE PRECISION   C( * ), RWORK( * )
  160: *     ..
  161: *
  162: *  =====================================================================
  163: *
  164: *     .. Local Scalars ..
  165:       LOGICAL            NOTRANS
  166:       INTEGER            KASE, I, J
  167:       DOUBLE PRECISION   AINVNM, ANORM, TMP
  168:       COMPLEX*16         ZDUM
  169: *     ..
  170: *     .. Local Arrays ..
  171:       INTEGER            ISAVE( 3 )
  172: *     ..
  173: *     .. External Functions ..
  174:       LOGICAL            LSAME
  175:       EXTERNAL           LSAME
  176: *     ..
  177: *     .. External Subroutines ..
  178:       EXTERNAL           ZLACN2, ZGETRS, XERBLA
  179: *     ..
  180: *     .. Intrinsic Functions ..
  181:       INTRINSIC          ABS, MAX, REAL, DIMAG
  182: *     ..
  183: *     .. Statement Functions ..
  184:       DOUBLE PRECISION   CABS1
  185: *     ..
  186: *     .. Statement Function Definitions ..
  187:       CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
  188: *     ..
  189: *     .. Executable Statements ..
  190:       ZLA_GERCOND_C = 0.0D+0
  191: *
  192:       INFO = 0
  193:       NOTRANS = LSAME( TRANS, 'N' )
  194:       IF ( .NOT. NOTRANS .AND. .NOT. LSAME( TRANS, 'T' ) .AND. .NOT.
  195:      $     LSAME( TRANS, 'C' ) ) THEN
  196:          INFO = -1
  197:       ELSE IF( N.LT.0 ) THEN
  198:          INFO = -2
  199:       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
  200:          INFO = -4
  201:       ELSE IF( LDAF.LT.MAX( 1, N ) ) THEN
  202:          INFO = -6
  203:       END IF
  204:       IF( INFO.NE.0 ) THEN
  205:          CALL XERBLA( 'ZLA_GERCOND_C', -INFO )
  206:          RETURN
  207:       END IF
  208: *
  209: *     Compute norm of op(A)*op2(C).
  210: *
  211:       ANORM = 0.0D+0
  212:       IF ( NOTRANS ) THEN
  213:          DO I = 1, N
  214:             TMP = 0.0D+0
  215:             IF ( CAPPLY ) THEN
  216:                DO J = 1, N
  217:                   TMP = TMP + CABS1( A( I, J ) ) / C( J )
  218:                END DO
  219:             ELSE
  220:                DO J = 1, N
  221:                   TMP = TMP + CABS1( A( I, J ) )
  222:                END DO
  223:             END IF
  224:             RWORK( I ) = TMP
  225:             ANORM = MAX( ANORM, TMP )
  226:          END DO
  227:       ELSE
  228:          DO I = 1, N
  229:             TMP = 0.0D+0
  230:             IF ( CAPPLY ) THEN
  231:                DO J = 1, N
  232:                   TMP = TMP + CABS1( A( J, I ) ) / C( J )
  233:                END DO
  234:             ELSE
  235:                DO J = 1, N
  236:                   TMP = TMP + CABS1( A( J, I ) )
  237:                END DO
  238:             END IF
  239:             RWORK( I ) = TMP
  240:             ANORM = MAX( ANORM, TMP )
  241:          END DO
  242:       END IF
  243: *
  244: *     Quick return if possible.
  245: *
  246:       IF( N.EQ.0 ) THEN
  247:          ZLA_GERCOND_C = 1.0D+0
  248:          RETURN
  249:       ELSE IF( ANORM .EQ. 0.0D+0 ) THEN
  250:          RETURN
  251:       END IF
  252: *
  253: *     Estimate the norm of inv(op(A)).
  254: *
  255:       AINVNM = 0.0D+0
  256: *
  257:       KASE = 0
  258:    10 CONTINUE
  259:       CALL ZLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
  260:       IF( KASE.NE.0 ) THEN
  261:          IF( KASE.EQ.2 ) THEN
  262: *
  263: *           Multiply by R.
  264: *
  265:             DO I = 1, N
  266:                WORK( I ) = WORK( I ) * RWORK( I )
  267:             END DO
  268: *
  269:             IF (NOTRANS) THEN
  270:                CALL ZGETRS( 'No transpose', N, 1, AF, LDAF, IPIV,
  271:      $            WORK, N, INFO )
  272:             ELSE
  273:                CALL ZGETRS( 'Conjugate transpose', N, 1, AF, LDAF, IPIV,
  274:      $            WORK, N, INFO )
  275:             ENDIF
  276: *
  277: *           Multiply by inv(C).
  278: *
  279:             IF ( CAPPLY ) THEN
  280:                DO I = 1, N
  281:                   WORK( I ) = WORK( I ) * C( I )
  282:                END DO
  283:             END IF
  284:          ELSE
  285: *
  286: *           Multiply by inv(C**H).
  287: *
  288:             IF ( CAPPLY ) THEN
  289:                DO I = 1, N
  290:                   WORK( I ) = WORK( I ) * C( I )
  291:                END DO
  292:             END IF
  293: *
  294:             IF ( NOTRANS ) THEN
  295:                CALL ZGETRS( 'Conjugate transpose', N, 1, AF, LDAF, IPIV,
  296:      $            WORK, N, INFO )
  297:             ELSE
  298:                CALL ZGETRS( 'No transpose', N, 1, AF, LDAF, IPIV,
  299:      $            WORK, N, INFO )
  300:             END IF
  301: *
  302: *           Multiply by R.
  303: *
  304:             DO I = 1, N
  305:                WORK( I ) = WORK( I ) * RWORK( I )
  306:             END DO
  307:          END IF
  308:          GO TO 10
  309:       END IF
  310: *
  311: *     Compute the estimate of the reciprocal condition number.
  312: *
  313:       IF( AINVNM .NE. 0.0D+0 )
  314:      $   ZLA_GERCOND_C = 1.0D+0 / AINVNM
  315: *
  316:       RETURN
  317: *
  318:       END

CVSweb interface <joel.bertrand@systella.fr>