001:       SUBROUTINE CLACN2( N, V, X, EST, KASE, ISAVE )
002: *
003: *  -- LAPACK auxiliary routine (version 3.2) --
004: *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
005: *     November 2006
006: *
007: *     .. Scalar Arguments ..
008:       INTEGER            KASE, N
009:       REAL               EST
010: *     ..
011: *     .. Array Arguments ..
012:       INTEGER            ISAVE( 3 )
013:       COMPLEX            V( * ), X( * )
014: *     ..
015: *
016: *  Purpose
017: *  =======
018: *
019: *  CLACN2 estimates the 1-norm of a square, complex matrix A.
020: *  Reverse communication is used for evaluating matrix-vector products.
021: *
022: *  Arguments
023: *  =========
024: *
025: *  N      (input) INTEGER
026: *         The order of the matrix.  N >= 1.
027: *
028: *  V      (workspace) COMPLEX array, dimension (N)
029: *         On the final return, V = A*W,  where  EST = norm(V)/norm(W)
030: *         (W is not returned).
031: *
032: *  X      (input/output) COMPLEX array, dimension (N)
033: *         On an intermediate return, X should be overwritten by
034: *               A * X,   if KASE=1,
035: *               A' * X,  if KASE=2,
036: *         where A' is the conjugate transpose of A, and CLACN2 must be
037: *         re-called with all the other parameters unchanged.
038: *
039: *  EST    (input/output) REAL
040: *         On entry with KASE = 1 or 2 and ISAVE(1) = 3, EST should be
041: *         unchanged from the previous call to CLACN2.
042: *         On exit, EST is an estimate (a lower bound) for norm(A). 
043: *
044: *  KASE   (input/output) INTEGER
045: *         On the initial call to CLACN2, KASE should be 0.
046: *         On an intermediate return, KASE will be 1 or 2, indicating
047: *         whether X should be overwritten by A * X  or A' * X.
048: *         On the final return from CLACN2, KASE will again be 0.
049: *
050: *  ISAVE  (input/output) INTEGER array, dimension (3)
051: *         ISAVE is used to save variables between calls to SLACN2
052: *
053: *  Further Details
054: *  ======= =======
055: *
056: *  Contributed by Nick Higham, University of Manchester.
057: *  Originally named CONEST, dated March 16, 1988.
058: *
059: *  Reference: N.J. Higham, "FORTRAN codes for estimating the one-norm of
060: *  a real or complex matrix, with applications to condition estimation",
061: *  ACM Trans. Math. Soft., vol. 14, no. 4, pp. 381-396, December 1988.
062: *
063: *  Last modified:  April, 1999
064: *
065: *  This is a thread safe version of CLACON, which uses the array ISAVE
066: *  in place of a SAVE statement, as follows:
067: *
068: *     CLACON     CLACN2
069: *      JUMP     ISAVE(1)
070: *      J        ISAVE(2)
071: *      ITER     ISAVE(3)
072: *
073: *  =====================================================================
074: *
075: *     .. Parameters ..
076:       INTEGER              ITMAX
077:       PARAMETER          ( ITMAX = 5 )
078:       REAL                 ONE,         TWO
079:       PARAMETER          ( ONE = 1.0E0, TWO = 2.0E0 )
080:       COMPLEX              CZERO, CONE
081:       PARAMETER          ( CZERO = ( 0.0E0, 0.0E0 ),
082:      $                            CONE = ( 1.0E0, 0.0E0 ) )
083: *     ..
084: *     .. Local Scalars ..
085:       INTEGER            I, JLAST
086:       REAL               ABSXI, ALTSGN, ESTOLD, SAFMIN, TEMP
087: *     ..
088: *     .. External Functions ..
089:       INTEGER            ICMAX1
090:       REAL               SCSUM1, SLAMCH
091:       EXTERNAL           ICMAX1, SCSUM1, SLAMCH
092: *     ..
093: *     .. External Subroutines ..
094:       EXTERNAL           CCOPY
095: *     ..
096: *     .. Intrinsic Functions ..
097:       INTRINSIC          ABS, AIMAG, CMPLX, REAL
098: *     ..
099: *     .. Executable Statements ..
100: *
101:       SAFMIN = SLAMCH( 'Safe minimum' )
102:       IF( KASE.EQ.0 ) THEN
103:          DO 10 I = 1, N
104:             X( I ) = CMPLX( ONE / REAL( N ) )
105:    10    CONTINUE
106:          KASE = 1
107:          ISAVE( 1 ) = 1
108:          RETURN
109:       END IF
110: *
111:       GO TO ( 20, 40, 70, 90, 120 )ISAVE( 1 )
112: *
113: *     ................ ENTRY   (ISAVE( 1 ) = 1)
114: *     FIRST ITERATION.  X HAS BEEN OVERWRITTEN BY A*X.
115: *
116:    20 CONTINUE
117:       IF( N.EQ.1 ) THEN
118:          V( 1 ) = X( 1 )
119:          EST = ABS( V( 1 ) )
120: *        ... QUIT
121:          GO TO 130
122:       END IF
123:       EST = SCSUM1( N, X, 1 )
124: *
125:       DO 30 I = 1, N
126:          ABSXI = ABS( X( I ) )
127:          IF( ABSXI.GT.SAFMIN ) THEN
128:             X( I ) = CMPLX( REAL( X( I ) ) / ABSXI,
129:      $               AIMAG( X( I ) ) / ABSXI )
130:          ELSE
131:             X( I ) = CONE
132:          END IF
133:    30 CONTINUE
134:       KASE = 2
135:       ISAVE( 1 ) = 2
136:       RETURN
137: *
138: *     ................ ENTRY   (ISAVE( 1 ) = 2)
139: *     FIRST ITERATION.  X HAS BEEN OVERWRITTEN BY CTRANS(A)*X.
140: *
141:    40 CONTINUE
142:       ISAVE( 2 ) = ICMAX1( N, X, 1 )
143:       ISAVE( 3 ) = 2
144: *
145: *     MAIN LOOP - ITERATIONS 2,3,...,ITMAX.
146: *
147:    50 CONTINUE
148:       DO 60 I = 1, N
149:          X( I ) = CZERO
150:    60 CONTINUE
151:       X( ISAVE( 2 ) ) = CONE
152:       KASE = 1
153:       ISAVE( 1 ) = 3
154:       RETURN
155: *
156: *     ................ ENTRY   (ISAVE( 1 ) = 3)
157: *     X HAS BEEN OVERWRITTEN BY A*X.
158: *
159:    70 CONTINUE
160:       CALL CCOPY( N, X, 1, V, 1 )
161:       ESTOLD = EST
162:       EST = SCSUM1( N, V, 1 )
163: *
164: *     TEST FOR CYCLING.
165:       IF( EST.LE.ESTOLD )
166:      $   GO TO 100
167: *
168:       DO 80 I = 1, N
169:          ABSXI = ABS( X( I ) )
170:          IF( ABSXI.GT.SAFMIN ) THEN
171:             X( I ) = CMPLX( REAL( X( I ) ) / ABSXI,
172:      $               AIMAG( X( I ) ) / ABSXI )
173:          ELSE
174:             X( I ) = CONE
175:          END IF
176:    80 CONTINUE
177:       KASE = 2
178:       ISAVE( 1 ) = 4
179:       RETURN
180: *
181: *     ................ ENTRY   (ISAVE( 1 ) = 4)
182: *     X HAS BEEN OVERWRITTEN BY CTRANS(A)*X.
183: *
184:    90 CONTINUE
185:       JLAST = ISAVE( 2 )
186:       ISAVE( 2 ) = ICMAX1( N, X, 1 )
187:       IF( ( ABS( X( JLAST ) ).NE.ABS( X( ISAVE( 2 ) ) ) ) .AND.
188:      $    ( ISAVE( 3 ).LT.ITMAX ) ) THEN
189:          ISAVE( 3 ) = ISAVE( 3 ) + 1
190:          GO TO 50
191:       END IF
192: *
193: *     ITERATION COMPLETE.  FINAL STAGE.
194: *
195:   100 CONTINUE
196:       ALTSGN = ONE
197:       DO 110 I = 1, N
198:          X( I ) = CMPLX( ALTSGN*( ONE + REAL( I-1 ) / REAL( N-1 ) ) )
199:          ALTSGN = -ALTSGN
200:   110 CONTINUE
201:       KASE = 1
202:       ISAVE( 1 ) = 5
203:       RETURN
204: *
205: *     ................ ENTRY   (ISAVE( 1 ) = 5)
206: *     X HAS BEEN OVERWRITTEN BY A*X.
207: *
208:   120 CONTINUE
209:       TEMP = TWO*( SCSUM1( N, X, 1 ) / REAL( 3*N ) )
210:       IF( TEMP.GT.EST ) THEN
211:          CALL CCOPY( N, X, 1, V, 1 )
212:          EST = TEMP
213:       END IF
214: *
215:   130 CONTINUE
216:       KASE = 0
217:       RETURN
218: *
219: *     End of CLACN2
220: *
221:       END
222: