1       DOUBLE PRECISION FUNCTION ZLA_SYRCOND_X( UPLO, N, A, LDA, AF,
  2      $                                         LDAF, IPIV, X, INFO,
  3      $                                         WORK, RWORK )
  4 *
  5 *     -- LAPACK routine (version 3.2.1)                                 --
  6 *     -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and --
  7 *     -- Jason Riedy of Univ. of California Berkeley.                 --
  8 *     -- April 2009                                                   --
  9 *
 10 *     -- LAPACK is a software package provided by Univ. of Tennessee, --
 11 *     -- Univ. of California Berkeley and NAG Ltd.                    --
 12 *
 13       IMPLICIT NONE
 14 *     ..
 15 *     .. Scalar Arguments ..
 16       CHARACTER          UPLO
 17       INTEGER            N, LDA, LDAF, INFO
 18 *     ..
 19 *     .. Array Arguments ..
 20       INTEGER            IPIV( * )
 21       COMPLEX*16         A( LDA, * ), AF( LDAF, * ), WORK( * ), X( * )
 22       DOUBLE PRECISION   RWORK( * )
 23 *     ..
 24 *
 25 *  Purpose
 26 *  =======
 27 *
 28 *     ZLA_SYRCOND_X Computes the infinity norm condition number of
 29 *     op(A) * diag(X) where X is a COMPLEX*16 vector.
 30 *
 31 *  Arguments
 32 *  =========
 33 *
 34 *     UPLO    (input) CHARACTER*1
 35 *       = 'U':  Upper triangle of A is stored;
 36 *       = 'L':  Lower triangle of A is stored.
 37 *
 38 *     N       (input) INTEGER
 39 *     The number of linear equations, i.e., the order of the
 40 *     matrix A.  N >= 0.
 41 *
 42 *     A       (input) COMPLEX*16 array, dimension (LDA,N)
 43 *     On entry, the N-by-N matrix A.
 44 *
 45 *     LDA     (input) INTEGER
 46 *     The leading dimension of the array A.  LDA >= max(1,N).
 47 *
 48 *     AF      (input) COMPLEX*16 array, dimension (LDAF,N)
 49 *     The block diagonal matrix D and the multipliers used to
 50 *     obtain the factor U or L as computed by ZSYTRF.
 51 *
 52 *     LDAF    (input) INTEGER
 53 *     The leading dimension of the array AF.  LDAF >= max(1,N).
 54 *
 55 *     IPIV    (input) INTEGER array, dimension (N)
 56 *     Details of the interchanges and the block structure of D
 57 *     as determined by ZSYTRF.
 58 *
 59 *     X       (input) COMPLEX*16 array, dimension (N)
 60 *     The vector X in the formula op(A) * diag(X).
 61 *
 62 *     INFO    (output) INTEGER
 63 *       = 0:  Successful exit.
 64 *     i > 0:  The ith argument is invalid.
 65 *
 66 *     WORK    (input) COMPLEX*16 array, dimension (2*N).
 67 *     Workspace.
 68 *
 69 *     RWORK   (input) DOUBLE PRECISION array, dimension (N).
 70 *     Workspace.
 71 *
 72 *  =====================================================================
 73 *
 74 *     .. Local Scalars ..
 75       INTEGER            KASE
 76       DOUBLE PRECISION   AINVNM, ANORM, TMP
 77       INTEGER            I, J
 78       LOGICAL            UP
 79       COMPLEX*16         ZDUM
 80 *     ..
 81 *     .. Local Arrays ..
 82       INTEGER            ISAVE( 3 )
 83 *     ..
 84 *     .. External Functions ..
 85       LOGICAL            LSAME
 86       EXTERNAL           LSAME
 87 *     ..
 88 *     .. External Subroutines ..
 89       EXTERNAL           ZLACN2, ZSYTRS, XERBLA
 90 *     ..
 91 *     .. Intrinsic Functions ..
 92       INTRINSIC          ABSMAX
 93 *     ..
 94 *     .. Statement Functions ..
 95       DOUBLE PRECISION   CABS1
 96 *     ..
 97 *     .. Statement Function Definitions ..
 98       CABS1( ZDUM ) = ABSDBLE( ZDUM ) ) + ABSDIMAG( ZDUM ) )
 99 *     ..
100 *     .. Executable Statements ..
101 *
102       ZLA_SYRCOND_X = 0.0D+0
103 *
104       INFO = 0
105       IF( N.LT.0 ) THEN
106          INFO = -2
107       END IF
108       IF( INFO.NE.0 ) THEN
109          CALL XERBLA( 'ZLA_SYRCOND_X'-INFO )
110          RETURN
111       END IF
112       UP = .FALSE.
113       IF ( LSAME( UPLO, 'U' ) ) UP = .TRUE.
114 *
115 *     Compute norm of op(A)*op2(C).
116 *
117       ANORM = 0.0D+0
118       IF ( UP ) THEN
119          DO I = 1, N
120             TMP = 0.0D+0
121             DO J = 1, I
122                TMP = TMP + CABS1( A( J, I ) * X( J ) )
123             END DO
124             DO J = I+1, N
125                TMP = TMP + CABS1( A( I, J ) * X( J ) )
126             END DO
127             RWORK( I ) = TMP
128             ANORM = MAX( ANORM, TMP )
129          END DO
130       ELSE
131          DO I = 1, N
132             TMP = 0.0D+0
133             DO J = 1, I
134                TMP = TMP + CABS1( A( I, J ) * X( J ) )
135             END DO
136             DO J = I+1, N
137                TMP = TMP + CABS1( A( J, I ) * X( J ) )
138             END DO
139             RWORK( I ) = TMP
140             ANORM = MAX( ANORM, TMP )
141          END DO
142       END IF
143 *
144 *     Quick return if possible.
145 *
146       IF( N.EQ.0 ) THEN
147          ZLA_SYRCOND_X = 1.0D+0
148          RETURN
149       ELSE IF( ANORM .EQ. 0.0D+0 ) THEN
150          RETURN
151       END IF
152 *
153 *     Estimate the norm of inv(op(A)).
154 *
155       AINVNM = 0.0D+0
156 *
157       KASE = 0
158    10 CONTINUE
159       CALL ZLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
160       IF( KASE.NE.0 ) THEN
161          IF( KASE.EQ.2 ) THEN
162 *
163 *           Multiply by R.
164 *
165             DO I = 1, N
166                WORK( I ) = WORK( I ) * RWORK( I )
167             END DO
168 *
169             IF ( UP ) THEN
170                CALL ZSYTRS( 'U', N, 1, AF, LDAF, IPIV,
171      $            WORK, N, INFO )
172             ELSE
173                CALL ZSYTRS( 'L', N, 1, AF, LDAF, IPIV,
174      $            WORK, N, INFO )
175             ENDIF
176 *
177 *           Multiply by inv(X).
178 *
179             DO I = 1, N
180                WORK( I ) = WORK( I ) / X( I )
181             END DO
182          ELSE
183 *
184 *           Multiply by inv(X**T).
185 *
186             DO I = 1, N
187                WORK( I ) = WORK( I ) / X( I )
188             END DO
189 *
190             IF ( UP ) THEN
191                CALL ZSYTRS( 'U', N, 1, AF, LDAF, IPIV,
192      $            WORK, N, INFO )
193             ELSE
194                CALL ZSYTRS( 'L', N, 1, AF, LDAF, IPIV,
195      $            WORK, N, INFO )
196             END IF
197 *
198 *           Multiply by R.
199 *
200             DO I = 1, N
201                WORK( I ) = WORK( I ) * RWORK( I )
202             END DO
203          END IF
204          GO TO 10
205       END IF
206 *
207 *     Compute the estimate of the reciprocal condition number.
208 *
209       IF( AINVNM .NE. 0.0D+0 )
210      $   ZLA_SYRCOND_X = 1.0D+0 / AINVNM
211 *
212       RETURN
213 *
214       END