LAPACK 3.12.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches
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
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 la_herpvgrw
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
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...
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48