Return to dpttrf.f CVS log

Up to [local] / rpl / lapack / lapack

Annotation of rpl/lapack/lapack/dpttrf.f, revision 1.1

1.1     ! bertrand    1:       SUBROUTINE DPTTRF( N, D, E, INFO )
        !             2: *
        !             3: *  -- LAPACK routine (version 3.2) --
        !             4: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
        !             5: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
        !             6: *     November 2006
        !             7: *
        !             8: *     .. Scalar Arguments ..
        !             9:       INTEGER            INFO, N
        !            10: *     ..
        !            11: *     .. Array Arguments ..
        !            12:       DOUBLE PRECISION   D( * ), E( * )
        !            13: *     ..
        !            14: *
        !            15: *  Purpose
        !            16: *  =======
        !            17: *
        !            18: *  DPTTRF computes the L*D*L' factorization of a real symmetric
        !            19: *  positive definite tridiagonal matrix A.  The factorization may also
        !            20: *  be regarded as having the form A = U'*D*U.
        !            21: *
        !            22: *  Arguments
        !            23: *  =========
        !            24: *
        !            25: *  N       (input) INTEGER
        !            26: *          The order of the matrix A.  N >= 0.
        !            27: *
        !            28: *  D       (input/output) DOUBLE PRECISION array, dimension (N)
        !            29: *          On entry, the n diagonal elements of the tridiagonal matrix
        !            30: *          A.  On exit, the n diagonal elements of the diagonal matrix
        !            31: *          D from the L*D*L' factorization of A.
        !            32: *
        !            33: *  E       (input/output) DOUBLE PRECISION array, dimension (N-1)
        !            34: *          On entry, the (n-1) subdiagonal elements of the tridiagonal
        !            35: *          matrix A.  On exit, the (n-1) subdiagonal elements of the
        !            36: *          unit bidiagonal factor L from the L*D*L' factorization of A.
        !            37: *          E can also be regarded as the superdiagonal of the unit
        !            38: *          bidiagonal factor U from the U'*D*U factorization of A.
        !            39: *
        !            40: *  INFO    (output) INTEGER
        !            41: *          = 0: successful exit
        !            42: *          < 0: if INFO = -k, the k-th argument had an illegal value
        !            43: *          > 0: if INFO = k, the leading minor of order k is not
        !            44: *               positive definite; if k < N, the factorization could not
        !            45: *               be completed, while if k = N, the factorization was
        !            46: *               completed, but D(N) <= 0.
        !            47: *
        !            48: *  =====================================================================
        !            49: *
        !            50: *     .. Parameters ..
        !            51:       DOUBLE PRECISION   ZERO
        !            52:       PARAMETER          ( ZERO = 0.0D+0 )
        !            53: *     ..
        !            54: *     .. Local Scalars ..
        !            55:       INTEGER            I, I4
        !            56:       DOUBLE PRECISION   EI
        !            57: *     ..
        !            58: *     .. External Subroutines ..
        !            59:       EXTERNAL           XERBLA
        !            60: *     ..
        !            61: *     .. Intrinsic Functions ..
        !            62:       INTRINSIC          MOD
        !            63: *     ..
        !            64: *     .. Executable Statements ..
        !            65: *
        !            66: *     Test the input parameters.
        !            67: *
        !            68:       INFO = 0
        !            69:       IF( N.LT.0 ) THEN
        !            70:          INFO = -1
        !            71:          CALL XERBLA( 'DPTTRF', -INFO )
        !            72:          RETURN
        !            73:       END IF
        !            74: *
        !            75: *     Quick return if possible
        !            76: *
        !            77:       IF( N.EQ.0 )
        !            78:      $   RETURN
        !            79: *
        !            80: *     Compute the L*D*L' (or U'*D*U) factorization of A.
        !            81: *
        !            82:       I4 = MOD( N-1, 4 )
        !            83:       DO 10 I = 1, I4
        !            84:          IF( D( I ).LE.ZERO ) THEN
        !            85:             INFO = I
        !            86:             GO TO 30
        !            87:          END IF
        !            88:          EI = E( I )
        !            89:          E( I ) = EI / D( I )
        !            90:          D( I+1 ) = D( I+1 ) - E( I )*EI
        !            91:    10 CONTINUE
        !            92: *
        !            93:       DO 20 I = I4 + 1, N - 4, 4
        !            94: *
        !            95: *        Drop out of the loop if d(i) <= 0: the matrix is not positive
        !            96: *        definite.
        !            97: *
        !            98:          IF( D( I ).LE.ZERO ) THEN
        !            99:             INFO = I
        !           100:             GO TO 30
        !           101:          END IF
        !           102: *
        !           103: *        Solve for e(i) and d(i+1).
        !           104: *
        !           105:          EI = E( I )
        !           106:          E( I ) = EI / D( I )
        !           107:          D( I+1 ) = D( I+1 ) - E( I )*EI
        !           108: *
        !           109:          IF( D( I+1 ).LE.ZERO ) THEN
        !           110:             INFO = I + 1
        !           111:             GO TO 30
        !           112:          END IF
        !           113: *
        !           114: *        Solve for e(i+1) and d(i+2).
        !           115: *
        !           116:          EI = E( I+1 )
        !           117:          E( I+1 ) = EI / D( I+1 )
        !           118:          D( I+2 ) = D( I+2 ) - E( I+1 )*EI
        !           119: *
        !           120:          IF( D( I+2 ).LE.ZERO ) THEN
        !           121:             INFO = I + 2
        !           122:             GO TO 30
        !           123:          END IF
        !           124: *
        !           125: *        Solve for e(i+2) and d(i+3).
        !           126: *
        !           127:          EI = E( I+2 )
        !           128:          E( I+2 ) = EI / D( I+2 )
        !           129:          D( I+3 ) = D( I+3 ) - E( I+2 )*EI
        !           130: *
        !           131:          IF( D( I+3 ).LE.ZERO ) THEN
        !           132:             INFO = I + 3
        !           133:             GO TO 30
        !           134:          END IF
        !           135: *
        !           136: *        Solve for e(i+3) and d(i+4).
        !           137: *
        !           138:          EI = E( I+3 )
        !           139:          E( I+3 ) = EI / D( I+3 )
        !           140:          D( I+4 ) = D( I+4 ) - E( I+3 )*EI
        !           141:    20 CONTINUE
        !           142: *
        !           143: *     Check d(n) for positive definiteness.
        !           144: *
        !           145:       IF( D( N ).LE.ZERO )
        !           146:      $   INFO = N
        !           147: *
        !           148:    30 CONTINUE
        !           149:       RETURN
        !           150: *
        !           151: *     End of DPTTRF
        !           152: *
        !           153:       END