SRC/cla_porpvgrw.f

   1 *> \brief \b CLA_PORPVGRW computes the reciprocal pivot growth factor norm(A)/norm(U) for a symmetric or Hermitian positive-definite matrix.
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download CLA_PORPVGRW + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cla_porpvgrw.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cla_porpvgrw.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cla_porpvgrw.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       REAL FUNCTION CLA_PORPVGRW( UPLO, NCOLS, A, LDA, AF, LDAF, WORK )
  22 *
  23 *       .. Scalar Arguments ..
  24 *       CHARACTER*1        UPLO
  25 *       INTEGER            NCOLS, LDA, LDAF
  26 *       ..
  27 *       .. Array Arguments ..
  28 *       COMPLEX            A( LDA, * ), AF( LDAF, * )
  29 *       REAL               WORK( * )
  30 *       ..
  31 *
  32 *
  33 *> \par Purpose:
  34 *  =============
  35 *>
  36 *> \verbatim
  37 *>
  38 *>
  39 *> CLA_PORPVGRW computes the reciprocal pivot growth factor
  40 *> norm(A)/norm(U). The "max absolute element" norm is used. If this is
  41 *> much less than 1, the stability of the LU factorization of the
  42 *> (equilibrated) matrix A could be poor. This also means that the
  43 *> solution X, estimated condition numbers, and error bounds could be
  44 *> unreliable.
  45 *> \endverbatim
  46 *
  47 *  Arguments:
  48 *  ==========
  49 *
  50 *> \param[in] UPLO
  51 *> \verbatim
  52 *>          UPLO is CHARACTER*1
  53 *>       = 'U':  Upper triangle of A is stored;
  54 *>       = 'L':  Lower triangle of A is stored.
  55 *> \endverbatim
  56 *>
  57 *> \param[in] NCOLS
  58 *> \verbatim
  59 *>          NCOLS is INTEGER
  60 *>     The number of columns of the matrix A. NCOLS >= 0.
  61 *> \endverbatim
  62 *>
  63 *> \param[in] A
  64 *> \verbatim
  65 *>          A is COMPLEX array, dimension (LDA,N)
  66 *>     On entry, the N-by-N matrix A.
  67 *> \endverbatim
  68 *>
  69 *> \param[in] LDA
  70 *> \verbatim
  71 *>          LDA is INTEGER
  72 *>     The leading dimension of the array A.  LDA >= max(1,N).
  73 *> \endverbatim
  74 *>
  75 *> \param[in] AF
  76 *> \verbatim
  77 *>          AF is COMPLEX array, dimension (LDAF,N)
  78 *>     The triangular factor U or L from the Cholesky factorization
  79 *>     A = U**T*U or A = L*L**T, as computed by CPOTRF.
  80 *> \endverbatim
  81 *>
  82 *> \param[in] LDAF
  83 *> \verbatim
  84 *>          LDAF is INTEGER
  85 *>     The leading dimension of the array AF.  LDAF >= max(1,N).
  86 *> \endverbatim
  87 *>
  88 *> \param[in] WORK
  89 *> \verbatim
  90 *>          WORK is REAL array, dimension (2*N)
  91 *> \endverbatim
  92 *
  93 *  Authors:
  94 *  ========
  95 *
  96 *> \author Univ. of Tennessee
  97 *> \author Univ. of California Berkeley
  98 *> \author Univ. of Colorado Denver
  99 *> \author NAG Ltd.
 100 *
 101 *> \date June 2016
 102 *
 103 *> \ingroup complexPOcomputational
 104 *
 105 *  =====================================================================
 106       REAL FUNCTION CLA_PORPVGRW( UPLO, NCOLS, A, LDA, AF, LDAF, WORK )
 107 *
 108 *  -- LAPACK computational routine (version 3.6.1) --
 109 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 110 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 111 *     June 2016
 112 *
 113 *     .. Scalar Arguments ..
 114       CHARACTER*1        UPLO
 115       INTEGER            NCOLS, LDA, LDAF
 116 *     ..
 117 *     .. Array Arguments ..
 118       COMPLEX            A( LDA, * ), AF( LDAF, * )
 119       REAL               WORK( * )
 120 *     ..
 121 *
 122 *  =====================================================================
 123 *
 124 *     .. Local Scalars ..
 125       INTEGER            I, J
 126       REAL               AMAX, UMAX, RPVGRW
 127       LOGICAL            UPPER
 128       COMPLEX            ZDUM
 129 *     ..
 130 *     .. External Functions ..
 131       EXTERNAL           LSAME, CLASET
 132       LOGICAL            LSAME
 133 *     ..
 134 *     .. Intrinsic Functions ..
 135       INTRINSIC          ABS, MAX, MIN, REAL, AIMAG
 136 *     ..
 137 *     .. Statement Functions ..
 138       REAL               CABS1
 139 *     ..
 140 *     .. Statement Function Definitions ..
 141       CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) )
 142 *     ..
 143 *     .. Executable Statements ..
 144       UPPER = LSAME( 'Upper', UPLO )
 145 *
 146 *     SPOTRF will have factored only the NCOLSxNCOLS leading minor, so
 147 *     we restrict the growth search to that minor and use only the first
 148 *     2*NCOLS workspace entries.
 149 *
 150       RPVGRW = 1.0
 151       DO I = 1, 2*NCOLS
 152          WORK( I ) = 0.0
 153       END DO
 154 *
 155 *     Find the max magnitude entry of each column.
 156 *
 157       IF ( UPPER ) THEN
 158          DO J = 1, NCOLS
 159             DO I = 1, J
 160                WORK( NCOLS+J ) =
 161      $              MAX( CABS1( A( I, J ) ), WORK( NCOLS+J ) )
 162             END DO
 163          END DO
 164       ELSE
 165          DO J = 1, NCOLS
 166             DO I = J, NCOLS
 167                WORK( NCOLS+J ) =
 168      $              MAX( CABS1( A( I, J ) ), WORK( NCOLS+J ) )
 169             END DO
 170          END DO
 171       END IF
 172 *
 173 *     Now find the max magnitude entry of each column of the factor in
 174 *     AF.  No pivoting, so no permutations.
 175 *
 176       IF ( LSAME( 'Upper', UPLO ) ) THEN
 177          DO J = 1, NCOLS
 178             DO I = 1, J
 179                WORK( J ) = MAX( CABS1( AF( I, J ) ), WORK( J ) )
 180             END DO
 181          END DO
 182       ELSE
 183          DO J = 1, NCOLS
 184             DO I = J, NCOLS
 185                WORK( J ) = MAX( CABS1( AF( I, J ) ), WORK( J ) )
 186             END DO
 187          END DO
 188       END IF
 189 *
 190 *     Compute the *inverse* of the max element growth factor.  Dividing
 191 *     by zero would imply the largest entry of the factor's column is
 192 *     zero.  Than can happen when either the column of A is zero or
 193 *     massive pivots made the factor underflow to zero.  Neither counts
 194 *     as growth in itself, so simply ignore terms with zero
 195 *     denominators.
 196 *
 197       IF ( LSAME( 'Upper', UPLO ) ) THEN
 198          DO I = 1, NCOLS
 199             UMAX = WORK( I )
 200             AMAX = WORK( NCOLS+I )
 201             IF ( UMAX /= 0.0 ) THEN
 202                RPVGRW = MIN( AMAX / UMAX, RPVGRW )
 203             END IF
 204          END DO
 205       ELSE
 206          DO I = 1, NCOLS
 207             UMAX = WORK( I )
 208             AMAX = WORK( NCOLS+I )
 209             IF ( UMAX /= 0.0 ) THEN
 210                RPVGRW = MIN( AMAX / UMAX, RPVGRW )
 211             END IF
 212          END DO
 213       END IF
 214
 215       CLA_PORPVGRW = RPVGRW
 216       END