File:  [local] / rpl / lapack / lapack / zla_syrcond_x.f
Revision 1.4: download - view: text, annotated - select for diffs - revision graph
Tue Dec 21 13:53:48 2010 UTC (13 years, 4 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_0, rpl-4_0_24, rpl-4_0_22, rpl-4_0_21, rpl-4_0_20, rpl-4_0, HEAD
Mise à jour de lapack vers la version 3.3.0.

    1:       DOUBLE PRECISION FUNCTION ZLA_SYRCOND_X( UPLO, N, A, LDA, AF,
    2:      $                                         LDAF, IPIV, X, INFO,
    3:      $                                         WORK, RWORK )
    4: *
    5: *     -- LAPACK routine (version 3.2.1)                                 --
    6: *     -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and --
    7: *     -- Jason Riedy of Univ. of California Berkeley.                 --
    8: *     -- April 2009                                                   --
    9: *
   10: *     -- LAPACK is a software package provided by Univ. of Tennessee, --
   11: *     -- Univ. of California Berkeley and NAG Ltd.                    --
   12: *
   13:       IMPLICIT NONE
   14: *     ..
   15: *     .. Scalar Arguments ..
   16:       CHARACTER          UPLO
   17:       INTEGER            N, LDA, LDAF, INFO
   18: *     ..
   19: *     .. Array Arguments ..
   20:       INTEGER            IPIV( * )
   21:       COMPLEX*16         A( LDA, * ), AF( LDAF, * ), WORK( * ), X( * )
   22:       DOUBLE PRECISION   RWORK( * )
   23: *     ..
   24: *
   25: *  Purpose
   26: *  =======
   27: *
   28: *     ZLA_SYRCOND_X Computes the infinity norm condition number of
   29: *     op(A) * diag(X) where X is a COMPLEX*16 vector.
   30: *
   31: *  Arguments
   32: *  =========
   33: *
   34: *     UPLO    (input) CHARACTER*1
   35: *       = 'U':  Upper triangle of A is stored;
   36: *       = 'L':  Lower triangle of A is stored.
   37: *
   38: *     N       (input) INTEGER
   39: *     The number of linear equations, i.e., the order of the
   40: *     matrix A.  N >= 0.
   41: *
   42: *     A       (input) COMPLEX*16 array, dimension (LDA,N)
   43: *     On entry, the N-by-N matrix A.
   44: *
   45: *     LDA     (input) INTEGER
   46: *     The leading dimension of the array A.  LDA >= max(1,N).
   47: *
   48: *     AF      (input) COMPLEX*16 array, dimension (LDAF,N)
   49: *     The block diagonal matrix D and the multipliers used to
   50: *     obtain the factor U or L as computed by ZSYTRF.
   51: *
   52: *     LDAF    (input) INTEGER
   53: *     The leading dimension of the array AF.  LDAF >= max(1,N).
   54: *
   55: *     IPIV    (input) INTEGER array, dimension (N)
   56: *     Details of the interchanges and the block structure of D
   57: *     as determined by ZSYTRF.
   58: *
   59: *     X       (input) COMPLEX*16 array, dimension (N)
   60: *     The vector X in the formula op(A) * diag(X).
   61: *
   62: *     INFO    (output) INTEGER
   63: *       = 0:  Successful exit.
   64: *     i > 0:  The ith argument is invalid.
   65: *
   66: *     WORK    (input) COMPLEX*16 array, dimension (2*N).
   67: *     Workspace.
   68: *
   69: *     RWORK   (input) DOUBLE PRECISION array, dimension (N).
   70: *     Workspace.
   71: *
   72: *  =====================================================================
   73: *
   74: *     .. Local Scalars ..
   75:       INTEGER            KASE
   76:       DOUBLE PRECISION   AINVNM, ANORM, TMP
   77:       INTEGER            I, J
   78:       LOGICAL            UP
   79:       COMPLEX*16         ZDUM
   80: *     ..
   81: *     .. Local Arrays ..
   82:       INTEGER            ISAVE( 3 )
   83: *     ..
   84: *     .. External Functions ..
   85:       LOGICAL            LSAME
   86:       EXTERNAL           LSAME
   87: *     ..
   88: *     .. External Subroutines ..
   89:       EXTERNAL           ZLACN2, ZSYTRS, XERBLA
   90: *     ..
   91: *     .. Intrinsic Functions ..
   92:       INTRINSIC          ABS, MAX
   93: *     ..
   94: *     .. Statement Functions ..
   95:       DOUBLE PRECISION   CABS1
   96: *     ..
   97: *     .. Statement Function Definitions ..
   98:       CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
   99: *     ..
  100: *     .. Executable Statements ..
  101: *
  102:       ZLA_SYRCOND_X = 0.0D+0
  103: *
  104:       INFO = 0
  105:       IF( N.LT.0 ) THEN
  106:          INFO = -2
  107:       END IF
  108:       IF( INFO.NE.0 ) THEN
  109:          CALL XERBLA( 'ZLA_SYRCOND_X', -INFO )
  110:          RETURN
  111:       END IF
  112:       UP = .FALSE.
  113:       IF ( LSAME( UPLO, 'U' ) ) UP = .TRUE.
  114: *
  115: *     Compute norm of op(A)*op2(C).
  116: *
  117:       ANORM = 0.0D+0
  118:       IF ( UP ) THEN
  119:          DO I = 1, N
  120:             TMP = 0.0D+0
  121:             DO J = 1, I
  122:                TMP = TMP + CABS1( A( J, I ) * X( J ) )
  123:             END DO
  124:             DO J = I+1, N
  125:                TMP = TMP + CABS1( A( I, J ) * X( J ) )
  126:             END DO
  127:             RWORK( I ) = TMP
  128:             ANORM = MAX( ANORM, TMP )
  129:          END DO
  130:       ELSE
  131:          DO I = 1, N
  132:             TMP = 0.0D+0
  133:             DO J = 1, I
  134:                TMP = TMP + CABS1( A( I, J ) * X( J ) )
  135:             END DO
  136:             DO J = I+1, N
  137:                TMP = TMP + CABS1( A( J, I ) * X( J ) )
  138:             END DO
  139:             RWORK( I ) = TMP
  140:             ANORM = MAX( ANORM, TMP )
  141:          END DO
  142:       END IF
  143: *
  144: *     Quick return if possible.
  145: *
  146:       IF( N.EQ.0 ) THEN
  147:          ZLA_SYRCOND_X = 1.0D+0
  148:          RETURN
  149:       ELSE IF( ANORM .EQ. 0.0D+0 ) THEN
  150:          RETURN
  151:       END IF
  152: *
  153: *     Estimate the norm of inv(op(A)).
  154: *
  155:       AINVNM = 0.0D+0
  156: *
  157:       KASE = 0
  158:    10 CONTINUE
  159:       CALL ZLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
  160:       IF( KASE.NE.0 ) THEN
  161:          IF( KASE.EQ.2 ) THEN
  162: *
  163: *           Multiply by R.
  164: *
  165:             DO I = 1, N
  166:                WORK( I ) = WORK( I ) * RWORK( I )
  167:             END DO
  168: *
  169:             IF ( UP ) THEN
  170:                CALL ZSYTRS( 'U', N, 1, AF, LDAF, IPIV,
  171:      $            WORK, N, INFO )
  172:             ELSE
  173:                CALL ZSYTRS( 'L', N, 1, AF, LDAF, IPIV,
  174:      $            WORK, N, INFO )
  175:             ENDIF
  176: *
  177: *           Multiply by inv(X).
  178: *
  179:             DO I = 1, N
  180:                WORK( I ) = WORK( I ) / X( I )
  181:             END DO
  182:          ELSE
  183: *
  184: *           Multiply by inv(X').
  185: *
  186:             DO I = 1, N
  187:                WORK( I ) = WORK( I ) / X( I )
  188:             END DO
  189: *
  190:             IF ( UP ) THEN
  191:                CALL ZSYTRS( 'U', N, 1, AF, LDAF, IPIV,
  192:      $            WORK, N, INFO )
  193:             ELSE
  194:                CALL ZSYTRS( 'L', N, 1, AF, LDAF, IPIV,
  195:      $            WORK, N, INFO )
  196:             END IF
  197: *
  198: *           Multiply by R.
  199: *
  200:             DO I = 1, N
  201:                WORK( I ) = WORK( I ) * RWORK( I )
  202:             END DO
  203:          END IF
  204:          GO TO 10
  205:       END IF
  206: *
  207: *     Compute the estimate of the reciprocal condition number.
  208: *
  209:       IF( AINVNM .NE. 0.0D+0 )
  210:      $   ZLA_SYRCOND_X = 1.0D+0 / AINVNM
  211: *
  212:       RETURN
  213: *
  214:       END

CVSweb interface <joel.bertrand@systella.fr>