1       SUBROUTINE DPTTRF( N, D, E, INFO )
  2 *
  3 *  -- LAPACK routine (version 3.3.1) --
  4 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  5 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  6 *  -- April 2011                                                      --
  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**T factorization of a real symmetric
 19 *  positive definite tridiagonal matrix A.  The factorization may also
 20 *  be regarded as having the form A = U**T*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**T 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**T factorization of A.
 37 *          E can also be regarded as the superdiagonal of the unit
 38 *          bidiagonal factor U from the U**T*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**T (or U**T*D*U) factorization of A.
 81 *
 82       I4 = MOD( N-14 )
 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 - 44
 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