SRC/slaqsp.f

   1 *> \brief \b SLAQSP scales a symmetric/Hermitian matrix in packed storage, using scaling factors computed by sppequ.
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download SLAQSP + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/slaqsp.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/slaqsp.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/slaqsp.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE SLAQSP( UPLO, N, AP, S, SCOND, AMAX, EQUED )
  22 *
  23 *       .. Scalar Arguments ..
  24 *       CHARACTER          EQUED, UPLO
  25 *       INTEGER            N
  26 *       REAL               AMAX, SCOND
  27 *       ..
  28 *       .. Array Arguments ..
  29 *       REAL               AP( * ), S( * )
  30 *       ..
  31 *
  32 *
  33 *> \par Purpose:
  34 *  =============
  35 *>
  36 *> \verbatim
  37 *>
  38 *> SLAQSP equilibrates a symmetric matrix A using the scaling factors
  39 *> in the vector S.
  40 *> \endverbatim
  41 *
  42 *  Arguments:
  43 *  ==========
  44 *
  45 *> \param[in] UPLO
  46 *> \verbatim
  47 *>          UPLO is CHARACTER*1
  48 *>          Specifies whether the upper or lower triangular part of the
  49 *>          symmetric matrix A is stored.
  50 *>          = 'U':  Upper triangular
  51 *>          = 'L':  Lower triangular
  52 *> \endverbatim
  53 *>
  54 *> \param[in] N
  55 *> \verbatim
  56 *>          N is INTEGER
  57 *>          The order of the matrix A.  N >= 0.
  58 *> \endverbatim
  59 *>
  60 *> \param[in,out] AP
  61 *> \verbatim
  62 *>          AP is REAL array, dimension (N*(N+1)/2)
  63 *>          On entry, the upper or lower triangle of the symmetric matrix
  64 *>          A, packed columnwise in a linear array.  The j-th column of A
  65 *>          is stored in the array AP as follows:
  66 *>          if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
  67 *>          if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n.
  68 *>
  69 *>          On exit, the equilibrated matrix:  diag(S) * A * diag(S), in
  70 *>          the same storage format as A.
  71 *> \endverbatim
  72 *>
  73 *> \param[in] S
  74 *> \verbatim
  75 *>          S is REAL array, dimension (N)
  76 *>          The scale factors for A.
  77 *> \endverbatim
  78 *>
  79 *> \param[in] SCOND
  80 *> \verbatim
  81 *>          SCOND is REAL
  82 *>          Ratio of the smallest S(i) to the largest S(i).
  83 *> \endverbatim
  84 *>
  85 *> \param[in] AMAX
  86 *> \verbatim
  87 *>          AMAX is REAL
  88 *>          Absolute value of largest matrix entry.
  89 *> \endverbatim
  90 *>
  91 *> \param[out] EQUED
  92 *> \verbatim
  93 *>          EQUED is CHARACTER*1
  94 *>          Specifies whether or not equilibration was done.
  95 *>          = 'N':  No equilibration.
  96 *>          = 'Y':  Equilibration was done, i.e., A has been replaced by
  97 *>                  diag(S) * A * diag(S).
  98 *> \endverbatim
  99 *
 100 *> \par Internal Parameters:
 101 *  =========================
 102 *>
 103 *> \verbatim
 104 *>  THRESH is a threshold value used to decide if scaling should be done
 105 *>  based on the ratio of the scaling factors.  If SCOND < THRESH,
 106 *>  scaling is done.
 107 *>
 108 *>  LARGE and SMALL are threshold values used to decide if scaling should
 109 *>  be done based on the absolute size of the largest matrix element.
 110 *>  If AMAX > LARGE or AMAX < SMALL, scaling is done.
 111 *> \endverbatim
 112 *
 113 *  Authors:
 114 *  ========
 115 *
 116 *> \author Univ. of Tennessee
 117 *> \author Univ. of California Berkeley
 118 *> \author Univ. of Colorado Denver
 119 *> \author NAG Ltd.
 120 *
 121 *> \date September 2012
 122 *
 123 *> \ingroup realOTHERauxiliary
 124 *
 125 *  =====================================================================
 126       SUBROUTINE SLAQSP( UPLO, N, AP, S, SCOND, AMAX, EQUED )
 127 *
 128 *  -- LAPACK auxiliary routine (version 3.4.2) --
 129 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 130 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 131 *     September 2012
 132 *
 133 *     .. Scalar Arguments ..
 134       CHARACTER          EQUED, UPLO
 135       INTEGER            N
 136       REAL               AMAX, SCOND
 137 *     ..
 138 *     .. Array Arguments ..
 139       REAL               AP( * ), S( * )
 140 *     ..
 141 *
 142 *  =====================================================================
 143 *
 144 *     .. Parameters ..
 145       REAL               ONE, THRESH
 146       PARAMETER          ( ONE = 1.0E+0, THRESH = 0.1E+0 )
 147 *     ..
 148 *     .. Local Scalars ..
 149       INTEGER            I, J, JC
 150       REAL               CJ, LARGE, SMALL
 151 *     ..
 152 *     .. External Functions ..
 153       LOGICAL            LSAME
 154       REAL               SLAMCH
 155       EXTERNAL           LSAME, SLAMCH
 156 *     ..
 157 *     .. Executable Statements ..
 158 *
 159 *     Quick return if possible
 160 *
 161       IF( N.LE.0 ) THEN
 162          EQUED = 'N'
 163          RETURN
 164       END IF
 165 *
 166 *     Initialize LARGE and SMALL.
 167 *
 168       SMALL = SLAMCH( 'Safe minimum' ) / SLAMCH( 'Precision' )
 169       LARGE = ONE / SMALL
 170 *
 171       IF( SCOND.GE.THRESH .AND. AMAX.GE.SMALL .AND. AMAX.LE.LARGE ) THEN
 172 *
 173 *        No equilibration
 174 *
 175          EQUED = 'N'
 176       ELSE
 177 *
 178 *        Replace A by diag(S) * A * diag(S).
 179 *
 180          IF( LSAME( UPLO, 'U' ) ) THEN
 181 *
 182 *           Upper triangle of A is stored.
 183 *
 184             JC = 1
 185             DO 20 J = 1, N
 186                CJ = S( J )
 187                DO 10 I = 1, J
 188                   AP( JC+I-1 ) = CJ*S( I )*AP( JC+I-1 )
 189    10          CONTINUE
 190                JC = JC + J
 191    20       CONTINUE
 192          ELSE
 193 *
 194 *           Lower triangle of A is stored.
 195 *
 196             JC = 1
 197             DO 40 J = 1, N
 198                CJ = S( J )
 199                DO 30 I = J, N
 200                   AP( JC+I-J ) = CJ*S( I )*AP( JC+I-J )
 201    30          CONTINUE
 202                JC = JC + N - J + 1
 203    40       CONTINUE
 204          END IF
 205          EQUED = 'Y'
 206       END IF
 207 *
 208       RETURN
 209 *
 210 *     End of SLAQSP
 211 *
 212       END