![]() ![]() | ![]() |
Mise à jour de lapack vers la version 3.3.0.
1: SUBROUTINE ZLACN2( N, V, X, EST, KASE, ISAVE ) 2: * 3: * -- LAPACK auxiliary routine (version 3.2) -- 4: * -- LAPACK is a software package provided by Univ. of Tennessee, -- 5: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 6: * November 2006 7: * 8: * .. Scalar Arguments .. 9: INTEGER KASE, N 10: DOUBLE PRECISION EST 11: * .. 12: * .. Array Arguments .. 13: INTEGER ISAVE( 3 ) 14: COMPLEX*16 V( * ), X( * ) 15: * .. 16: * 17: * Purpose 18: * ======= 19: * 20: * ZLACN2 estimates the 1-norm of a square, complex matrix A. 21: * Reverse communication is used for evaluating matrix-vector products. 22: * 23: * Arguments 24: * ========= 25: * 26: * N (input) INTEGER 27: * The order of the matrix. N >= 1. 28: * 29: * V (workspace) COMPLEX*16 array, dimension (N) 30: * On the final return, V = A*W, where EST = norm(V)/norm(W) 31: * (W is not returned). 32: * 33: * X (input/output) COMPLEX*16 array, dimension (N) 34: * On an intermediate return, X should be overwritten by 35: * A * X, if KASE=1, 36: * A' * X, if KASE=2, 37: * where A' is the conjugate transpose of A, and ZLACN2 must be 38: * re-called with all the other parameters unchanged. 39: * 40: * EST (input/output) DOUBLE PRECISION 41: * On entry with KASE = 1 or 2 and ISAVE(1) = 3, EST should be 42: * unchanged from the previous call to ZLACN2. 43: * On exit, EST is an estimate (a lower bound) for norm(A). 44: * 45: * KASE (input/output) INTEGER 46: * On the initial call to ZLACN2, KASE should be 0. 47: * On an intermediate return, KASE will be 1 or 2, indicating 48: * whether X should be overwritten by A * X or A' * X. 49: * On the final return from ZLACN2, KASE will again be 0. 50: * 51: * ISAVE (input/output) INTEGER array, dimension (3) 52: * ISAVE is used to save variables between calls to ZLACN2 53: * 54: * Further Details 55: * ======= ======= 56: * 57: * Contributed by Nick Higham, University of Manchester. 58: * Originally named CONEST, dated March 16, 1988. 59: * 60: * Reference: N.J. Higham, "FORTRAN codes for estimating the one-norm of 61: * a real or complex matrix, with applications to condition estimation", 62: * ACM Trans. Math. Soft., vol. 14, no. 4, pp. 381-396, December 1988. 63: * 64: * Last modified: April, 1999 65: * 66: * This is a thread safe version of ZLACON, which uses the array ISAVE 67: * in place of a SAVE statement, as follows: 68: * 69: * ZLACON ZLACN2 70: * JUMP ISAVE(1) 71: * J ISAVE(2) 72: * ITER ISAVE(3) 73: * 74: * ===================================================================== 75: * 76: * .. Parameters .. 77: INTEGER ITMAX 78: PARAMETER ( ITMAX = 5 ) 79: DOUBLE PRECISION ONE, TWO 80: PARAMETER ( ONE = 1.0D0, TWO = 2.0D0 ) 81: COMPLEX*16 CZERO, CONE 82: PARAMETER ( CZERO = ( 0.0D0, 0.0D0 ), 83: $ CONE = ( 1.0D0, 0.0D0 ) ) 84: * .. 85: * .. Local Scalars .. 86: INTEGER I, JLAST 87: DOUBLE PRECISION ABSXI, ALTSGN, ESTOLD, SAFMIN, TEMP 88: * .. 89: * .. External Functions .. 90: INTEGER IZMAX1 91: DOUBLE PRECISION DLAMCH, DZSUM1 92: EXTERNAL IZMAX1, DLAMCH, DZSUM1 93: * .. 94: * .. External Subroutines .. 95: EXTERNAL ZCOPY 96: * .. 97: * .. Intrinsic Functions .. 98: INTRINSIC ABS, DBLE, DCMPLX, DIMAG 99: * .. 100: * .. Executable Statements .. 101: * 102: SAFMIN = DLAMCH( 'Safe minimum' ) 103: IF( KASE.EQ.0 ) THEN 104: DO 10 I = 1, N 105: X( I ) = DCMPLX( ONE / DBLE( N ) ) 106: 10 CONTINUE 107: KASE = 1 108: ISAVE( 1 ) = 1 109: RETURN 110: END IF 111: * 112: GO TO ( 20, 40, 70, 90, 120 )ISAVE( 1 ) 113: * 114: * ................ ENTRY (ISAVE( 1 ) = 1) 115: * FIRST ITERATION. X HAS BEEN OVERWRITTEN BY A*X. 116: * 117: 20 CONTINUE 118: IF( N.EQ.1 ) THEN 119: V( 1 ) = X( 1 ) 120: EST = ABS( V( 1 ) ) 121: * ... QUIT 122: GO TO 130 123: END IF 124: EST = DZSUM1( N, X, 1 ) 125: * 126: DO 30 I = 1, N 127: ABSXI = ABS( X( I ) ) 128: IF( ABSXI.GT.SAFMIN ) THEN 129: X( I ) = DCMPLX( DBLE( X( I ) ) / ABSXI, 130: $ DIMAG( X( I ) ) / ABSXI ) 131: ELSE 132: X( I ) = CONE 133: END IF 134: 30 CONTINUE 135: KASE = 2 136: ISAVE( 1 ) = 2 137: RETURN 138: * 139: * ................ ENTRY (ISAVE( 1 ) = 2) 140: * FIRST ITERATION. X HAS BEEN OVERWRITTEN BY CTRANS(A)*X. 141: * 142: 40 CONTINUE 143: ISAVE( 2 ) = IZMAX1( N, X, 1 ) 144: ISAVE( 3 ) = 2 145: * 146: * MAIN LOOP - ITERATIONS 2,3,...,ITMAX. 147: * 148: 50 CONTINUE 149: DO 60 I = 1, N 150: X( I ) = CZERO 151: 60 CONTINUE 152: X( ISAVE( 2 ) ) = CONE 153: KASE = 1 154: ISAVE( 1 ) = 3 155: RETURN 156: * 157: * ................ ENTRY (ISAVE( 1 ) = 3) 158: * X HAS BEEN OVERWRITTEN BY A*X. 159: * 160: 70 CONTINUE 161: CALL ZCOPY( N, X, 1, V, 1 ) 162: ESTOLD = EST 163: EST = DZSUM1( N, V, 1 ) 164: * 165: * TEST FOR CYCLING. 166: IF( EST.LE.ESTOLD ) 167: $ GO TO 100 168: * 169: DO 80 I = 1, N 170: ABSXI = ABS( X( I ) ) 171: IF( ABSXI.GT.SAFMIN ) THEN 172: X( I ) = DCMPLX( DBLE( X( I ) ) / ABSXI, 173: $ DIMAG( X( I ) ) / ABSXI ) 174: ELSE 175: X( I ) = CONE 176: END IF 177: 80 CONTINUE 178: KASE = 2 179: ISAVE( 1 ) = 4 180: RETURN 181: * 182: * ................ ENTRY (ISAVE( 1 ) = 4) 183: * X HAS BEEN OVERWRITTEN BY CTRANS(A)*X. 184: * 185: 90 CONTINUE 186: JLAST = ISAVE( 2 ) 187: ISAVE( 2 ) = IZMAX1( N, X, 1 ) 188: IF( ( ABS( X( JLAST ) ).NE.ABS( X( ISAVE( 2 ) ) ) ) .AND. 189: $ ( ISAVE( 3 ).LT.ITMAX ) ) THEN 190: ISAVE( 3 ) = ISAVE( 3 ) + 1 191: GO TO 50 192: END IF 193: * 194: * ITERATION COMPLETE. FINAL STAGE. 195: * 196: 100 CONTINUE 197: ALTSGN = ONE 198: DO 110 I = 1, N 199: X( I ) = DCMPLX( ALTSGN*( ONE+DBLE( I-1 ) / DBLE( N-1 ) ) ) 200: ALTSGN = -ALTSGN 201: 110 CONTINUE 202: KASE = 1 203: ISAVE( 1 ) = 5 204: RETURN 205: * 206: * ................ ENTRY (ISAVE( 1 ) = 5) 207: * X HAS BEEN OVERWRITTEN BY A*X. 208: * 209: 120 CONTINUE 210: TEMP = TWO*( DZSUM1( N, X, 1 ) / DBLE( 3*N ) ) 211: IF( TEMP.GT.EST ) THEN 212: CALL ZCOPY( N, X, 1, V, 1 ) 213: EST = TEMP 214: END IF 215: * 216: 130 CONTINUE 217: KASE = 0 218: RETURN 219: * 220: * End of ZLACN2 221: * 222: END