LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ zla_herpvgrw()

double precision function zla_herpvgrw ( character*1 uplo,
integer n,
integer info,
complex*16, dimension( lda, * ) a,
integer lda,
complex*16, dimension( ldaf, * ) af,
integer ldaf,
integer, dimension( * ) ipiv,
double precision, dimension( * ) work )

ZLA_HERPVGRW

Download ZLA_HERPVGRW + dependencies [TGZ] [ZIP] [TXT]

Purpose:
!>
!>
!> ZLA_HERPVGRW computes the reciprocal pivot growth factor
!> norm(A)/norm(U). The  norm is used. If this is
!> much less than 1, the stability of the LU factorization of the
!> (equilibrated) matrix A could be poor. This also means that the
!> solution X, estimated condition numbers, and error bounds could be
!> unreliable.
!> 
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>       = 'U':  Upper triangle of A is stored;
!>       = 'L':  Lower triangle of A is stored.
!> 
[in]N
!>          N is INTEGER
!>     The number of linear equations, i.e., the order of the
!>     matrix A.  N >= 0.
!> 
[in]INFO
!>          INFO is INTEGER
!>     The value of INFO returned from ZHETRF, .i.e., the pivot in
!>     column INFO is exactly 0.
!> 
[in]A
!>          A is COMPLEX*16 array, dimension (LDA,N)
!>     On entry, the N-by-N matrix A.
!> 
[in]LDA
!>          LDA is INTEGER
!>     The leading dimension of the array A.  LDA >= max(1,N).
!> 
[in]AF
!>          AF is COMPLEX*16 array, dimension (LDAF,N)
!>     The block diagonal matrix D and the multipliers used to
!>     obtain the factor U or L as computed by ZHETRF.
!> 
[in]LDAF
!>          LDAF is INTEGER
!>     The leading dimension of the array AF.  LDAF >= max(1,N).
!> 
[in]IPIV
!>          IPIV is INTEGER array, dimension (N)
!>     Details of the interchanges and the block structure of D
!>     as determined by ZHETRF.
!> 
[out]WORK
!>          WORK is DOUBLE PRECISION array, dimension (2*N)
!> 
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 119 of file zla_herpvgrw.f.

122*
123* -- LAPACK computational routine --
124* -- LAPACK is a software package provided by Univ. of Tennessee, --
125* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
126*
127* .. Scalar Arguments ..
128 CHARACTER*1 UPLO
129 INTEGER N, INFO, LDA, LDAF
130* ..
131* .. Array Arguments ..
132 INTEGER IPIV( * )
133 COMPLEX*16 A( LDA, * ), AF( LDAF, * )
134 DOUBLE PRECISION WORK( * )
135* ..
136*
137* =====================================================================
138*
139* .. Local Scalars ..
140 INTEGER NCOLS, I, J, K, KP
141 DOUBLE PRECISION AMAX, UMAX, RPVGRW, TMP
142 LOGICAL UPPER, LSAME
143 COMPLEX*16 ZDUM
144* ..
145* .. External Functions ..
146 EXTERNAL lsame
147* ..
148* .. Intrinsic Functions ..
149 INTRINSIC abs, real, dimag, max, min
150* ..
151* .. Statement Functions ..
152 DOUBLE PRECISION CABS1
153* ..
154* .. Statement Function Definitions ..
155 cabs1( zdum ) = abs( dble( zdum ) ) + abs( dimag( zdum ) )
156* ..
157* .. Executable Statements ..
158*
159 upper = lsame( 'Upper', uplo )
160 IF ( info.EQ.0 ) THEN
161 IF (upper) THEN
162 ncols = 1
163 ELSE
164 ncols = n
165 END IF
166 ELSE
167 ncols = info
168 END IF
169
170 rpvgrw = 1.0d+0
171 DO i = 1, 2*n
172 work( i ) = 0.0d+0
173 END DO
174*
175* Find the max magnitude entry of each column of A. Compute the max
176* for all N columns so we can apply the pivot permutation while
177* looping below. Assume a full factorization is the common case.
178*
179 IF ( upper ) THEN
180 DO j = 1, n
181 DO i = 1, j
182 work( n+i ) = max( cabs1( a( i,j ) ), work( n+i ) )
183 work( n+j ) = max( cabs1( a( i,j ) ), work( n+j ) )
184 END DO
185 END DO
186 ELSE
187 DO j = 1, n
188 DO i = j, n
189 work( n+i ) = max( cabs1( a( i, j ) ), work( n+i ) )
190 work( n+j ) = max( cabs1( a( i, j ) ), work( n+j ) )
191 END DO
192 END DO
193 END IF
194*
195* Now find the max magnitude entry of each column of U or L. Also
196* permute the magnitudes of A above so they're in the same order as
197* the factor.
198*
199* The iteration orders and permutations were copied from zsytrs.
200* Calls to SSWAP would be severe overkill.
201*
202 IF ( upper ) THEN
203 k = n
204 DO WHILE ( k .LT. ncols .AND. k.GT.0 )
205 IF ( ipiv( k ).GT.0 ) THEN
206! 1x1 pivot
207 kp = ipiv( k )
208 IF ( kp .NE. k ) THEN
209 tmp = work( n+k )
210 work( n+k ) = work( n+kp )
211 work( n+kp ) = tmp
212 END IF
213 DO i = 1, k
214 work( k ) = max( cabs1( af( i, k ) ), work( k ) )
215 END DO
216 k = k - 1
217 ELSE
218! 2x2 pivot
219 kp = -ipiv( k )
220 tmp = work( n+k-1 )
221 work( n+k-1 ) = work( n+kp )
222 work( n+kp ) = tmp
223 DO i = 1, k-1
224 work( k ) = max( cabs1( af( i, k ) ), work( k ) )
225 work( k-1 ) =
226 $ max( cabs1( af( i, k-1 ) ), work( k-1 ) )
227 END DO
228 work( k ) = max( cabs1( af( k, k ) ), work( k ) )
229 k = k - 2
230 END IF
231 END DO
232 k = ncols
233 DO WHILE ( k .LE. n )
234 IF ( ipiv( k ).GT.0 ) THEN
235 kp = ipiv( k )
236 IF ( kp .NE. k ) THEN
237 tmp = work( n+k )
238 work( n+k ) = work( n+kp )
239 work( n+kp ) = tmp
240 END IF
241 k = k + 1
242 ELSE
243 kp = -ipiv( k )
244 tmp = work( n+k )
245 work( n+k ) = work( n+kp )
246 work( n+kp ) = tmp
247 k = k + 2
248 END IF
249 END DO
250 ELSE
251 k = 1
252 DO WHILE ( k .LE. ncols )
253 IF ( ipiv( k ).GT.0 ) THEN
254! 1x1 pivot
255 kp = ipiv( k )
256 IF ( kp .NE. k ) THEN
257 tmp = work( n+k )
258 work( n+k ) = work( n+kp )
259 work( n+kp ) = tmp
260 END IF
261 DO i = k, n
262 work( k ) = max( cabs1( af( i, k ) ), work( k ) )
263 END DO
264 k = k + 1
265 ELSE
266! 2x2 pivot
267 kp = -ipiv( k )
268 tmp = work( n+k+1 )
269 work( n+k+1 ) = work( n+kp )
270 work( n+kp ) = tmp
271 DO i = k+1, n
272 work( k ) = max( cabs1( af( i, k ) ), work( k ) )
273 work( k+1 ) =
274 $ max( cabs1( af( i, k+1 ) ) , work( k+1 ) )
275 END DO
276 work(k) = max( cabs1( af( k, k ) ), work( k ) )
277 k = k + 2
278 END IF
279 END DO
280 k = ncols
281 DO WHILE ( k .GE. 1 )
282 IF ( ipiv( k ).GT.0 ) THEN
283 kp = ipiv( k )
284 IF ( kp .NE. k ) THEN
285 tmp = work( n+k )
286 work( n+k ) = work( n+kp )
287 work( n+kp ) = tmp
288 END IF
289 k = k - 1
290 ELSE
291 kp = -ipiv( k )
292 tmp = work( n+k )
293 work( n+k ) = work( n+kp )
294 work( n+kp ) = tmp
295 k = k - 2
296 ENDIF
297 END DO
298 END IF
299*
300* Compute the *inverse* of the max element growth factor. Dividing
301* by zero would imply the largest entry of the factor's column is
302* zero. Than can happen when either the column of A is zero or
303* massive pivots made the factor underflow to zero. Neither counts
304* as growth in itself, so simply ignore terms with zero
305* denominators.
306*
307 IF ( upper ) THEN
308 DO i = ncols, n
309 umax = work( i )
310 amax = work( n+i )
311 IF ( umax /= 0.0d+0 ) THEN
312 rpvgrw = min( amax / umax, rpvgrw )
313 END IF
314 END DO
315 ELSE
316 DO i = 1, ncols
317 umax = work( i )
318 amax = work( n+i )
319 IF ( umax /= 0.0d+0 ) THEN
320 rpvgrw = min( amax / umax, rpvgrw )
321 END IF
322 END DO
323 END IF
324
325 zla_herpvgrw = rpvgrw
326*
327* End of ZLA_HERPVGRW
328*
double precision function zla_herpvgrw(uplo, n, info, a, lda, af, ldaf, ipiv, work)
ZLA_HERPVGRW
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
Here is the call graph for this function:
Here is the caller graph for this function: