LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
clacon.f
Go to the documentation of this file.
1 *> \brief \b CLACON 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 CLACON + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/clacon.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/clacon.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/clacon.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE CLACON( N, V, X, EST, KASE )
22 *
23 * .. Scalar Arguments ..
24 * INTEGER KASE, N
25 * REAL EST
26 * ..
27 * .. Array Arguments ..
28 * COMPLEX V( N ), X( N )
29 * ..
30 *
31 *
32 *> \par Purpose:
33 * =============
34 *>
35 *> \verbatim
36 *>
37 *> CLACON estimates the 1-norm of a square, complex matrix A.
38 *> Reverse communication is used for evaluating matrix-vector products.
39 *> \endverbatim
40 *
41 * Arguments:
42 * ==========
43 *
44 *> \param[in] N
45 *> \verbatim
46 *> N is INTEGER
47 *> The order of the matrix. N >= 1.
48 *> \endverbatim
49 *>
50 *> \param[out] V
51 *> \verbatim
52 *> V is COMPLEX array, dimension (N)
53 *> On the final return, V = A*W, where EST = norm(V)/norm(W)
54 *> (W is not returned).
55 *> \endverbatim
56 *>
57 *> \param[in,out] X
58 *> \verbatim
59 *> X is COMPLEX array, dimension (N)
60 *> On an intermediate return, X should be overwritten by
61 *> A * X, if KASE=1,
62 *> A**H * X, if KASE=2,
63 *> where A**H is the conjugate transpose of A, and CLACON must be
64 *> re-called with all the other parameters unchanged.
65 *> \endverbatim
66 *>
67 *> \param[in,out] EST
68 *> \verbatim
69 *> EST is REAL
70 *> On entry with KASE = 1 or 2 and JUMP = 3, EST should be
71 *> unchanged from the previous call to CLACON.
72 *> On exit, EST is an estimate (a lower bound) for norm(A).
73 *> \endverbatim
74 *>
75 *> \param[in,out] KASE
76 *> \verbatim
77 *> KASE is INTEGER
78 *> On the initial call to CLACON, KASE should be 0.
79 *> On an intermediate return, KASE will be 1 or 2, indicating
80 *> whether X should be overwritten by A * X or A**H * X.
81 *> On the final return from CLACON, KASE will again be 0.
82 *> \endverbatim
83 *
84 * Authors:
85 * ========
86 *
87 *> \author Univ. of Tennessee
88 *> \author Univ. of California Berkeley
89 *> \author Univ. of Colorado Denver
90 *> \author NAG Ltd.
91 *
92 *> \date September 2012
93 *
94 *> \ingroup complexOTHERauxiliary
95 *
96 *> \par Further Details:
97 * =====================
98 *>
99 *> Originally named CONEST, dated March 16, 1988. \n
100 *> Last modified: April, 1999
101 *
102 *> \par Contributors:
103 * ==================
104 *>
105 *> Nick Higham, University of Manchester
106 *
107 *> \par References:
108 * ================
109 *>
110 *> N.J. Higham, "FORTRAN codes for estimating the one-norm of
111 *> a real or complex matrix, with applications to condition estimation",
112 *> ACM Trans. Math. Soft., vol. 14, no. 4, pp. 381-396, December 1988.
113 *>
114 * =====================================================================
115  SUBROUTINE clacon( N, V, X, EST, KASE )
116 *
117 * -- LAPACK auxiliary routine (version 3.4.2) --
118 * -- LAPACK is a software package provided by Univ. of Tennessee, --
119 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
120 * September 2012
121 *
122 * .. Scalar Arguments ..
123  INTEGER kase, n
124  REAL est
125 * ..
126 * .. Array Arguments ..
127  COMPLEX v( n ), x( n )
128 * ..
129 *
130 * =====================================================================
131 *
132 * .. Parameters ..
133  INTEGER itmax
134  parameter( itmax = 5 )
135  REAL one, two
136  parameter( one = 1.0e0, two = 2.0e0 )
137  COMPLEX czero, cone
138  parameter( czero = ( 0.0e0, 0.0e0 ),
139  $ cone = ( 1.0e0, 0.0e0 ) )
140 * ..
141 * .. Local Scalars ..
142  INTEGER i, iter, j, jlast, jump
143  REAL absxi, altsgn, estold, safmin, temp
144 * ..
145 * .. External Functions ..
146  INTEGER icmax1
147  REAL scsum1, slamch
148  EXTERNAL icmax1, scsum1, slamch
149 * ..
150 * .. External Subroutines ..
151  EXTERNAL ccopy
152 * ..
153 * .. Intrinsic Functions ..
154  INTRINSIC abs, aimag, cmplx, real
155 * ..
156 * .. Save statement ..
157  SAVE
158 * ..
159 * .. Executable Statements ..
160 *
161  safmin = slamch( 'Safe minimum' )
162  IF( kase.EQ.0 ) THEN
163  DO 10 i = 1, n
164  x( i ) = cmplx( one / REAL( N ) )
165  10 continue
166  kase = 1
167  jump = 1
168  return
169  END IF
170 *
171  go to( 20, 40, 70, 90, 120 )jump
172 *
173 * ................ ENTRY (JUMP = 1)
174 * FIRST ITERATION. X HAS BEEN OVERWRITTEN BY A*X.
175 *
176  20 continue
177  IF( n.EQ.1 ) THEN
178  v( 1 ) = x( 1 )
179  est = abs( v( 1 ) )
180 * ... QUIT
181  go to 130
182  END IF
183  est = scsum1( n, x, 1 )
184 *
185  DO 30 i = 1, n
186  absxi = abs( x( i ) )
187  IF( absxi.GT.safmin ) THEN
188  x( i ) = cmplx( REAL( X( I ) ) / absxi,
189  $ aimag( x( i ) ) / absxi )
190  ELSE
191  x( i ) = cone
192  END IF
193  30 continue
194  kase = 2
195  jump = 2
196  return
197 *
198 * ................ ENTRY (JUMP = 2)
199 * FIRST ITERATION. X HAS BEEN OVERWRITTEN BY CTRANS(A)*X.
200 *
201  40 continue
202  j = icmax1( n, x, 1 )
203  iter = 2
204 *
205 * MAIN LOOP - ITERATIONS 2,3,...,ITMAX.
206 *
207  50 continue
208  DO 60 i = 1, n
209  x( i ) = czero
210  60 continue
211  x( j ) = cone
212  kase = 1
213  jump = 3
214  return
215 *
216 * ................ ENTRY (JUMP = 3)
217 * X HAS BEEN OVERWRITTEN BY A*X.
218 *
219  70 continue
220  CALL ccopy( n, x, 1, v, 1 )
221  estold = est
222  est = scsum1( n, v, 1 )
223 *
224 * TEST FOR CYCLING.
225  IF( est.LE.estold )
226  $ go to 100
227 *
228  DO 80 i = 1, n
229  absxi = abs( x( i ) )
230  IF( absxi.GT.safmin ) THEN
231  x( i ) = cmplx( REAL( X( I ) ) / absxi,
232  $ aimag( x( i ) ) / absxi )
233  ELSE
234  x( i ) = cone
235  END IF
236  80 continue
237  kase = 2
238  jump = 4
239  return
240 *
241 * ................ ENTRY (JUMP = 4)
242 * X HAS BEEN OVERWRITTEN BY CTRANS(A)*X.
243 *
244  90 continue
245  jlast = j
246  j = icmax1( n, x, 1 )
247  IF( ( abs( x( jlast ) ).NE.abs( x( j ) ) ) .AND.
248  $ ( iter.LT.itmax ) ) THEN
249  iter = iter + 1
250  go to 50
251  END IF
252 *
253 * ITERATION COMPLETE. FINAL STAGE.
254 *
255  100 continue
256  altsgn = one
257  DO 110 i = 1, n
258  x( i ) = cmplx( altsgn*( one+REAL( I-1 ) / REAL( N-1 ) ) )
259  altsgn = -altsgn
260  110 continue
261  kase = 1
262  jump = 5
263  return
264 *
265 * ................ ENTRY (JUMP = 5)
266 * X HAS BEEN OVERWRITTEN BY A*X.
267 *
268  120 continue
269  temp = two*( scsum1( n, x, 1 ) / REAL( 3*N ) )
270  IF( temp.GT.est ) THEN
271  CALL ccopy( n, x, 1, v, 1 )
272  est = temp
273  END IF
274 *
275  130 continue
276  kase = 0
277  return
278 *
279 * End of CLACON
280 *
281  END