001:       SUBROUTINE SLACN2( N, V, X, ISGN, 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            ISGN( * ), ISAVE( 3 )
013:       REAL               V( * ), X( * )
014: *     ..
015: *
016: *  Purpose
017: *  =======
018: *
019: *  SLACN2 estimates the 1-norm of a square, real 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) REAL 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) REAL 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: *         and SLACN2 must be re-called with all the other parameters
037: *         unchanged.
038: *
039: *  ISGN   (workspace) INTEGER array, dimension (N)
040: *
041: *  EST    (input/output) REAL
042: *         On entry with KASE = 1 or 2 and ISAVE(1) = 3, EST should be
043: *         unchanged from the previous call to SLACN2.
044: *         On exit, EST is an estimate (a lower bound) for norm(A). 
045: *
046: *  KASE   (input/output) INTEGER
047: *         On the initial call to SLACN2, KASE should be 0.
048: *         On an intermediate return, KASE will be 1 or 2, indicating
049: *         whether X should be overwritten by A * X  or A' * X.
050: *         On the final return from SLACN2, KASE will again be 0.
051: *
052: *  ISAVE  (input/output) INTEGER array, dimension (3)
053: *         ISAVE is used to save variables between calls to SLACN2
054: *
055: *  Further Details
056: *  ======= =======
057: *
058: *  Contributed by Nick Higham, University of Manchester.
059: *  Originally named SONEST, dated March 16, 1988.
060: *
061: *  Reference: N.J. Higham, "FORTRAN codes for estimating the one-norm of
062: *  a real or complex matrix, with applications to condition estimation",
063: *  ACM Trans. Math. Soft., vol. 14, no. 4, pp. 381-396, December 1988.
064: *
065: *  This is a thread safe version of SLACON, which uses the array ISAVE
066: *  in place of a SAVE statement, as follows:
067: *
068: *     SLACON     SLACN2
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               ZERO, ONE, TWO
079:       PARAMETER          ( ZERO = 0.0E+0, ONE = 1.0E+0, TWO = 2.0E+0 )
080: *     ..
081: *     .. Local Scalars ..
082:       INTEGER            I, JLAST
083:       REAL               ALTSGN, ESTOLD, TEMP
084: *     ..
085: *     .. External Functions ..
086:       INTEGER            ISAMAX
087:       REAL               SASUM
088:       EXTERNAL           ISAMAX, SASUM
089: *     ..
090: *     .. External Subroutines ..
091:       EXTERNAL           SCOPY
092: *     ..
093: *     .. Intrinsic Functions ..
094:       INTRINSIC          ABS, NINT, REAL, SIGN
095: *     ..
096: *     .. Executable Statements ..
097: *
098:       IF( KASE.EQ.0 ) THEN
099:          DO 10 I = 1, N
100:             X( I ) = ONE / REAL( N )
101:    10    CONTINUE
102:          KASE = 1
103:          ISAVE( 1 ) = 1
104:          RETURN
105:       END IF
106: *
107:       GO TO ( 20, 40, 70, 110, 140 )ISAVE( 1 )
108: *
109: *     ................ ENTRY   (ISAVE( 1 ) = 1)
110: *     FIRST ITERATION.  X HAS BEEN OVERWRITTEN BY A*X.
111: *
112:    20 CONTINUE
113:       IF( N.EQ.1 ) THEN
114:          V( 1 ) = X( 1 )
115:          EST = ABS( V( 1 ) )
116: *        ... QUIT
117:          GO TO 150
118:       END IF
119:       EST = SASUM( N, X, 1 )
120: *
121:       DO 30 I = 1, N
122:          X( I ) = SIGN( ONE, X( I ) )
123:          ISGN( I ) = NINT( X( I ) )
124:    30 CONTINUE
125:       KASE = 2
126:       ISAVE( 1 ) = 2
127:       RETURN
128: *
129: *     ................ ENTRY   (ISAVE( 1 ) = 2)
130: *     FIRST ITERATION.  X HAS BEEN OVERWRITTEN BY TRANSPOSE(A)*X.
131: *
132:    40 CONTINUE
133:       ISAVE( 2 ) = ISAMAX( N, X, 1 )
134:       ISAVE( 3 ) = 2
135: *
136: *     MAIN LOOP - ITERATIONS 2,3,...,ITMAX.
137: *
138:    50 CONTINUE
139:       DO 60 I = 1, N
140:          X( I ) = ZERO
141:    60 CONTINUE
142:       X( ISAVE( 2 ) ) = ONE
143:       KASE = 1
144:       ISAVE( 1 ) = 3
145:       RETURN
146: *
147: *     ................ ENTRY   (ISAVE( 1 ) = 3)
148: *     X HAS BEEN OVERWRITTEN BY A*X.
149: *
150:    70 CONTINUE
151:       CALL SCOPY( N, X, 1, V, 1 )
152:       ESTOLD = EST
153:       EST = SASUM( N, V, 1 )
154:       DO 80 I = 1, N
155:          IF( NINT( SIGN( ONE, X( I ) ) ).NE.ISGN( I ) )
156:      $      GO TO 90
157:    80 CONTINUE
158: *     REPEATED SIGN VECTOR DETECTED, HENCE ALGORITHM HAS CONVERGED.
159:       GO TO 120
160: *
161:    90 CONTINUE
162: *     TEST FOR CYCLING.
163:       IF( EST.LE.ESTOLD )
164:      $   GO TO 120
165: *
166:       DO 100 I = 1, N
167:          X( I ) = SIGN( ONE, X( I ) )
168:          ISGN( I ) = NINT( X( I ) )
169:   100 CONTINUE
170:       KASE = 2
171:       ISAVE( 1 ) = 4
172:       RETURN
173: *
174: *     ................ ENTRY   (ISAVE( 1 ) = 4)
175: *     X HAS BEEN OVERWRITTEN BY TRANSPOSE(A)*X.
176: *
177:   110 CONTINUE
178:       JLAST = ISAVE( 2 )
179:       ISAVE( 2 ) = ISAMAX( N, X, 1 )
180:       IF( ( X( JLAST ).NE.ABS( X( ISAVE( 2 ) ) ) ) .AND.
181:      $    ( ISAVE( 3 ).LT.ITMAX ) ) THEN
182:          ISAVE( 3 ) = ISAVE( 3 ) + 1
183:          GO TO 50
184:       END IF
185: *
186: *     ITERATION COMPLETE.  FINAL STAGE.
187: *
188:   120 CONTINUE
189:       ALTSGN = ONE
190:       DO 130 I = 1, N
191:          X( I ) = ALTSGN*( ONE+REAL( I-1 ) / REAL( N-1 ) )
192:          ALTSGN = -ALTSGN
193:   130 CONTINUE
194:       KASE = 1
195:       ISAVE( 1 ) = 5
196:       RETURN
197: *
198: *     ................ ENTRY   (ISAVE( 1 ) = 5)
199: *     X HAS BEEN OVERWRITTEN BY A*X.
200: *
201:   140 CONTINUE
202:       TEMP = TWO*( SASUM( N, X, 1 ) / REAL( 3*N ) )
203:       IF( TEMP.GT.EST ) THEN
204:          CALL SCOPY( N, X, 1, V, 1 )
205:          EST = TEMP
206:       END IF
207: *
208:   150 CONTINUE
209:       KASE = 0
210:       RETURN
211: *
212: *     End of SLACN2
213: *
214:       END
215: