Annotation of rpl/lapack/lapack/dpotrf2.f, revision 1.3
1.1 bertrand 1: *> \brief \b DPOTRF2
2: *
3: * =========== DOCUMENTATION ===========
4: *
1.3 ! bertrand 5: * Online html documentation available at
! 6: * http://www.netlib.org/lapack/explore-html/
1.1 bertrand 7: *
8: * Definition:
9: * ===========
10: *
11: * RECURSIVE SUBROUTINE DPOTRF2( UPLO, N, A, LDA, INFO )
1.3 ! bertrand 12: *
1.1 bertrand 13: * .. Scalar Arguments ..
14: * CHARACTER UPLO
15: * INTEGER INFO, LDA, N
16: * ..
17: * .. Array Arguments ..
18: * REAL A( LDA, * )
19: * ..
1.3 ! bertrand 20: *
1.1 bertrand 21: *
22: *> \par Purpose:
23: * =============
24: *>
25: *> \verbatim
26: *>
27: *> DPOTRF2 computes the Cholesky factorization of a real symmetric
28: *> positive definite matrix A using the recursive algorithm.
29: *>
30: *> The factorization has the form
31: *> A = U**T * U, if UPLO = 'U', or
32: *> A = L * L**T, if UPLO = 'L',
33: *> where U is an upper triangular matrix and L is lower triangular.
34: *>
35: *> This is the recursive version of the algorithm. It divides
36: *> the matrix into four submatrices:
37: *>
38: *> [ A11 | A12 ] where A11 is n1 by n1 and A22 is n2 by n2
39: *> A = [ -----|----- ] with n1 = n/2
40: *> [ A21 | A22 ] n2 = n-n1
41: *>
42: *> The subroutine calls itself to factor A11. Update and scale A21
43: *> or A12, update A22 then calls itself to factor A22.
1.3 ! bertrand 44: *>
1.1 bertrand 45: *> \endverbatim
46: *
47: * Arguments:
48: * ==========
49: *
50: *> \param[in] UPLO
51: *> \verbatim
52: *> UPLO is CHARACTER*1
53: *> = 'U': Upper triangle of A is stored;
54: *> = 'L': Lower triangle of A is stored.
55: *> \endverbatim
56: *>
57: *> \param[in] N
58: *> \verbatim
59: *> N is INTEGER
60: *> The order of the matrix A. N >= 0.
61: *> \endverbatim
62: *>
63: *> \param[in,out] A
64: *> \verbatim
65: *> A is DOUBLE PRECISION array, dimension (LDA,N)
66: *> On entry, the symmetric matrix A. If UPLO = 'U', the leading
67: *> N-by-N upper triangular part of A contains the upper
68: *> triangular part of the matrix A, and the strictly lower
69: *> triangular part of A is not referenced. If UPLO = 'L', the
70: *> leading N-by-N lower triangular part of A contains the lower
71: *> triangular part of the matrix A, and the strictly upper
72: *> triangular part of A is not referenced.
73: *>
74: *> On exit, if INFO = 0, the factor U or L from the Cholesky
75: *> factorization A = U**T*U or A = L*L**T.
76: *> \endverbatim
77: *>
78: *> \param[in] LDA
79: *> \verbatim
80: *> LDA is INTEGER
81: *> The leading dimension of the array A. LDA >= max(1,N).
82: *> \endverbatim
83: *>
84: *> \param[out] INFO
85: *> \verbatim
86: *> INFO is INTEGER
87: *> = 0: successful exit
88: *> < 0: if INFO = -i, the i-th argument had an illegal value
89: *> > 0: if INFO = i, the leading minor of order i is not
90: *> positive definite, and the factorization could not be
91: *> completed.
92: *> \endverbatim
93: *
94: * Authors:
95: * ========
96: *
1.3 ! bertrand 97: *> \author Univ. of Tennessee
! 98: *> \author Univ. of California Berkeley
! 99: *> \author Univ. of Colorado Denver
! 100: *> \author NAG Ltd.
1.1 bertrand 101: *
1.3 ! bertrand 102: *> \date December 2016
1.1 bertrand 103: *
104: *> \ingroup doublePOcomputational
105: *
106: * =====================================================================
107: RECURSIVE SUBROUTINE DPOTRF2( UPLO, N, A, LDA, INFO )
108: *
1.3 ! bertrand 109: * -- LAPACK computational routine (version 3.7.0) --
1.1 bertrand 110: * -- LAPACK is a software package provided by Univ. of Tennessee, --
111: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.3 ! bertrand 112: * December 2016
1.1 bertrand 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 ..
1.3 ! bertrand 129: LOGICAL UPPER
1.1 bertrand 130: INTEGER N1, N2, IINFO
131: * ..
132: * .. External Functions ..
133: LOGICAL LSAME, DISNAN
134: EXTERNAL LSAME, DISNAN
135: * ..
136: * .. External Subroutines ..
137: EXTERNAL DSYRK, DTRSM, XERBLA
138: * ..
139: * .. Intrinsic Functions ..
140: INTRINSIC MAX, SQRT
141: * ..
142: * .. Executable Statements ..
143: *
144: * Test the input parameters
145: *
146: INFO = 0
147: UPPER = LSAME( UPLO, 'U' )
148: IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
149: INFO = -1
150: ELSE IF( N.LT.0 ) THEN
151: INFO = -2
152: ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
153: INFO = -4
154: END IF
155: IF( INFO.NE.0 ) THEN
156: CALL XERBLA( 'DPOTRF2', -INFO )
157: RETURN
158: END IF
159: *
160: * Quick return if possible
161: *
162: IF( N.EQ.0 )
163: $ RETURN
164: *
165: * N=1 case
166: *
167: IF( N.EQ.1 ) THEN
168: *
169: * Test for non-positive-definiteness
170: *
171: IF( A( 1, 1 ).LE.ZERO.OR.DISNAN( A( 1, 1 ) ) ) THEN
172: INFO = 1
173: RETURN
174: END IF
175: *
176: * Factor
177: *
178: A( 1, 1 ) = SQRT( A( 1, 1 ) )
179: *
180: * Use recursive code
181: *
182: ELSE
183: N1 = N/2
184: N2 = N-N1
185: *
186: * Factor A11
187: *
188: CALL DPOTRF2( UPLO, N1, A( 1, 1 ), LDA, IINFO )
189: IF ( IINFO.NE.0 ) THEN
190: INFO = IINFO
191: RETURN
1.3 ! bertrand 192: END IF
1.1 bertrand 193: *
194: * Compute the Cholesky factorization A = U**T*U
195: *
196: IF( UPPER ) THEN
197: *
198: * Update and scale A12
199: *
200: CALL DTRSM( 'L', 'U', 'T', 'N', N1, N2, ONE,
1.3 ! bertrand 201: $ A( 1, 1 ), LDA, A( 1, N1+1 ), LDA )
1.1 bertrand 202: *
203: * Update and factor A22
1.3 ! bertrand 204: *
1.1 bertrand 205: CALL DSYRK( UPLO, 'T', N2, N1, -ONE, A( 1, N1+1 ), LDA,
206: $ ONE, A( N1+1, N1+1 ), LDA )
207: CALL DPOTRF2( UPLO, N2, A( N1+1, N1+1 ), LDA, IINFO )
208: IF ( IINFO.NE.0 ) THEN
209: INFO = IINFO + N1
210: RETURN
211: END IF
212: *
213: * Compute the Cholesky factorization A = L*L**T
214: *
215: ELSE
216: *
217: * Update and scale A21
218: *
1.3 ! bertrand 219: CALL DTRSM( 'R', 'L', 'T', 'N', N2, N1, ONE,
1.1 bertrand 220: $ A( 1, 1 ), LDA, A( N1+1, 1 ), LDA )
221: *
222: * Update and factor A22
223: *
224: CALL DSYRK( UPLO, 'N', N2, N1, -ONE, A( N1+1, 1 ), LDA,
225: $ ONE, A( N1+1, N1+1 ), LDA )
226: CALL DPOTRF2( UPLO, N2, A( N1+1, N1+1 ), LDA, IINFO )
227: IF ( IINFO.NE.0 ) THEN
228: INFO = IINFO + N1
229: RETURN
230: END IF
231: END IF
232: END IF
233: RETURN
234: *
235: * End of DPOTRF2
236: *
237: END
CVSweb interface <joel.bertrand@systella.fr>