Annotation of rpl/lapack/lapack/zla_hercond_c.f, revision 1.1
1.1 ! bertrand 1: DOUBLE PRECISION FUNCTION ZLA_HERCOND_C( UPLO, N, A, LDA, AF,
! 2: $ LDAF, IPIV, C, CAPPLY,
! 3: $ INFO, 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: LOGICAL CAPPLY
! 18: INTEGER N, LDA, LDAF, INFO
! 19: * ..
! 20: * .. Array Arguments ..
! 21: INTEGER IPIV( * )
! 22: COMPLEX*16 A( LDA, * ), AF( LDAF, * ), WORK( * )
! 23: DOUBLE PRECISION C ( * ), RWORK( * )
! 24: * ..
! 25: *
! 26: * Purpose
! 27: * =======
! 28: *
! 29: * ZLA_HERCOND_C computes the infinity norm condition number of
! 30: * op(A) * inv(diag(C)) where C is a DOUBLE PRECISION vector.
! 31: *
! 32: * Arguments
! 33: * =========
! 34: *
! 35: * UPLO (input) CHARACTER*1
! 36: * = 'U': Upper triangle of A is stored;
! 37: * = 'L': Lower triangle of A is stored.
! 38: *
! 39: * N (input) INTEGER
! 40: * The number of linear equations, i.e., the order of the
! 41: * matrix A. N >= 0.
! 42: *
! 43: * A (input) COMPLEX*16 array, dimension (LDA,N)
! 44: * On entry, the N-by-N matrix A
! 45: *
! 46: * LDA (input) INTEGER
! 47: * The leading dimension of the array A. LDA >= max(1,N).
! 48: *
! 49: * AF (input) COMPLEX*16 array, dimension (LDAF,N)
! 50: * The block diagonal matrix D and the multipliers used to
! 51: * obtain the factor U or L as computed by ZHETRF.
! 52: *
! 53: * LDAF (input) INTEGER
! 54: * The leading dimension of the array AF. LDAF >= max(1,N).
! 55: *
! 56: * IPIV (input) INTEGER array, dimension (N)
! 57: * Details of the interchanges and the block structure of D
! 58: * as determined by CHETRF.
! 59: *
! 60: * C (input) DOUBLE PRECISION array, dimension (N)
! 61: * The vector C in the formula op(A) * inv(diag(C)).
! 62: *
! 63: * CAPPLY (input) LOGICAL
! 64: * If .TRUE. then access the vector C in the formula above.
! 65: *
! 66: * INFO (output) INTEGER
! 67: * = 0: Successful exit.
! 68: * i > 0: The ith argument is invalid.
! 69: *
! 70: * WORK (input) COMPLEX*16 array, dimension (2*N).
! 71: * Workspace.
! 72: *
! 73: * RWORK (input) DOUBLE PRECISION array, dimension (N).
! 74: * Workspace.
! 75: *
! 76: * =====================================================================
! 77: *
! 78: * .. Local Scalars ..
! 79: INTEGER KASE, I, J
! 80: DOUBLE PRECISION AINVNM, ANORM, TMP
! 81: LOGICAL UP
! 82: COMPLEX*16 ZDUM
! 83: * ..
! 84: * .. Local Arrays ..
! 85: INTEGER ISAVE( 3 )
! 86: * ..
! 87: * .. External Functions ..
! 88: LOGICAL LSAME
! 89: EXTERNAL LSAME
! 90: * ..
! 91: * .. External Subroutines ..
! 92: EXTERNAL ZLACN2, ZHETRS, XERBLA
! 93: * ..
! 94: * .. Intrinsic Functions ..
! 95: INTRINSIC ABS, MAX
! 96: * ..
! 97: * .. Statement Functions ..
! 98: DOUBLE PRECISION CABS1
! 99: * ..
! 100: * .. Statement Function Definitions ..
! 101: CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
! 102: * ..
! 103: * .. Executable Statements ..
! 104: *
! 105: ZLA_HERCOND_C = 0.0D+0
! 106: *
! 107: INFO = 0
! 108: IF( N.LT.0 ) THEN
! 109: INFO = -2
! 110: END IF
! 111: IF( INFO.NE.0 ) THEN
! 112: CALL XERBLA( 'ZLA_HERCOND_C', -INFO )
! 113: RETURN
! 114: END IF
! 115: UP = .FALSE.
! 116: IF ( LSAME( UPLO, 'U' ) ) UP = .TRUE.
! 117: *
! 118: * Compute norm of op(A)*op2(C).
! 119: *
! 120: ANORM = 0.0D+0
! 121: IF ( UP ) THEN
! 122: DO I = 1, N
! 123: TMP = 0.0D+0
! 124: IF ( CAPPLY ) THEN
! 125: DO J = 1, I
! 126: TMP = TMP + CABS1( A( J, I ) ) / C( J )
! 127: END DO
! 128: DO J = I+1, N
! 129: TMP = TMP + CABS1( A( I, J ) ) / C( J )
! 130: END DO
! 131: ELSE
! 132: DO J = 1, I
! 133: TMP = TMP + CABS1( A( J, I ) )
! 134: END DO
! 135: DO J = I+1, N
! 136: TMP = TMP + CABS1( A( I, J ) )
! 137: END DO
! 138: END IF
! 139: RWORK( I ) = TMP
! 140: ANORM = MAX( ANORM, TMP )
! 141: END DO
! 142: ELSE
! 143: DO I = 1, N
! 144: TMP = 0.0D+0
! 145: IF ( CAPPLY ) THEN
! 146: DO J = 1, I
! 147: TMP = TMP + CABS1( A( I, J ) ) / C( J )
! 148: END DO
! 149: DO J = I+1, N
! 150: TMP = TMP + CABS1( A( J, I ) ) / C( J )
! 151: END DO
! 152: ELSE
! 153: DO J = 1, I
! 154: TMP = TMP + CABS1( A( I, J ) )
! 155: END DO
! 156: DO J = I+1, N
! 157: TMP = TMP + CABS1( A( J, I ) )
! 158: END DO
! 159: END IF
! 160: RWORK( I ) = TMP
! 161: ANORM = MAX( ANORM, TMP )
! 162: END DO
! 163: END IF
! 164: *
! 165: * Quick return if possible.
! 166: *
! 167: IF( N.EQ.0 ) THEN
! 168: ZLA_HERCOND_C = 1.0D+0
! 169: RETURN
! 170: ELSE IF( ANORM .EQ. 0.0D+0 ) THEN
! 171: RETURN
! 172: END IF
! 173: *
! 174: * Estimate the norm of inv(op(A)).
! 175: *
! 176: AINVNM = 0.0D+0
! 177: *
! 178: KASE = 0
! 179: 10 CONTINUE
! 180: CALL ZLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
! 181: IF( KASE.NE.0 ) THEN
! 182: IF( KASE.EQ.2 ) THEN
! 183: *
! 184: * Multiply by R.
! 185: *
! 186: DO I = 1, N
! 187: WORK( I ) = WORK( I ) * RWORK( I )
! 188: END DO
! 189: *
! 190: IF ( UP ) THEN
! 191: CALL ZHETRS( 'U', N, 1, AF, LDAF, IPIV,
! 192: $ WORK, N, INFO )
! 193: ELSE
! 194: CALL ZHETRS( 'L', N, 1, AF, LDAF, IPIV,
! 195: $ WORK, N, INFO )
! 196: ENDIF
! 197: *
! 198: * Multiply by inv(C).
! 199: *
! 200: IF ( CAPPLY ) THEN
! 201: DO I = 1, N
! 202: WORK( I ) = WORK( I ) * C( I )
! 203: END DO
! 204: END IF
! 205: ELSE
! 206: *
! 207: * Multiply by inv(C').
! 208: *
! 209: IF ( CAPPLY ) THEN
! 210: DO I = 1, N
! 211: WORK( I ) = WORK( I ) * C( I )
! 212: END DO
! 213: END IF
! 214: *
! 215: IF ( UP ) THEN
! 216: CALL ZHETRS( 'U', N, 1, AF, LDAF, IPIV,
! 217: $ WORK, N, INFO )
! 218: ELSE
! 219: CALL ZHETRS( 'L', N, 1, AF, LDAF, IPIV,
! 220: $ WORK, N, INFO )
! 221: END IF
! 222: *
! 223: * Multiply by R.
! 224: *
! 225: DO I = 1, N
! 226: WORK( I ) = WORK( I ) * RWORK( I )
! 227: END DO
! 228: END IF
! 229: GO TO 10
! 230: END IF
! 231: *
! 232: * Compute the estimate of the reciprocal condition number.
! 233: *
! 234: IF( AINVNM .NE. 0.0D+0 )
! 235: $ ZLA_HERCOND_C = 1.0D+0 / AINVNM
! 236: *
! 237: RETURN
! 238: *
! 239: END
CVSweb interface <joel.bertrand@systella.fr>