![]() ![]() | ![]() |
Mise à jour de lapack vers la version 3.3.0.
1: SUBROUTINE ZGTCON( NORM, N, DL, D, DU, DU2, IPIV, ANORM, RCOND, 2: $ WORK, INFO ) 3: * 4: * -- LAPACK routine (version 3.2) -- 5: * -- LAPACK is a software package provided by Univ. of Tennessee, -- 6: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 7: * November 2006 8: * 9: * Modified to call ZLACN2 in place of ZLACON, 10 Feb 03, SJH. 10: * 11: * .. Scalar Arguments .. 12: CHARACTER NORM 13: INTEGER INFO, N 14: DOUBLE PRECISION ANORM, RCOND 15: * .. 16: * .. Array Arguments .. 17: INTEGER IPIV( * ) 18: COMPLEX*16 D( * ), DL( * ), DU( * ), DU2( * ), WORK( * ) 19: * .. 20: * 21: * Purpose 22: * ======= 23: * 24: * ZGTCON estimates the reciprocal of the condition number of a complex 25: * tridiagonal matrix A using the LU factorization as computed by 26: * ZGTTRF. 27: * 28: * An estimate is obtained for norm(inv(A)), and the reciprocal of the 29: * condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))). 30: * 31: * Arguments 32: * ========= 33: * 34: * NORM (input) CHARACTER*1 35: * Specifies whether the 1-norm condition number or the 36: * infinity-norm condition number is required: 37: * = '1' or 'O': 1-norm; 38: * = 'I': Infinity-norm. 39: * 40: * N (input) INTEGER 41: * The order of the matrix A. N >= 0. 42: * 43: * DL (input) COMPLEX*16 array, dimension (N-1) 44: * The (n-1) multipliers that define the matrix L from the 45: * LU factorization of A as computed by ZGTTRF. 46: * 47: * D (input) COMPLEX*16 array, dimension (N) 48: * The n diagonal elements of the upper triangular matrix U from 49: * the LU factorization of A. 50: * 51: * DU (input) COMPLEX*16 array, dimension (N-1) 52: * The (n-1) elements of the first superdiagonal of U. 53: * 54: * DU2 (input) COMPLEX*16 array, dimension (N-2) 55: * The (n-2) elements of the second superdiagonal of U. 56: * 57: * IPIV (input) INTEGER array, dimension (N) 58: * The pivot indices; for 1 <= i <= n, row i of the matrix was 59: * interchanged with row IPIV(i). IPIV(i) will always be either 60: * i or i+1; IPIV(i) = i indicates a row interchange was not 61: * required. 62: * 63: * ANORM (input) DOUBLE PRECISION 64: * If NORM = '1' or 'O', the 1-norm of the original matrix A. 65: * If NORM = 'I', the infinity-norm of the original matrix A. 66: * 67: * RCOND (output) DOUBLE PRECISION 68: * The reciprocal of the condition number of the matrix A, 69: * computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is an 70: * estimate of the 1-norm of inv(A) computed in this routine. 71: * 72: * WORK (workspace) COMPLEX*16 array, dimension (2*N) 73: * 74: * INFO (output) INTEGER 75: * = 0: successful exit 76: * < 0: if INFO = -i, the i-th argument had an illegal value 77: * 78: * ===================================================================== 79: * 80: * .. Parameters .. 81: DOUBLE PRECISION ONE, ZERO 82: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 83: * .. 84: * .. Local Scalars .. 85: LOGICAL ONENRM 86: INTEGER I, KASE, KASE1 87: DOUBLE PRECISION AINVNM 88: * .. 89: * .. Local Arrays .. 90: INTEGER ISAVE( 3 ) 91: * .. 92: * .. External Functions .. 93: LOGICAL LSAME 94: EXTERNAL LSAME 95: * .. 96: * .. External Subroutines .. 97: EXTERNAL XERBLA, ZGTTRS, ZLACN2 98: * .. 99: * .. Intrinsic Functions .. 100: INTRINSIC DCMPLX 101: * .. 102: * .. Executable Statements .. 103: * 104: * Test the input arguments. 105: * 106: INFO = 0 107: ONENRM = NORM.EQ.'1' .OR. LSAME( NORM, 'O' ) 108: IF( .NOT.ONENRM .AND. .NOT.LSAME( NORM, 'I' ) ) THEN 109: INFO = -1 110: ELSE IF( N.LT.0 ) THEN 111: INFO = -2 112: ELSE IF( ANORM.LT.ZERO ) THEN 113: INFO = -8 114: END IF 115: IF( INFO.NE.0 ) THEN 116: CALL XERBLA( 'ZGTCON', -INFO ) 117: RETURN 118: END IF 119: * 120: * Quick return if possible 121: * 122: RCOND = ZERO 123: IF( N.EQ.0 ) THEN 124: RCOND = ONE 125: RETURN 126: ELSE IF( ANORM.EQ.ZERO ) THEN 127: RETURN 128: END IF 129: * 130: * Check that D(1:N) is non-zero. 131: * 132: DO 10 I = 1, N 133: IF( D( I ).EQ.DCMPLX( ZERO ) ) 134: $ RETURN 135: 10 CONTINUE 136: * 137: AINVNM = ZERO 138: IF( ONENRM ) THEN 139: KASE1 = 1 140: ELSE 141: KASE1 = 2 142: END IF 143: KASE = 0 144: 20 CONTINUE 145: CALL ZLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE ) 146: IF( KASE.NE.0 ) THEN 147: IF( KASE.EQ.KASE1 ) THEN 148: * 149: * Multiply by inv(U)*inv(L). 150: * 151: CALL ZGTTRS( 'No transpose', N, 1, DL, D, DU, DU2, IPIV, 152: $ WORK, N, INFO ) 153: ELSE 154: * 155: * Multiply by inv(L')*inv(U'). 156: * 157: CALL ZGTTRS( 'Conjugate transpose', N, 1, DL, D, DU, DU2, 158: $ IPIV, WORK, N, INFO ) 159: END IF 160: GO TO 20 161: END IF 162: * 163: * Compute the estimate of the reciprocal condition number. 164: * 165: IF( AINVNM.NE.ZERO ) 166: $ RCOND = ( ONE / AINVNM ) / ANORM 167: * 168: RETURN 169: * 170: * End of ZGTCON 171: * 172: END