LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
zla_syrpvgrw.f
Go to the documentation of this file.
1 *> \brief \b ZLA_SYRPVGRW computes the reciprocal pivot growth factor norm(A)/norm(U) for a symmetric indefinite matrix.
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download ZLA_SYRPVGRW + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zla_syrpvgrw.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zla_syrpvgrw.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zla_syrpvgrw.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * DOUBLE PRECISION FUNCTION ZLA_SYRPVGRW( UPLO, N, INFO, A, LDA, AF,
22 * LDAF, IPIV, WORK )
23 *
24 * .. Scalar Arguments ..
25 * CHARACTER*1 UPLO
26 * INTEGER N, INFO, LDA, LDAF
27 * ..
28 * .. Array Arguments ..
29 * COMPLEX*16 A( LDA, * ), AF( LDAF, * )
30 * DOUBLE PRECISION WORK( * )
31 * INTEGER IPIV( * )
32 * ..
33 *
34 *
35 *> \par Purpose:
36 * =============
37 *>
38 *> \verbatim
39 *>
40 *>
41 *> ZLA_SYRPVGRW computes the reciprocal pivot growth factor
42 *> norm(A)/norm(U). The "max absolute element" norm is used. If this is
43 *> much less than 1, the stability of the LU factorization of the
44 *> (equilibrated) matrix A could be poor. This also means that the
45 *> solution X, estimated condition numbers, and error bounds could be
46 *> unreliable.
47 *> \endverbatim
48 *
49 * Arguments:
50 * ==========
51 *
52 *> \param[in] UPLO
53 *> \verbatim
54 *> UPLO is CHARACTER*1
55 *> = 'U': Upper triangle of A is stored;
56 *> = 'L': Lower triangle of A is stored.
57 *> \endverbatim
58 *>
59 *> \param[in] N
60 *> \verbatim
61 *> N is INTEGER
62 *> The number of linear equations, i.e., the order of the
63 *> matrix A. N >= 0.
64 *> \endverbatim
65 *>
66 *> \param[in] INFO
67 *> \verbatim
68 *> INFO is INTEGER
69 *> The value of INFO returned from ZSYTRF, .i.e., the pivot in
70 *> column INFO is exactly 0.
71 *> \endverbatim
72 *>
73 *> \param[in] A
74 *> \verbatim
75 *> A is COMPLEX*16 array, dimension (LDA,N)
76 *> On entry, the N-by-N matrix A.
77 *> \endverbatim
78 *>
79 *> \param[in] LDA
80 *> \verbatim
81 *> LDA is INTEGER
82 *> The leading dimension of the array A. LDA >= max(1,N).
83 *> \endverbatim
84 *>
85 *> \param[in] AF
86 *> \verbatim
87 *> AF is COMPLEX*16 array, dimension (LDAF,N)
88 *> The block diagonal matrix D and the multipliers used to
89 *> obtain the factor U or L as computed by ZSYTRF.
90 *> \endverbatim
91 *>
92 *> \param[in] LDAF
93 *> \verbatim
94 *> LDAF is INTEGER
95 *> The leading dimension of the array AF. LDAF >= max(1,N).
96 *> \endverbatim
97 *>
98 *> \param[in] IPIV
99 *> \verbatim
100 *> IPIV is INTEGER array, dimension (N)
101 *> Details of the interchanges and the block structure of D
102 *> as determined by ZSYTRF.
103 *> \endverbatim
104 *>
105 *> \param[out] WORK
106 *> \verbatim
107 *> WORK is DOUBLE PRECISION array, dimension (2*N)
108 *> \endverbatim
109 *
110 * Authors:
111 * ========
112 *
113 *> \author Univ. of Tennessee
114 *> \author Univ. of California Berkeley
115 *> \author Univ. of Colorado Denver
116 *> \author NAG Ltd.
117 *
118 *> \ingroup complex16SYcomputational
119 *
120 * =====================================================================
121  DOUBLE PRECISION FUNCTION zla_syrpvgrw( UPLO, N, INFO, A, LDA, AF,
122  $ LDAF, IPIV, WORK )
123 *
124 * -- LAPACK computational routine --
125 * -- LAPACK is a software package provided by Univ. of Tennessee, --
126 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
127 *
128 * .. Scalar Arguments ..
129  CHARACTER*1 uplo
130  INTEGER n, info, lda, ldaf
131 * ..
132 * .. Array Arguments ..
133  COMPLEX*16 a( lda, * ), af( ldaf, * )
134  DOUBLE PRECISION work( * )
135  INTEGER ipiv( * )
136 * ..
137 *
138 * =====================================================================
139 *
140 * .. Local Scalars ..
141  INTEGER ncols, i, j, k, kp
142  DOUBLE PRECISION amax, umax, rpvgrw, tmp
143  LOGICAL upper
144  COMPLEX*16 zdum
145 * ..
146 * .. Intrinsic Functions ..
147  INTRINSIC abs, real, dimag, max, min
148 * ..
149 * .. External Subroutines ..
150  EXTERNAL lsame
151  LOGICAL lsame
152 * ..
153 * .. Statement Functions ..
154  DOUBLE PRECISION cabs1
155 * ..
156 * .. Statement Function Definitions ..
157  cabs1( zdum ) = abs( dble( zdum ) ) + abs( dimag( zdum ) )
158 * ..
159 * .. Executable Statements ..
160 *
161  upper = lsame( 'Upper', uplo )
162  IF ( info.EQ.0 ) THEN
163  IF ( upper ) THEN
164  ncols = 1
165  ELSE
166  ncols = n
167  END IF
168  ELSE
169  ncols = info
170  END IF
171 
172  rpvgrw = 1.0d+0
173  DO i = 1, 2*n
174  work( i ) = 0.0d+0
175  END DO
176 *
177 * Find the max magnitude entry of each column of A. Compute the max
178 * for all N columns so we can apply the pivot permutation while
179 * looping below. Assume a full factorization is the common case.
180 *
181  IF ( upper ) THEN
182  DO j = 1, n
183  DO i = 1, j
184  work( n+i ) = max( cabs1( a( i, j ) ), work( n+i ) )
185  work( n+j ) = max( cabs1( a( i, j ) ), work( n+j ) )
186  END DO
187  END DO
188  ELSE
189  DO j = 1, n
190  DO i = j, n
191  work( n+i ) = max( cabs1( a( i, j ) ), work( n+i ) )
192  work( n+j ) = max( cabs1( a( i, j ) ), work( n+j ) )
193  END DO
194  END DO
195  END IF
196 *
197 * Now find the max magnitude entry of each column of U or L. Also
198 * permute the magnitudes of A above so they're in the same order as
199 * the factor.
200 *
201 * The iteration orders and permutations were copied from zsytrs.
202 * Calls to SSWAP would be severe overkill.
203 *
204  IF ( upper ) THEN
205  k = n
206  DO WHILE ( k .LT. ncols .AND. k.GT.0 )
207  IF ( ipiv( k ).GT.0 ) THEN
208 ! 1x1 pivot
209  kp = ipiv( k )
210  IF ( kp .NE. k ) THEN
211  tmp = work( n+k )
212  work( n+k ) = work( n+kp )
213  work( n+kp ) = tmp
214  END IF
215  DO i = 1, k
216  work( k ) = max( cabs1( af( i, k ) ), work( k ) )
217  END DO
218  k = k - 1
219  ELSE
220 ! 2x2 pivot
221  kp = -ipiv( k )
222  tmp = work( n+k-1 )
223  work( n+k-1 ) = work( n+kp )
224  work( n+kp ) = tmp
225  DO i = 1, k-1
226  work( k ) = max( cabs1( af( i, k ) ), work( k ) )
227  work( k-1 ) =
228  $ max( cabs1( af( i, k-1 ) ), work( k-1 ) )
229  END DO
230  work( k ) = max( cabs1( af( k, k ) ), work( k ) )
231  k = k - 2
232  END IF
233  END DO
234  k = ncols
235  DO WHILE ( k .LE. n )
236  IF ( ipiv( k ).GT.0 ) THEN
237  kp = ipiv( k )
238  IF ( kp .NE. k ) THEN
239  tmp = work( n+k )
240  work( n+k ) = work( n+kp )
241  work( n+kp ) = tmp
242  END IF
243  k = k + 1
244  ELSE
245  kp = -ipiv( k )
246  tmp = work( n+k )
247  work( n+k ) = work( n+kp )
248  work( n+kp ) = tmp
249  k = k + 2
250  END IF
251  END DO
252  ELSE
253  k = 1
254  DO WHILE ( k .LE. ncols )
255  IF ( ipiv( k ).GT.0 ) THEN
256 ! 1x1 pivot
257  kp = ipiv( k )
258  IF ( kp .NE. k ) THEN
259  tmp = work( n+k )
260  work( n+k ) = work( n+kp )
261  work( n+kp ) = tmp
262  END IF
263  DO i = k, n
264  work( k ) = max( cabs1( af( i, k ) ), work( k ) )
265  END DO
266  k = k + 1
267  ELSE
268 ! 2x2 pivot
269  kp = -ipiv( k )
270  tmp = work( n+k+1 )
271  work( n+k+1 ) = work( n+kp )
272  work( n+kp ) = tmp
273  DO i = k+1, n
274  work( k ) = max( cabs1( af( i, k ) ), work( k ) )
275  work( k+1 ) =
276  $ max( cabs1( af( i, k+1 ) ), work( k+1 ) )
277  END DO
278  work( k ) = max( cabs1( af( k, k ) ), work( k ) )
279  k = k + 2
280  END IF
281  END DO
282  k = ncols
283  DO WHILE ( k .GE. 1 )
284  IF ( ipiv( k ).GT.0 ) THEN
285  kp = ipiv( k )
286  IF ( kp .NE. k ) THEN
287  tmp = work( n+k )
288  work( n+k ) = work( n+kp )
289  work( n+kp ) = tmp
290  END IF
291  k = k - 1
292  ELSE
293  kp = -ipiv( k )
294  tmp = work( n+k )
295  work( n+k ) = work( n+kp )
296  work( n+kp ) = tmp
297  k = k - 2
298  ENDIF
299  END DO
300  END IF
301 *
302 * Compute the *inverse* of the max element growth factor. Dividing
303 * by zero would imply the largest entry of the factor's column is
304 * zero. Than can happen when either the column of A is zero or
305 * massive pivots made the factor underflow to zero. Neither counts
306 * as growth in itself, so simply ignore terms with zero
307 * denominators.
308 *
309  IF ( upper ) THEN
310  DO i = ncols, n
311  umax = work( i )
312  amax = work( n+i )
313  IF ( umax /= 0.0d+0 ) THEN
314  rpvgrw = min( amax / umax, rpvgrw )
315  END IF
316  END DO
317  ELSE
318  DO i = 1, ncols
319  umax = work( i )
320  amax = work( n+i )
321  IF ( umax /= 0.0d+0 ) THEN
322  rpvgrw = min( amax / umax, rpvgrw )
323  END IF
324  END DO
325  END IF
326 
327  zla_syrpvgrw = rpvgrw
328 *
329 * End of ZLA_SYRPVGRW
330 *
331  END
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
double precision function zla_syrpvgrw(UPLO, N, INFO, A, LDA, AF, LDAF, IPIV, WORK)
ZLA_SYRPVGRW computes the reciprocal pivot growth factor norm(A)/norm(U) for a symmetric indefinite m...
Definition: zla_syrpvgrw.f:123