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