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