Annotation of rpl/lapack/lapack/zsycon.f, revision 1.8

1.1       bertrand    1:       SUBROUTINE ZSYCON( UPLO, N, A, LDA, IPIV, ANORM, RCOND, WORK,
                      2:      $                   INFO )
                      3: *
1.8     ! bertrand    4: *  -- LAPACK routine (version 3.3.1) --
1.1       bertrand    5: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
                      6: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.8     ! bertrand    7: *  -- April 2011                                                      --
1.1       bertrand    8: *
                      9: *     Modified to call ZLACN2 in place of ZLACON, 10 Feb 03, SJH.
                     10: *
                     11: *     .. Scalar Arguments ..
                     12:       CHARACTER          UPLO
                     13:       INTEGER            INFO, LDA, N
                     14:       DOUBLE PRECISION   ANORM, RCOND
                     15: *     ..
                     16: *     .. Array Arguments ..
                     17:       INTEGER            IPIV( * )
                     18:       COMPLEX*16         A( LDA, * ), WORK( * )
                     19: *     ..
                     20: *
                     21: *  Purpose
                     22: *  =======
                     23: *
                     24: *  ZSYCON estimates the reciprocal of the condition number (in the
                     25: *  1-norm) of a complex symmetric matrix A using the factorization
                     26: *  A = U*D*U**T or A = L*D*L**T computed by ZSYTRF.
                     27: *
                     28: *  An estimate is obtained for norm(inv(A)), and the reciprocal of the
                     29: *  condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))).
                     30: *
                     31: *  Arguments
                     32: *  =========
                     33: *
                     34: *  UPLO    (input) CHARACTER*1
                     35: *          Specifies whether the details of the factorization are stored
                     36: *          as an upper or lower triangular matrix.
                     37: *          = 'U':  Upper triangular, form is A = U*D*U**T;
                     38: *          = 'L':  Lower triangular, form is A = L*D*L**T.
                     39: *
                     40: *  N       (input) INTEGER
                     41: *          The order of the matrix A.  N >= 0.
                     42: *
                     43: *  A       (input) COMPLEX*16 array, dimension (LDA,N)
                     44: *          The block diagonal matrix D and the multipliers used to
                     45: *          obtain the factor U or L as computed by ZSYTRF.
                     46: *
                     47: *  LDA     (input) INTEGER
                     48: *          The leading dimension of the array A.  LDA >= max(1,N).
                     49: *
                     50: *  IPIV    (input) INTEGER array, dimension (N)
                     51: *          Details of the interchanges and the block structure of D
                     52: *          as determined by ZSYTRF.
                     53: *
                     54: *  ANORM   (input) DOUBLE PRECISION
                     55: *          The 1-norm of the original matrix A.
                     56: *
                     57: *  RCOND   (output) DOUBLE PRECISION
                     58: *          The reciprocal of the condition number of the matrix A,
                     59: *          computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is an
                     60: *          estimate of the 1-norm of inv(A) computed in this routine.
                     61: *
                     62: *  WORK    (workspace) COMPLEX*16 array, dimension (2*N)
                     63: *
                     64: *  INFO    (output) INTEGER
                     65: *          = 0:  successful exit
                     66: *          < 0:  if INFO = -i, the i-th argument had an illegal value
                     67: *
                     68: *  =====================================================================
                     69: *
                     70: *     .. Parameters ..
                     71:       DOUBLE PRECISION   ONE, ZERO
                     72:       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
                     73: *     ..
                     74: *     .. Local Scalars ..
                     75:       LOGICAL            UPPER
                     76:       INTEGER            I, KASE
                     77:       DOUBLE PRECISION   AINVNM
                     78: *     ..
                     79: *     .. Local Arrays ..
                     80:       INTEGER            ISAVE( 3 )
                     81: *     ..
                     82: *     .. External Functions ..
                     83:       LOGICAL            LSAME
                     84:       EXTERNAL           LSAME
                     85: *     ..
                     86: *     .. External Subroutines ..
                     87:       EXTERNAL           XERBLA, ZLACN2, ZSYTRS
                     88: *     ..
                     89: *     .. Intrinsic Functions ..
                     90:       INTRINSIC          MAX
                     91: *     ..
                     92: *     .. Executable Statements ..
                     93: *
                     94: *     Test the input parameters.
                     95: *
                     96:       INFO = 0
                     97:       UPPER = LSAME( UPLO, 'U' )
                     98:       IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
                     99:          INFO = -1
                    100:       ELSE IF( N.LT.0 ) THEN
                    101:          INFO = -2
                    102:       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
                    103:          INFO = -4
                    104:       ELSE IF( ANORM.LT.ZERO ) THEN
                    105:          INFO = -6
                    106:       END IF
                    107:       IF( INFO.NE.0 ) THEN
                    108:          CALL XERBLA( 'ZSYCON', -INFO )
                    109:          RETURN
                    110:       END IF
                    111: *
                    112: *     Quick return if possible
                    113: *
                    114:       RCOND = ZERO
                    115:       IF( N.EQ.0 ) THEN
                    116:          RCOND = ONE
                    117:          RETURN
                    118:       ELSE IF( ANORM.LE.ZERO ) THEN
                    119:          RETURN
                    120:       END IF
                    121: *
                    122: *     Check that the diagonal matrix D is nonsingular.
                    123: *
                    124:       IF( UPPER ) THEN
                    125: *
                    126: *        Upper triangular storage: examine D from bottom to top
                    127: *
                    128:          DO 10 I = N, 1, -1
                    129:             IF( IPIV( I ).GT.0 .AND. A( I, I ).EQ.ZERO )
                    130:      $         RETURN
                    131:    10    CONTINUE
                    132:       ELSE
                    133: *
                    134: *        Lower triangular storage: examine D from top to bottom.
                    135: *
                    136:          DO 20 I = 1, N
                    137:             IF( IPIV( I ).GT.0 .AND. A( I, I ).EQ.ZERO )
                    138:      $         RETURN
                    139:    20    CONTINUE
                    140:       END IF
                    141: *
                    142: *     Estimate the 1-norm of the inverse.
                    143: *
                    144:       KASE = 0
                    145:    30 CONTINUE
                    146:       CALL ZLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
                    147:       IF( KASE.NE.0 ) THEN
                    148: *
1.8     ! bertrand  149: *        Multiply by inv(L*D*L**T) or inv(U*D*U**T).
1.1       bertrand  150: *
                    151:          CALL ZSYTRS( UPLO, N, 1, A, LDA, IPIV, WORK, N, INFO )
                    152:          GO TO 30
                    153:       END IF
                    154: *
                    155: *     Compute the estimate of the reciprocal condition number.
                    156: *
                    157:       IF( AINVNM.NE.ZERO )
                    158:      $   RCOND = ( ONE / AINVNM ) / ANORM
                    159: *
                    160:       RETURN
                    161: *
                    162: *     End of ZSYCON
                    163: *
                    164:       END

CVSweb interface <joel.bertrand@systella.fr>