 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE ZLACN2( N, V, X, EST, KASE, ISAVE ) 00002 * 00003 * -- LAPACK auxiliary routine (version 3.2) -- 00004 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00005 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00006 * November 2006 00007 * 00008 * .. Scalar Arguments .. 00009 INTEGER KASE, N 00010 DOUBLE PRECISION EST 00011 * .. 00012 * .. Array Arguments .. 00013 INTEGER ISAVE( 3 ) 00014 COMPLEX*16 V( * ), X( * ) 00015 * .. 00016 * 00017 * Purpose 00018 * ======= 00019 * 00020 * ZLACN2 estimates the 1-norm of a square, complex matrix A. 00021 * Reverse communication is used for evaluating matrix-vector products. 00022 * 00023 * Arguments 00024 * ========= 00025 * 00026 * N (input) INTEGER 00027 * The order of the matrix. N >= 1. 00028 * 00029 * V (workspace) COMPLEX*16 array, dimension (N) 00030 * On the final return, V = A*W, where EST = norm(V)/norm(W) 00031 * (W is not returned). 00032 * 00033 * X (input/output) COMPLEX*16 array, dimension (N) 00034 * On an intermediate return, X should be overwritten by 00035 * A * X, if KASE=1, 00036 * A**H * X, if KASE=2, 00037 * where A**H is the conjugate transpose of A, and ZLACN2 must be 00038 * re-called with all the other parameters unchanged. 00039 * 00040 * EST (input/output) DOUBLE PRECISION 00041 * On entry with KASE = 1 or 2 and ISAVE(1) = 3, EST should be 00042 * unchanged from the previous call to ZLACN2. 00043 * On exit, EST is an estimate (a lower bound) for norm(A). 00044 * 00045 * KASE (input/output) INTEGER 00046 * On the initial call to ZLACN2, KASE should be 0. 00047 * On an intermediate return, KASE will be 1 or 2, indicating 00048 * whether X should be overwritten by A * X or A**H * X. 00049 * On the final return from ZLACN2, KASE will again be 0. 00050 * 00051 * ISAVE (input/output) INTEGER array, dimension (3) 00052 * ISAVE is used to save variables between calls to ZLACN2 00053 * 00054 * Further Details 00055 * ======= ======= 00056 * 00057 * Contributed by Nick Higham, University of Manchester. 00058 * Originally named CONEST, dated March 16, 1988. 00059 * 00060 * Reference: N.J. Higham, "FORTRAN codes for estimating the one-norm of 00061 * a real or complex matrix, with applications to condition estimation", 00062 * ACM Trans. Math. Soft., vol. 14, no. 4, pp. 381-396, December 1988. 00063 * 00064 * Last modified: April, 1999 00065 * 00066 * This is a thread safe version of ZLACON, which uses the array ISAVE 00067 * in place of a SAVE statement, as follows: 00068 * 00069 * ZLACON ZLACN2 00070 * JUMP ISAVE(1) 00071 * J ISAVE(2) 00072 * ITER ISAVE(3) 00073 * 00074 * ===================================================================== 00075 * 00076 * .. Parameters .. 00077 INTEGER ITMAX 00078 PARAMETER ( ITMAX = 5 ) 00079 DOUBLE PRECISION ONE, TWO 00080 PARAMETER ( ONE = 1.0D0, TWO = 2.0D0 ) 00081 COMPLEX*16 CZERO, CONE 00082 PARAMETER ( CZERO = ( 0.0D0, 0.0D0 ), 00083 $ CONE = ( 1.0D0, 0.0D0 ) ) 00084 * .. 00085 * .. Local Scalars .. 00086 INTEGER I, JLAST 00087 DOUBLE PRECISION ABSXI, ALTSGN, ESTOLD, SAFMIN, TEMP 00088 * .. 00089 * .. External Functions .. 00090 INTEGER IZMAX1 00091 DOUBLE PRECISION DLAMCH, DZSUM1 00092 EXTERNAL IZMAX1, DLAMCH, DZSUM1 00093 * .. 00094 * .. External Subroutines .. 00095 EXTERNAL ZCOPY 00096 * .. 00097 * .. Intrinsic Functions .. 00098 INTRINSIC ABS, DBLE, DCMPLX, DIMAG 00099 * .. 00100 * .. Executable Statements .. 00101 * 00102 SAFMIN = DLAMCH( 'Safe minimum' ) 00103 IF( KASE.EQ.0 ) THEN 00104 DO 10 I = 1, N 00105 X( I ) = DCMPLX( ONE / DBLE( N ) ) 00106 10 CONTINUE 00107 KASE = 1 00108 ISAVE( 1 ) = 1 00109 RETURN 00110 END IF 00111 * 00112 GO TO ( 20, 40, 70, 90, 120 )ISAVE( 1 ) 00113 * 00114 * ................ ENTRY (ISAVE( 1 ) = 1) 00115 * FIRST ITERATION. X HAS BEEN OVERWRITTEN BY A*X. 00116 * 00117 20 CONTINUE 00118 IF( N.EQ.1 ) THEN 00119 V( 1 ) = X( 1 ) 00120 EST = ABS( V( 1 ) ) 00121 * ... QUIT 00122 GO TO 130 00123 END IF 00124 EST = DZSUM1( N, X, 1 ) 00125 * 00126 DO 30 I = 1, N 00127 ABSXI = ABS( X( I ) ) 00128 IF( ABSXI.GT.SAFMIN ) THEN 00129 X( I ) = DCMPLX( DBLE( X( I ) ) / ABSXI, 00130 $ DIMAG( X( I ) ) / ABSXI ) 00131 ELSE 00132 X( I ) = CONE 00133 END IF 00134 30 CONTINUE 00135 KASE = 2 00136 ISAVE( 1 ) = 2 00137 RETURN 00138 * 00139 * ................ ENTRY (ISAVE( 1 ) = 2) 00140 * FIRST ITERATION. X HAS BEEN OVERWRITTEN BY CTRANS(A)*X. 00141 * 00142 40 CONTINUE 00143 ISAVE( 2 ) = IZMAX1( N, X, 1 ) 00144 ISAVE( 3 ) = 2 00145 * 00146 * MAIN LOOP - ITERATIONS 2,3,...,ITMAX. 00147 * 00148 50 CONTINUE 00149 DO 60 I = 1, N 00150 X( I ) = CZERO 00151 60 CONTINUE 00152 X( ISAVE( 2 ) ) = CONE 00153 KASE = 1 00154 ISAVE( 1 ) = 3 00155 RETURN 00156 * 00157 * ................ ENTRY (ISAVE( 1 ) = 3) 00158 * X HAS BEEN OVERWRITTEN BY A*X. 00159 * 00160 70 CONTINUE 00161 CALL ZCOPY( N, X, 1, V, 1 ) 00162 ESTOLD = EST 00163 EST = DZSUM1( N, V, 1 ) 00164 * 00165 * TEST FOR CYCLING. 00166 IF( EST.LE.ESTOLD ) 00167 $ GO TO 100 00168 * 00169 DO 80 I = 1, N 00170 ABSXI = ABS( X( I ) ) 00171 IF( ABSXI.GT.SAFMIN ) THEN 00172 X( I ) = DCMPLX( DBLE( X( I ) ) / ABSXI, 00173 $ DIMAG( X( I ) ) / ABSXI ) 00174 ELSE 00175 X( I ) = CONE 00176 END IF 00177 80 CONTINUE 00178 KASE = 2 00179 ISAVE( 1 ) = 4 00180 RETURN 00181 * 00182 * ................ ENTRY (ISAVE( 1 ) = 4) 00183 * X HAS BEEN OVERWRITTEN BY CTRANS(A)*X. 00184 * 00185 90 CONTINUE 00186 JLAST = ISAVE( 2 ) 00187 ISAVE( 2 ) = IZMAX1( N, X, 1 ) 00188 IF( ( ABS( X( JLAST ) ).NE.ABS( X( ISAVE( 2 ) ) ) ) .AND. 00189 $ ( ISAVE( 3 ).LT.ITMAX ) ) THEN 00190 ISAVE( 3 ) = ISAVE( 3 ) + 1 00191 GO TO 50 00192 END IF 00193 * 00194 * ITERATION COMPLETE. FINAL STAGE. 00195 * 00196 100 CONTINUE 00197 ALTSGN = ONE 00198 DO 110 I = 1, N 00199 X( I ) = DCMPLX( ALTSGN*( ONE+DBLE( I-1 ) / DBLE( N-1 ) ) ) 00200 ALTSGN = -ALTSGN 00201 110 CONTINUE 00202 KASE = 1 00203 ISAVE( 1 ) = 5 00204 RETURN 00205 * 00206 * ................ ENTRY (ISAVE( 1 ) = 5) 00207 * X HAS BEEN OVERWRITTEN BY A*X. 00208 * 00209 120 CONTINUE 00210 TEMP = TWO*( DZSUM1( N, X, 1 ) / DBLE( 3*N ) ) 00211 IF( TEMP.GT.EST ) THEN 00212 CALL ZCOPY( N, X, 1, V, 1 ) 00213 EST = TEMP 00214 END IF 00215 * 00216 130 CONTINUE 00217 KASE = 0 00218 RETURN 00219 * 00220 * End of ZLACN2 00221 * 00222 END
1.7.3 
