SRC/cla_syrcond_x.f

   1 *> \brief \b CLA_SYRCOND_X computes the infinity norm condition number of op(A)*diag(x) for symmetric indefinite matrices.
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download CLA_SYRCOND_X + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cla_syrcond_x.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cla_syrcond_x.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cla_syrcond_x.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       REAL FUNCTION CLA_SYRCOND_X( UPLO, N, A, LDA, AF, LDAF, IPIV, X,
  22 *                                    INFO, WORK, RWORK )
  23 *
  24 *       .. Scalar Arguments ..
  25 *       CHARACTER          UPLO
  26 *       INTEGER            N, LDA, LDAF, INFO
  27 *       ..
  28 *       .. Array Arguments ..
  29 *       INTEGER            IPIV( * )
  30 *       COMPLEX            A( LDA, * ), AF( LDAF, * ), WORK( * ), X( * )
  31 *       REAL               RWORK( * )
  32 *       ..
  33 *
  34 *
  35 *> \par Purpose:
  36 *  =============
  37 *>
  38 *> \verbatim
  39 *>
  40 *>    CLA_SYRCOND_X Computes the infinity norm condition number of
  41 *>    op(A) * diag(X) where X is a COMPLEX vector.
  42 *> \endverbatim
  43 *
  44 *  Arguments:
  45 *  ==========
  46 *
  47 *> \param[in] UPLO
  48 *> \verbatim
  49 *>          UPLO is CHARACTER*1
  50 *>       = 'U':  Upper triangle of A is stored;
  51 *>       = 'L':  Lower triangle of A is stored.
  52 *> \endverbatim
  53 *>
  54 *> \param[in] N
  55 *> \verbatim
  56 *>          N is INTEGER
  57 *>     The number of linear equations, i.e., the order of the
  58 *>     matrix A.  N >= 0.
  59 *> \endverbatim
  60 *>
  61 *> \param[in] A
  62 *> \verbatim
  63 *>          A is COMPLEX array, dimension (LDA,N)
  64 *>     On entry, the N-by-N matrix A.
  65 *> \endverbatim
  66 *>
  67 *> \param[in] LDA
  68 *> \verbatim
  69 *>          LDA is INTEGER
  70 *>     The leading dimension of the array A.  LDA >= max(1,N).
  71 *> \endverbatim
  72 *>
  73 *> \param[in] AF
  74 *> \verbatim
  75 *>          AF is COMPLEX array, dimension (LDAF,N)
  76 *>     The block diagonal matrix D and the multipliers used to
  77 *>     obtain the factor U or L as computed by CSYTRF.
  78 *> \endverbatim
  79 *>
  80 *> \param[in] LDAF
  81 *> \verbatim
  82 *>          LDAF is INTEGER
  83 *>     The leading dimension of the array AF.  LDAF >= max(1,N).
  84 *> \endverbatim
  85 *>
  86 *> \param[in] IPIV
  87 *> \verbatim
  88 *>          IPIV is INTEGER array, dimension (N)
  89 *>     Details of the interchanges and the block structure of D
  90 *>     as determined by CSYTRF.
  91 *> \endverbatim
  92 *>
  93 *> \param[in] X
  94 *> \verbatim
  95 *>          X is COMPLEX array, dimension (N)
  96 *>     The vector X in the formula op(A) * diag(X).
  97 *> \endverbatim
  98 *>
  99 *> \param[out] INFO
 100 *> \verbatim
 101 *>          INFO is INTEGER
 102 *>       = 0:  Successful exit.
 103 *>     i > 0:  The ith argument is invalid.
 104 *> \endverbatim
 105 *>
 106 *> \param[in] WORK
 107 *> \verbatim
 108 *>          WORK is COMPLEX array, dimension (2*N).
 109 *>     Workspace.
 110 *> \endverbatim
 111 *>
 112 *> \param[in] RWORK
 113 *> \verbatim
 114 *>          RWORK is REAL array, dimension (N).
 115 *>     Workspace.
 116 *> \endverbatim
 117 *
 118 *  Authors:
 119 *  ========
 120 *
 121 *> \author Univ. of Tennessee
 122 *> \author Univ. of California Berkeley
 123 *> \author Univ. of Colorado Denver
 124 *> \author NAG Ltd.
 125 *
 126 *> \date September 2012
 127 *
 128 *> \ingroup complexSYcomputational
 129 *
 130 *  =====================================================================
 131       REAL FUNCTION CLA_SYRCOND_X( UPLO, N, A, LDA, AF, LDAF, IPIV, X,
 132      $                             INFO, WORK, RWORK )
 133 *
 134 *  -- LAPACK computational routine (version 3.4.2) --
 135 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 136 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 137 *     September 2012
 138 *
 139 *     .. Scalar Arguments ..
 140       CHARACTER          UPLO
 141       INTEGER            N, LDA, LDAF, INFO
 142 *     ..
 143 *     .. Array Arguments ..
 144       INTEGER            IPIV( * )
 145       COMPLEX            A( LDA, * ), AF( LDAF, * ), WORK( * ), X( * )
 146       REAL               RWORK( * )
 147 *     ..
 148 *
 149 *  =====================================================================
 150 *
 151 *     .. Local Scalars ..
 152       INTEGER            KASE
 153       REAL               AINVNM, ANORM, TMP
 154       INTEGER            I, J
 155       LOGICAL            UP, UPPER
 156       COMPLEX            ZDUM
 157 *     ..
 158 *     .. Local Arrays ..
 159       INTEGER            ISAVE( 3 )
 160 *     ..
 161 *     .. External Functions ..
 162       LOGICAL            LSAME
 163       EXTERNAL           LSAME
 164 *     ..
 165 *     .. External Subroutines ..
 166       EXTERNAL           CLACN2, CSYTRS, XERBLA
 167 *     ..
 168 *     .. Intrinsic Functions ..
 169       INTRINSIC          ABS, MAX
 170 *     ..
 171 *     .. Statement Functions ..
 172       REAL               CABS1
 173 *     ..
 174 *     .. Statement Function Definitions ..
 175       CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) )
 176 *     ..
 177 *     .. Executable Statements ..
 178 *
 179       CLA_SYRCOND_X = 0.0E+0
 180 *
 181       INFO = 0
 182       UPPER = LSAME( UPLO, 'U' )
 183       IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
 184          INFO = -1
 185       ELSE IF ( N.LT.0 ) THEN
 186          INFO = -2
 187       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
 188          INFO = -4
 189       ELSE IF( LDAF.LT.MAX( 1, N ) ) THEN
 190          INFO = -6
 191       END IF
 192       IF( INFO.NE.0 ) THEN
 193          CALL XERBLA( 'CLA_SYRCOND_X', -INFO )
 194          RETURN
 195       END IF
 196       UP = .FALSE.
 197       IF ( LSAME( UPLO, 'U' ) ) UP = .TRUE.
 198 *
 199 *     Compute norm of op(A)*op2(C).
 200 *
 201       ANORM = 0.0
 202       IF ( UP ) THEN
 203          DO I = 1, N
 204             TMP = 0.0E+0
 205             DO J = 1, I
 206                TMP = TMP + CABS1( A( J, I ) * X( J ) )
 207             END DO
 208             DO J = I+1, N
 209                TMP = TMP + CABS1( A( I, J ) * X( J ) )
 210             END DO
 211             RWORK( I ) = TMP
 212             ANORM = MAX( ANORM, TMP )
 213          END DO
 214       ELSE
 215          DO I = 1, N
 216             TMP = 0.0E+0
 217             DO J = 1, I
 218                TMP = TMP + CABS1( A( I, J ) * X( J ) )
 219             END DO
 220             DO J = I+1, N
 221                TMP = TMP + CABS1( A( J, I ) * X( J ) )
 222             END DO
 223             RWORK( I ) = TMP
 224             ANORM = MAX( ANORM, TMP )
 225          END DO
 226       END IF
 227 *
 228 *     Quick return if possible.
 229 *
 230       IF( N.EQ.0 ) THEN
 231          CLA_SYRCOND_X = 1.0E+0
 232          RETURN
 233       ELSE IF( ANORM .EQ. 0.0E+0 ) THEN
 234          RETURN
 235       END IF
 236 *
 237 *     Estimate the norm of inv(op(A)).
 238 *
 239       AINVNM = 0.0E+0
 240 *
 241       KASE = 0
 242    10 CONTINUE
 243       CALL CLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
 244       IF( KASE.NE.0 ) THEN
 245          IF( KASE.EQ.2 ) THEN
 246 *
 247 *           Multiply by R.
 248 *
 249             DO I = 1, N
 250                WORK( I ) = WORK( I ) * RWORK( I )
 251             END DO
 252 *
 253             IF ( UP ) THEN
 254                CALL CSYTRS( 'U', N, 1, AF, LDAF, IPIV,
 255      $            WORK, N, INFO )
 256             ELSE
 257                CALL CSYTRS( 'L', N, 1, AF, LDAF, IPIV,
 258      $            WORK, N, INFO )
 259             ENDIF
 260 *
 261 *           Multiply by inv(X).
 262 *
 263             DO I = 1, N
 264                WORK( I ) = WORK( I ) / X( I )
 265             END DO
 266          ELSE
 267 *
 268 *           Multiply by inv(X**T).
 269 *
 270             DO I = 1, N
 271                WORK( I ) = WORK( I ) / X( I )
 272             END DO
 273 *
 274             IF ( UP ) THEN
 275                CALL CSYTRS( 'U', N, 1, AF, LDAF, IPIV,
 276      $            WORK, N, INFO )
 277             ELSE
 278                CALL CSYTRS( 'L', N, 1, AF, LDAF, IPIV,
 279      $            WORK, N, INFO )
 280             END IF
 281 *
 282 *           Multiply by R.
 283 *
 284             DO I = 1, N
 285                WORK( I ) = WORK( I ) * RWORK( I )
 286             END DO
 287          END IF
 288          GO TO 10
 289       END IF
 290 *
 291 *     Compute the estimate of the reciprocal condition number.
 292 *
 293       IF( AINVNM .NE. 0.0E+0 )
 294      $   CLA_SYRCOND_X = 1.0E+0 / AINVNM
 295 *
 296       RETURN
 297 *
 298       END