1: *> \brief \b ZLACN2 estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vector products.
2: *
3: * =========== DOCUMENTATION ===========
4: *
5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
7: *
8: *> \htmlonly
9: *> Download ZLACN2 + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zlacn2.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zlacn2.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zlacn2.f">
15: *> [TXT]</a>
16: *> \endhtmlonly
17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE ZLACN2( N, V, X, EST, KASE, ISAVE )
22: *
23: * .. Scalar Arguments ..
24: * INTEGER KASE, N
25: * DOUBLE PRECISION EST
26: * ..
27: * .. Array Arguments ..
28: * INTEGER ISAVE( 3 )
29: * COMPLEX*16 V( * ), X( * )
30: * ..
31: *
32: *
33: *> \par Purpose:
34: * =============
35: *>
36: *> \verbatim
37: *>
38: *> ZLACN2 estimates the 1-norm of a square, complex matrix A.
39: *> Reverse communication is used for evaluating matrix-vector products.
40: *> \endverbatim
41: *
42: * Arguments:
43: * ==========
44: *
45: *> \param[in] N
46: *> \verbatim
47: *> N is INTEGER
48: *> The order of the matrix. N >= 1.
49: *> \endverbatim
50: *>
51: *> \param[out] V
52: *> \verbatim
53: *> V is COMPLEX*16 array, dimension (N)
54: *> On the final return, V = A*W, where EST = norm(V)/norm(W)
55: *> (W is not returned).
56: *> \endverbatim
57: *>
58: *> \param[in,out] X
59: *> \verbatim
60: *> X is COMPLEX*16 array, dimension (N)
61: *> On an intermediate return, X should be overwritten by
62: *> A * X, if KASE=1,
63: *> A**H * X, if KASE=2,
64: *> where A**H is the conjugate transpose of A, and ZLACN2 must be
65: *> re-called with all the other parameters unchanged.
66: *> \endverbatim
67: *>
68: *> \param[in,out] EST
69: *> \verbatim
70: *> EST is DOUBLE PRECISION
71: *> On entry with KASE = 1 or 2 and ISAVE(1) = 3, EST should be
72: *> unchanged from the previous call to ZLACN2.
73: *> On exit, EST is an estimate (a lower bound) for norm(A).
74: *> \endverbatim
75: *>
76: *> \param[in,out] KASE
77: *> \verbatim
78: *> KASE is INTEGER
79: *> On the initial call to ZLACN2, KASE should be 0.
80: *> On an intermediate return, KASE will be 1 or 2, indicating
81: *> whether X should be overwritten by A * X or A**H * X.
82: *> On the final return from ZLACN2, KASE will again be 0.
83: *> \endverbatim
84: *>
85: *> \param[in,out] ISAVE
86: *> \verbatim
87: *> ISAVE is INTEGER array, dimension (3)
88: *> ISAVE is used to save variables between calls to ZLACN2
89: *> \endverbatim
90: *
91: * Authors:
92: * ========
93: *
94: *> \author Univ. of Tennessee
95: *> \author Univ. of California Berkeley
96: *> \author Univ. of Colorado Denver
97: *> \author NAG Ltd.
98: *
99: *> \ingroup complex16OTHERauxiliary
100: *
101: *> \par Further Details:
102: * =====================
103: *>
104: *> \verbatim
105: *>
106: *> Originally named CONEST, dated March 16, 1988.
107: *>
108: *> Last modified: April, 1999
109: *>
110: *> This is a thread safe version of ZLACON, which uses the array ISAVE
111: *> in place of a SAVE statement, as follows:
112: *>
113: *> ZLACON ZLACN2
114: *> JUMP ISAVE(1)
115: *> J ISAVE(2)
116: *> ITER ISAVE(3)
117: *> \endverbatim
118: *
119: *> \par Contributors:
120: * ==================
121: *>
122: *> Nick Higham, University of Manchester
123: *
124: *> \par References:
125: * ================
126: *>
127: *> N.J. Higham, "FORTRAN codes for estimating the one-norm of
128: *> a real or complex matrix, with applications to condition estimation",
129: *> ACM Trans. Math. Soft., vol. 14, no. 4, pp. 381-396, December 1988.
130: *>
131: * =====================================================================
132: SUBROUTINE ZLACN2( N, V, X, EST, KASE, ISAVE )
133: *
134: * -- LAPACK auxiliary routine --
135: * -- LAPACK is a software package provided by Univ. of Tennessee, --
136: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
137: *
138: * .. Scalar Arguments ..
139: INTEGER KASE, N
140: DOUBLE PRECISION EST
141: * ..
142: * .. Array Arguments ..
143: INTEGER ISAVE( 3 )
144: COMPLEX*16 V( * ), X( * )
145: * ..
146: *
147: * =====================================================================
148: *
149: * .. Parameters ..
150: INTEGER ITMAX
151: PARAMETER ( ITMAX = 5 )
152: DOUBLE PRECISION ONE, TWO
153: PARAMETER ( ONE = 1.0D0, TWO = 2.0D0 )
154: COMPLEX*16 CZERO, CONE
155: PARAMETER ( CZERO = ( 0.0D0, 0.0D0 ),
156: $ CONE = ( 1.0D0, 0.0D0 ) )
157: * ..
158: * .. Local Scalars ..
159: INTEGER I, JLAST
160: DOUBLE PRECISION ABSXI, ALTSGN, ESTOLD, SAFMIN, TEMP
161: * ..
162: * .. External Functions ..
163: INTEGER IZMAX1
164: DOUBLE PRECISION DLAMCH, DZSUM1
165: EXTERNAL IZMAX1, DLAMCH, DZSUM1
166: * ..
167: * .. External Subroutines ..
168: EXTERNAL ZCOPY
169: * ..
170: * .. Intrinsic Functions ..
171: INTRINSIC ABS, DBLE, DCMPLX, DIMAG
172: * ..
173: * .. Executable Statements ..
174: *
175: SAFMIN = DLAMCH( 'Safe minimum' )
176: IF( KASE.EQ.0 ) THEN
177: DO 10 I = 1, N
178: X( I ) = DCMPLX( ONE / DBLE( N ) )
179: 10 CONTINUE
180: KASE = 1
181: ISAVE( 1 ) = 1
182: RETURN
183: END IF
184: *
185: GO TO ( 20, 40, 70, 90, 120 )ISAVE( 1 )
186: *
187: * ................ ENTRY (ISAVE( 1 ) = 1)
188: * FIRST ITERATION. X HAS BEEN OVERWRITTEN BY A*X.
189: *
190: 20 CONTINUE
191: IF( N.EQ.1 ) THEN
192: V( 1 ) = X( 1 )
193: EST = ABS( V( 1 ) )
194: * ... QUIT
195: GO TO 130
196: END IF
197: EST = DZSUM1( N, X, 1 )
198: *
199: DO 30 I = 1, N
200: ABSXI = ABS( X( I ) )
201: IF( ABSXI.GT.SAFMIN ) THEN
202: X( I ) = DCMPLX( DBLE( X( I ) ) / ABSXI,
203: $ DIMAG( X( I ) ) / ABSXI )
204: ELSE
205: X( I ) = CONE
206: END IF
207: 30 CONTINUE
208: KASE = 2
209: ISAVE( 1 ) = 2
210: RETURN
211: *
212: * ................ ENTRY (ISAVE( 1 ) = 2)
213: * FIRST ITERATION. X HAS BEEN OVERWRITTEN BY CTRANS(A)*X.
214: *
215: 40 CONTINUE
216: ISAVE( 2 ) = IZMAX1( N, X, 1 )
217: ISAVE( 3 ) = 2
218: *
219: * MAIN LOOP - ITERATIONS 2,3,...,ITMAX.
220: *
221: 50 CONTINUE
222: DO 60 I = 1, N
223: X( I ) = CZERO
224: 60 CONTINUE
225: X( ISAVE( 2 ) ) = CONE
226: KASE = 1
227: ISAVE( 1 ) = 3
228: RETURN
229: *
230: * ................ ENTRY (ISAVE( 1 ) = 3)
231: * X HAS BEEN OVERWRITTEN BY A*X.
232: *
233: 70 CONTINUE
234: CALL ZCOPY( N, X, 1, V, 1 )
235: ESTOLD = EST
236: EST = DZSUM1( N, V, 1 )
237: *
238: * TEST FOR CYCLING.
239: IF( EST.LE.ESTOLD )
240: $ GO TO 100
241: *
242: DO 80 I = 1, N
243: ABSXI = ABS( X( I ) )
244: IF( ABSXI.GT.SAFMIN ) THEN
245: X( I ) = DCMPLX( DBLE( X( I ) ) / ABSXI,
246: $ DIMAG( X( I ) ) / ABSXI )
247: ELSE
248: X( I ) = CONE
249: END IF
250: 80 CONTINUE
251: KASE = 2
252: ISAVE( 1 ) = 4
253: RETURN
254: *
255: * ................ ENTRY (ISAVE( 1 ) = 4)
256: * X HAS BEEN OVERWRITTEN BY CTRANS(A)*X.
257: *
258: 90 CONTINUE
259: JLAST = ISAVE( 2 )
260: ISAVE( 2 ) = IZMAX1( N, X, 1 )
261: IF( ( ABS( X( JLAST ) ).NE.ABS( X( ISAVE( 2 ) ) ) ) .AND.
262: $ ( ISAVE( 3 ).LT.ITMAX ) ) THEN
263: ISAVE( 3 ) = ISAVE( 3 ) + 1
264: GO TO 50
265: END IF
266: *
267: * ITERATION COMPLETE. FINAL STAGE.
268: *
269: 100 CONTINUE
270: ALTSGN = ONE
271: DO 110 I = 1, N
272: X( I ) = DCMPLX( ALTSGN*( ONE+DBLE( I-1 ) / DBLE( N-1 ) ) )
273: ALTSGN = -ALTSGN
274: 110 CONTINUE
275: KASE = 1
276: ISAVE( 1 ) = 5
277: RETURN
278: *
279: * ................ ENTRY (ISAVE( 1 ) = 5)
280: * X HAS BEEN OVERWRITTEN BY A*X.
281: *
282: 120 CONTINUE
283: TEMP = TWO*( DZSUM1( N, X, 1 ) / DBLE( 3*N ) )
284: IF( TEMP.GT.EST ) THEN
285: CALL ZCOPY( N, X, 1, V, 1 )
286: EST = TEMP
287: END IF
288: *
289: 130 CONTINUE
290: KASE = 0
291: RETURN
292: *
293: * End of ZLACN2
294: *
295: END
CVSweb interface <joel.bertrand@systella.fr>