Annotation of rpl/lapack/lapack/dpotf2.f, revision 1.19
1.12 bertrand 1: *> \brief \b DPOTF2 computes the Cholesky factorization of a symmetric/Hermitian positive definite matrix (unblocked algorithm).
1.9 bertrand 2: *
3: * =========== DOCUMENTATION ===========
4: *
1.16 bertrand 5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
1.9 bertrand 7: *
8: *> \htmlonly
1.16 bertrand 9: *> Download DPOTF2 + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dpotf2.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dpotf2.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dpotf2.f">
1.9 bertrand 15: *> [TXT]</a>
1.16 bertrand 16: *> \endhtmlonly
1.9 bertrand 17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE DPOTF2( UPLO, N, A, LDA, INFO )
1.16 bertrand 22: *
1.9 bertrand 23: * .. Scalar Arguments ..
24: * CHARACTER UPLO
25: * INTEGER INFO, LDA, N
26: * ..
27: * .. Array Arguments ..
28: * DOUBLE PRECISION A( LDA, * )
29: * ..
1.16 bertrand 30: *
1.9 bertrand 31: *
32: *> \par Purpose:
33: * =============
34: *>
35: *> \verbatim
36: *>
37: *> DPOTF2 computes the Cholesky factorization of a real symmetric
38: *> positive definite matrix A.
39: *>
40: *> The factorization has the form
41: *> A = U**T * U , if UPLO = 'U', or
42: *> A = L * L**T, if UPLO = 'L',
43: *> where U is an upper triangular matrix and L is lower triangular.
44: *>
45: *> This is the unblocked version of the algorithm, calling Level 2 BLAS.
46: *> \endverbatim
47: *
48: * Arguments:
49: * ==========
50: *
51: *> \param[in] UPLO
52: *> \verbatim
53: *> UPLO is CHARACTER*1
54: *> Specifies whether the upper or lower triangular part of the
55: *> symmetric matrix A is stored.
56: *> = 'U': Upper triangular
57: *> = 'L': Lower triangular
58: *> \endverbatim
59: *>
60: *> \param[in] N
61: *> \verbatim
62: *> N is INTEGER
63: *> The order of the matrix A. N >= 0.
64: *> \endverbatim
65: *>
66: *> \param[in,out] A
67: *> \verbatim
68: *> A is DOUBLE PRECISION array, dimension (LDA,N)
69: *> On entry, the symmetric matrix A. If UPLO = 'U', the leading
70: *> n by n upper triangular part of A contains the upper
71: *> triangular part of the matrix A, and the strictly lower
72: *> triangular part of A is not referenced. If UPLO = 'L', the
73: *> leading n by n lower triangular part of A contains the lower
74: *> triangular part of the matrix A, and the strictly upper
75: *> triangular part of A is not referenced.
76: *>
77: *> On exit, if INFO = 0, the factor U or L from the Cholesky
78: *> factorization A = U**T *U or A = L*L**T.
79: *> \endverbatim
80: *>
81: *> \param[in] LDA
82: *> \verbatim
83: *> LDA is INTEGER
84: *> The leading dimension of the array A. LDA >= max(1,N).
85: *> \endverbatim
86: *>
87: *> \param[out] INFO
88: *> \verbatim
89: *> INFO is INTEGER
90: *> = 0: successful exit
91: *> < 0: if INFO = -k, the k-th argument had an illegal value
92: *> > 0: if INFO = k, the leading minor of order k is not
93: *> positive definite, and the factorization could not be
94: *> completed.
95: *> \endverbatim
96: *
97: * Authors:
98: * ========
99: *
1.16 bertrand 100: *> \author Univ. of Tennessee
101: *> \author Univ. of California Berkeley
102: *> \author Univ. of Colorado Denver
103: *> \author NAG Ltd.
1.9 bertrand 104: *
105: *> \ingroup doublePOcomputational
106: *
107: * =====================================================================
1.1 bertrand 108: SUBROUTINE DPOTF2( UPLO, N, A, LDA, INFO )
109: *
1.19 ! bertrand 110: * -- LAPACK computational routine --
1.1 bertrand 111: * -- LAPACK is a software package provided by Univ. of Tennessee, --
112: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
113: *
114: * .. Scalar Arguments ..
115: CHARACTER UPLO
116: INTEGER INFO, LDA, N
117: * ..
118: * .. Array Arguments ..
119: DOUBLE PRECISION A( LDA, * )
120: * ..
121: *
122: * =====================================================================
123: *
124: * .. Parameters ..
125: DOUBLE PRECISION ONE, ZERO
126: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
127: * ..
128: * .. Local Scalars ..
129: LOGICAL UPPER
130: INTEGER J
131: DOUBLE PRECISION AJJ
132: * ..
133: * .. External Functions ..
134: LOGICAL LSAME, DISNAN
135: DOUBLE PRECISION DDOT
136: EXTERNAL LSAME, DDOT, DISNAN
137: * ..
138: * .. External Subroutines ..
139: EXTERNAL DGEMV, DSCAL, XERBLA
140: * ..
141: * .. Intrinsic Functions ..
142: INTRINSIC MAX, SQRT
143: * ..
144: * .. Executable Statements ..
145: *
146: * Test the input parameters.
147: *
148: INFO = 0
149: UPPER = LSAME( UPLO, 'U' )
150: IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
151: INFO = -1
152: ELSE IF( N.LT.0 ) THEN
153: INFO = -2
154: ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
155: INFO = -4
156: END IF
157: IF( INFO.NE.0 ) THEN
158: CALL XERBLA( 'DPOTF2', -INFO )
159: RETURN
160: END IF
161: *
162: * Quick return if possible
163: *
164: IF( N.EQ.0 )
165: $ RETURN
166: *
167: IF( UPPER ) THEN
168: *
1.8 bertrand 169: * Compute the Cholesky factorization A = U**T *U.
1.1 bertrand 170: *
171: DO 10 J = 1, N
172: *
173: * Compute U(J,J) and test for non-positive-definiteness.
174: *
175: AJJ = A( J, J ) - DDOT( J-1, A( 1, J ), 1, A( 1, J ), 1 )
176: IF( AJJ.LE.ZERO.OR.DISNAN( AJJ ) ) THEN
177: A( J, J ) = AJJ
178: GO TO 30
179: END IF
180: AJJ = SQRT( AJJ )
181: A( J, J ) = AJJ
182: *
183: * Compute elements J+1:N of row J.
184: *
185: IF( J.LT.N ) THEN
186: CALL DGEMV( 'Transpose', J-1, N-J, -ONE, A( 1, J+1 ),
187: $ LDA, A( 1, J ), 1, ONE, A( J, J+1 ), LDA )
188: CALL DSCAL( N-J, ONE / AJJ, A( J, J+1 ), LDA )
189: END IF
190: 10 CONTINUE
191: ELSE
192: *
1.8 bertrand 193: * Compute the Cholesky factorization A = L*L**T.
1.1 bertrand 194: *
195: DO 20 J = 1, N
196: *
197: * Compute L(J,J) and test for non-positive-definiteness.
198: *
199: AJJ = A( J, J ) - DDOT( J-1, A( J, 1 ), LDA, A( J, 1 ),
200: $ LDA )
201: IF( AJJ.LE.ZERO.OR.DISNAN( AJJ ) ) THEN
202: A( J, J ) = AJJ
203: GO TO 30
204: END IF
205: AJJ = SQRT( AJJ )
206: A( J, J ) = AJJ
207: *
208: * Compute elements J+1:N of column J.
209: *
210: IF( J.LT.N ) THEN
211: CALL DGEMV( 'No transpose', N-J, J-1, -ONE, A( J+1, 1 ),
212: $ LDA, A( J, 1 ), LDA, ONE, A( J+1, J ), 1 )
213: CALL DSCAL( N-J, ONE / AJJ, A( J+1, J ), 1 )
214: END IF
215: 20 CONTINUE
216: END IF
217: GO TO 40
218: *
219: 30 CONTINUE
220: INFO = J
221: *
222: 40 CONTINUE
223: RETURN
224: *
225: * End of DPOTF2
226: *
227: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dpotf2.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack