LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
real function cla_syrpvgrw ( character*1  UPLO,
integer  N,
integer  INFO,
complex, dimension( lda, * )  A,
integer  LDA,
complex, dimension( ldaf, * )  AF,
integer  LDAF,
integer, dimension( * )  IPIV,
real, dimension( * )  WORK 
)

CLA_SYRPVGRW computes the reciprocal pivot growth factor norm(A)/norm(U) for a symmetric indefinite matrix.

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

Purpose: 
 CLA_SYRPVGRW computes the reciprocal pivot growth factor
 norm(A)/norm(U). The "max absolute element" 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 CSYTRF, .i.e., the pivot in
     column INFO is exactly 0.
[in]A
          A is COMPLEX 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 array, dimension (LDAF,N)
     The block diagonal matrix D and the multipliers used to
     obtain the factor U or L as computed by CSYTRF.
[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 CSYTRF.
[in]WORK
          WORK is REAL array, dimension (2*N)
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2015

Definition at line 125 of file cla_syrpvgrw.f.

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

Here is the call graph for this function:

Here is the caller graph for this function: