1: *> \brief \b DPTTRF
2: *
3: * =========== DOCUMENTATION ===========
4: *
5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
7: *
8: *> \htmlonly
9: *> Download DPTTRF + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dpttrf.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dpttrf.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dpttrf.f">
15: *> [TXT]</a>
16: *> \endhtmlonly
17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE DPTTRF( N, D, E, INFO )
22: *
23: * .. Scalar Arguments ..
24: * INTEGER INFO, N
25: * ..
26: * .. Array Arguments ..
27: * DOUBLE PRECISION D( * ), E( * )
28: * ..
29: *
30: *
31: *> \par Purpose:
32: * =============
33: *>
34: *> \verbatim
35: *>
36: *> DPTTRF computes the L*D*L**T factorization of a real symmetric
37: *> positive definite tridiagonal matrix A. The factorization may also
38: *> be regarded as having the form A = U**T*D*U.
39: *> \endverbatim
40: *
41: * Arguments:
42: * ==========
43: *
44: *> \param[in] N
45: *> \verbatim
46: *> N is INTEGER
47: *> The order of the matrix A. N >= 0.
48: *> \endverbatim
49: *>
50: *> \param[in,out] D
51: *> \verbatim
52: *> D is DOUBLE PRECISION array, dimension (N)
53: *> On entry, the n diagonal elements of the tridiagonal matrix
54: *> A. On exit, the n diagonal elements of the diagonal matrix
55: *> D from the L*D*L**T factorization of A.
56: *> \endverbatim
57: *>
58: *> \param[in,out] E
59: *> \verbatim
60: *> E is DOUBLE PRECISION array, dimension (N-1)
61: *> On entry, the (n-1) subdiagonal elements of the tridiagonal
62: *> matrix A. On exit, the (n-1) subdiagonal elements of the
63: *> unit bidiagonal factor L from the L*D*L**T factorization of A.
64: *> E can also be regarded as the superdiagonal of the unit
65: *> bidiagonal factor U from the U**T*D*U factorization of A.
66: *> \endverbatim
67: *>
68: *> \param[out] INFO
69: *> \verbatim
70: *> INFO is INTEGER
71: *> = 0: successful exit
72: *> < 0: if INFO = -k, the k-th argument had an illegal value
73: *> > 0: if INFO = k, the leading minor of order k is not
74: *> positive definite; if k < N, the factorization could not
75: *> be completed, while if k = N, the factorization was
76: *> completed, but D(N) <= 0.
77: *> \endverbatim
78: *
79: * Authors:
80: * ========
81: *
82: *> \author Univ. of Tennessee
83: *> \author Univ. of California Berkeley
84: *> \author Univ. of Colorado Denver
85: *> \author NAG Ltd.
86: *
87: *> \ingroup doublePTcomputational
88: *
89: * =====================================================================
90: SUBROUTINE DPTTRF( N, D, E, INFO )
91: *
92: * -- LAPACK computational routine --
93: * -- LAPACK is a software package provided by Univ. of Tennessee, --
94: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
95: *
96: * .. Scalar Arguments ..
97: INTEGER INFO, N
98: * ..
99: * .. Array Arguments ..
100: DOUBLE PRECISION D( * ), E( * )
101: * ..
102: *
103: * =====================================================================
104: *
105: * .. Parameters ..
106: DOUBLE PRECISION ZERO
107: PARAMETER ( ZERO = 0.0D+0 )
108: * ..
109: * .. Local Scalars ..
110: INTEGER I, I4
111: DOUBLE PRECISION EI
112: * ..
113: * .. External Subroutines ..
114: EXTERNAL XERBLA
115: * ..
116: * .. Intrinsic Functions ..
117: INTRINSIC MOD
118: * ..
119: * .. Executable Statements ..
120: *
121: * Test the input parameters.
122: *
123: INFO = 0
124: IF( N.LT.0 ) THEN
125: INFO = -1
126: CALL XERBLA( 'DPTTRF', -INFO )
127: RETURN
128: END IF
129: *
130: * Quick return if possible
131: *
132: IF( N.EQ.0 )
133: $ RETURN
134: *
135: * Compute the L*D*L**T (or U**T*D*U) factorization of A.
136: *
137: I4 = MOD( N-1, 4 )
138: DO 10 I = 1, I4
139: IF( D( I ).LE.ZERO ) THEN
140: INFO = I
141: GO TO 30
142: END IF
143: EI = E( I )
144: E( I ) = EI / D( I )
145: D( I+1 ) = D( I+1 ) - E( I )*EI
146: 10 CONTINUE
147: *
148: DO 20 I = I4 + 1, N - 4, 4
149: *
150: * Drop out of the loop if d(i) <= 0: the matrix is not positive
151: * definite.
152: *
153: IF( D( I ).LE.ZERO ) THEN
154: INFO = I
155: GO TO 30
156: END IF
157: *
158: * Solve for e(i) and d(i+1).
159: *
160: EI = E( I )
161: E( I ) = EI / D( I )
162: D( I+1 ) = D( I+1 ) - E( I )*EI
163: *
164: IF( D( I+1 ).LE.ZERO ) THEN
165: INFO = I + 1
166: GO TO 30
167: END IF
168: *
169: * Solve for e(i+1) and d(i+2).
170: *
171: EI = E( I+1 )
172: E( I+1 ) = EI / D( I+1 )
173: D( I+2 ) = D( I+2 ) - E( I+1 )*EI
174: *
175: IF( D( I+2 ).LE.ZERO ) THEN
176: INFO = I + 2
177: GO TO 30
178: END IF
179: *
180: * Solve for e(i+2) and d(i+3).
181: *
182: EI = E( I+2 )
183: E( I+2 ) = EI / D( I+2 )
184: D( I+3 ) = D( I+3 ) - E( I+2 )*EI
185: *
186: IF( D( I+3 ).LE.ZERO ) THEN
187: INFO = I + 3
188: GO TO 30
189: END IF
190: *
191: * Solve for e(i+3) and d(i+4).
192: *
193: EI = E( I+3 )
194: E( I+3 ) = EI / D( I+3 )
195: D( I+4 ) = D( I+4 ) - E( I+3 )*EI
196: 20 CONTINUE
197: *
198: * Check d(n) for positive definiteness.
199: *
200: IF( D( N ).LE.ZERO )
201: $ INFO = N
202: *
203: 30 CONTINUE
204: RETURN
205: *
206: * End of DPTTRF
207: *
208: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dpttrf.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack