LAPACK 3.3.1
Linear Algebra PACKage

zla_porcond_x.f

Go to the documentation of this file.
00001       DOUBLE PRECISION FUNCTION ZLA_PORCOND_X( UPLO, N, A, LDA, AF,
00002      $                                         LDAF, X, INFO, WORK,
00003      $                                         RWORK )
00004 *
00005 *     -- LAPACK routine (version 3.2.1)                                 --
00006 *     -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and --
00007 *     -- Jason Riedy of Univ. of California Berkeley.                 --
00008 *     -- April 2009                                                   --
00009 *
00010 *     -- LAPACK is a software package provided by Univ. of Tennessee, --
00011 *     -- Univ. of California Berkeley and NAG Ltd.                    --
00012 *
00013       IMPLICIT NONE
00014 *     ..
00015 *     .. Scalar Arguments ..
00016       CHARACTER          UPLO
00017       INTEGER            N, LDA, LDAF, INFO
00018 *     ..
00019 *     .. Array Arguments ..
00020       COMPLEX*16         A( LDA, * ), AF( LDAF, * ), WORK( * ), X( * )
00021       DOUBLE PRECISION   RWORK( * )
00022 *     ..
00023 *
00024 *  Purpose
00025 *  =======
00026 *
00027 *     ZLA_PORCOND_X Computes the infinity norm condition number of
00028 *     op(A) * diag(X) where X is a COMPLEX*16 vector.
00029 *
00030 *  Arguments
00031 *  =========
00032 *
00033 *     UPLO    (input) CHARACTER*1
00034 *       = 'U':  Upper triangle of A is stored;
00035 *       = 'L':  Lower triangle of A is stored.
00036 *
00037 *     N       (input) INTEGER
00038 *     The number of linear equations, i.e., the order of the
00039 *     matrix A.  N >= 0.
00040 *
00041 *     A       (input) COMPLEX*16 array, dimension (LDA,N)
00042 *     On entry, the N-by-N matrix A.
00043 *
00044 *     LDA     (input) INTEGER
00045 *     The leading dimension of the array A.  LDA >= max(1,N).
00046 *
00047 *     AF      (input) COMPLEX*16 array, dimension (LDAF,N)
00048 *     The triangular factor U or L from the Cholesky factorization
00049 *     A = U**H*U or A = L*L**H, as computed by ZPOTRF.
00050 *
00051 *     LDAF    (input) INTEGER
00052 *     The leading dimension of the array AF.  LDAF >= max(1,N).
00053 *
00054 *     X       (input) COMPLEX*16 array, dimension (N)
00055 *     The vector X in the formula op(A) * diag(X).
00056 *
00057 *     INFO    (output) INTEGER
00058 *       = 0:  Successful exit.
00059 *     i > 0:  The ith argument is invalid.
00060 *
00061 *     WORK    (input) COMPLEX*16 array, dimension (2*N).
00062 *     Workspace.
00063 *
00064 *     RWORK   (input) DOUBLE PRECISION array, dimension (N).
00065 *     Workspace.
00066 *
00067 *  =====================================================================
00068 *
00069 *     .. Local Scalars ..
00070       INTEGER            KASE, I, J
00071       DOUBLE PRECISION   AINVNM, ANORM, TMP
00072       LOGICAL            UP
00073       COMPLEX*16         ZDUM
00074 *     ..
00075 *     .. Local Arrays ..
00076       INTEGER            ISAVE( 3 )
00077 *     ..
00078 *     .. External Functions ..
00079       LOGICAL            LSAME
00080       EXTERNAL           LSAME
00081 *     ..
00082 *     .. External Subroutines ..
00083       EXTERNAL           ZLACN2, ZPOTRS, XERBLA
00084 *     ..
00085 *     .. Intrinsic Functions ..
00086       INTRINSIC          ABS, MAX, REAL, DIMAG
00087 *     ..
00088 *     .. Statement Functions ..
00089       DOUBLE PRECISION CABS1
00090 *     ..
00091 *     .. Statement Function Definitions ..
00092       CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
00093 *     ..
00094 *     .. Executable Statements ..
00095 *
00096       ZLA_PORCOND_X = 0.0D+0
00097 *
00098       INFO = 0
00099       IF( N.LT.0 ) THEN
00100          INFO = -2
00101       END IF
00102       IF( INFO.NE.0 ) THEN
00103          CALL XERBLA( 'ZLA_PORCOND_X', -INFO )
00104          RETURN
00105       END IF
00106       UP = .FALSE.
00107       IF ( LSAME( UPLO, 'U' ) ) UP = .TRUE.
00108 *
00109 *     Compute norm of op(A)*op2(C).
00110 *
00111       ANORM = 0.0D+0
00112       IF ( UP ) THEN
00113          DO I = 1, N
00114             TMP = 0.0D+0
00115             DO J = 1, I
00116                TMP = TMP + CABS1( A( J, I ) * X( J ) )
00117             END DO
00118             DO J = I+1, N
00119                TMP = TMP + CABS1( A( I, J ) * X( J ) )
00120             END DO
00121             RWORK( I ) = TMP
00122             ANORM = MAX( ANORM, TMP )
00123          END DO
00124       ELSE
00125          DO I = 1, N
00126             TMP = 0.0D+0
00127             DO J = 1, I
00128                TMP = TMP + CABS1( A( I, J ) * X( J ) )
00129             END DO
00130             DO J = I+1, N
00131                TMP = TMP + CABS1( A( J, I ) * X( J ) )
00132             END DO
00133             RWORK( I ) = TMP
00134             ANORM = MAX( ANORM, TMP )
00135          END DO
00136       END IF
00137 *
00138 *     Quick return if possible.
00139 *
00140       IF( N.EQ.0 ) THEN
00141          ZLA_PORCOND_X = 1.0D+0
00142          RETURN
00143       ELSE IF( ANORM .EQ. 0.0D+0 ) THEN
00144          RETURN
00145       END IF
00146 *
00147 *     Estimate the norm of inv(op(A)).
00148 *
00149       AINVNM = 0.0D+0
00150 *
00151       KASE = 0
00152    10 CONTINUE
00153       CALL ZLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
00154       IF( KASE.NE.0 ) THEN
00155          IF( KASE.EQ.2 ) THEN
00156 *
00157 *           Multiply by R.
00158 *
00159             DO I = 1, N
00160                WORK( I ) = WORK( I ) * RWORK( I )
00161             END DO
00162 *
00163             IF ( UP ) THEN
00164                CALL ZPOTRS( 'U', N, 1, AF, LDAF,
00165      $            WORK, N, INFO )
00166             ELSE
00167                CALL ZPOTRS( 'L', N, 1, AF, LDAF,
00168      $            WORK, N, INFO )
00169             ENDIF
00170 *
00171 *           Multiply by inv(X).
00172 *
00173             DO I = 1, N
00174                WORK( I ) = WORK( I ) / X( I )
00175             END DO
00176          ELSE
00177 *
00178 *           Multiply by inv(X**H).
00179 *
00180             DO I = 1, N
00181                WORK( I ) = WORK( I ) / X( I )
00182             END DO
00183 *
00184             IF ( UP ) THEN
00185                CALL ZPOTRS( 'U', N, 1, AF, LDAF,
00186      $            WORK, N, INFO )
00187             ELSE
00188                CALL ZPOTRS( 'L', N, 1, AF, LDAF,
00189      $            WORK, N, INFO )
00190             END IF
00191 *
00192 *           Multiply by R.
00193 *
00194             DO I = 1, N
00195                WORK( I ) = WORK( I ) * RWORK( I )
00196             END DO
00197          END IF
00198          GO TO 10
00199       END IF
00200 *
00201 *     Compute the estimate of the reciprocal condition number.
00202 *
00203       IF( AINVNM .NE. 0.0D+0 )
00204      $   ZLA_PORCOND_X = 1.0D+0 / AINVNM
00205 *
00206       RETURN
00207 *
00208       END
 All Files Functions