LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ dpotrf()

subroutine dpotrf ( character  UPLO,
integer  N,
double precision, dimension( lda, * )  A,
integer  LDA,
integer  INFO 
)

DPOTRF

DPOTRF VARIANT: top-looking block version of the algorithm, calling Level 3 BLAS.

Download DPOTRF + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 DPOTRF computes the Cholesky factorization of a real symmetric
 positive definite matrix A.

 The factorization has the form
    A = U**T * U,  if UPLO = 'U', or
    A = L  * L**T,  if UPLO = 'L',
 where U is an upper triangular matrix and L is lower triangular.

 This is the block version of the algorithm, calling Level 3 BLAS.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is DOUBLE PRECISION array, dimension (LDA,N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the leading
          N-by-N upper triangular part of A contains the upper
          triangular part of the matrix A, and the strictly lower
          triangular part of A is not referenced.  If UPLO = 'L', the
          leading N-by-N lower triangular part of A contains the lower
          triangular part of the matrix A, and the strictly upper
          triangular part of A is not referenced.

          On exit, if INFO = 0, the factor U or L from the Cholesky
          factorization A = U**T*U or A = L*L**T.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  if INFO = i, the leading minor of order i is not
                positive definite, and the factorization could not be
                completed.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Purpose:

 DPOTRF computes the Cholesky factorization of a real symmetric
 positive definite matrix A.

 The factorization has the form
    A = U**T * U,  if UPLO = 'U', or
    A = L  * L**T,  if UPLO = 'L',
 where U is an upper triangular matrix and L is lower triangular.

 This is the top-looking block version of the algorithm, calling Level 3 BLAS.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is DOUBLE PRECISION array, dimension (LDA,N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the leading
          N-by-N upper triangular part of A contains the upper
          triangular part of the matrix A, and the strictly lower
          triangular part of A is not referenced.  If UPLO = 'L', the
          leading N-by-N lower triangular part of A contains the lower
          triangular part of the matrix A, and the strictly upper
          triangular part of A is not referenced.
          On exit, if INFO = 0, the factor U or L from the Cholesky
          factorization A = U**T*U or A = L*L**T.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  if INFO = i, the leading minor of order i is not
                positive definite, and the factorization could not be
                completed.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
December 2016

Definition at line 106 of file dpotrf.f.

107 *
108 * -- LAPACK computational routine --
109 * -- LAPACK is a software package provided by Univ. of Tennessee, --
110 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
111 *
112 * .. Scalar Arguments ..
113  CHARACTER UPLO
114  INTEGER INFO, LDA, N
115 * ..
116 * .. Array Arguments ..
117  DOUBLE PRECISION A( LDA, * )
118 * ..
119 *
120 * =====================================================================
121 *
122 * .. Parameters ..
123  DOUBLE PRECISION ONE
124  parameter( one = 1.0d+0 )
125 * ..
126 * .. Local Scalars ..
127  LOGICAL UPPER
128  INTEGER J, JB, NB
129 * ..
130 * .. External Functions ..
131  LOGICAL LSAME
132  INTEGER ILAENV
133  EXTERNAL lsame, ilaenv
134 * ..
135 * .. External Subroutines ..
136  EXTERNAL dgemm, dpotrf2, dsyrk, dtrsm, xerbla
137 * ..
138 * .. Intrinsic Functions ..
139  INTRINSIC max, min
140 * ..
141 * .. Executable Statements ..
142 *
143 * Test the input parameters.
144 *
145  info = 0
146  upper = lsame( uplo, 'U' )
147  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
148  info = -1
149  ELSE IF( n.LT.0 ) THEN
150  info = -2
151  ELSE IF( lda.LT.max( 1, n ) ) THEN
152  info = -4
153  END IF
154  IF( info.NE.0 ) THEN
155  CALL xerbla( 'DPOTRF', -info )
156  RETURN
157  END IF
158 *
159 * Quick return if possible
160 *
161  IF( n.EQ.0 )
162  $ RETURN
163 *
164 * Determine the block size for this environment.
165 *
166  nb = ilaenv( 1, 'DPOTRF', uplo, n, -1, -1, -1 )
167  IF( nb.LE.1 .OR. nb.GE.n ) THEN
168 *
169 * Use unblocked code.
170 *
171  CALL dpotrf2( uplo, n, a, lda, info )
172  ELSE
173 *
174 * Use blocked code.
175 *
176  IF( upper ) THEN
177 *
178 * Compute the Cholesky factorization A = U**T*U.
179 *
180  DO 10 j = 1, n, nb
181 *
182 * Update and factorize the current diagonal block and test
183 * for non-positive-definiteness.
184 *
185  jb = min( nb, n-j+1 )
186  CALL dsyrk( 'Upper', 'Transpose', jb, j-1, -one,
187  $ a( 1, j ), lda, one, a( j, j ), lda )
188  CALL dpotrf2( 'Upper', jb, a( j, j ), lda, info )
189  IF( info.NE.0 )
190  $ GO TO 30
191  IF( j+jb.LE.n ) THEN
192 *
193 * Compute the current block row.
194 *
195  CALL dgemm( 'Transpose', 'No transpose', jb, n-j-jb+1,
196  $ j-1, -one, a( 1, j ), lda, a( 1, j+jb ),
197  $ lda, one, a( j, j+jb ), lda )
198  CALL dtrsm( 'Left', 'Upper', 'Transpose', 'Non-unit',
199  $ jb, n-j-jb+1, one, a( j, j ), lda,
200  $ a( j, j+jb ), lda )
201  END IF
202  10 CONTINUE
203 *
204  ELSE
205 *
206 * Compute the Cholesky factorization A = L*L**T.
207 *
208  DO 20 j = 1, n, nb
209 *
210 * Update and factorize the current diagonal block and test
211 * for non-positive-definiteness.
212 *
213  jb = min( nb, n-j+1 )
214  CALL dsyrk( 'Lower', 'No transpose', jb, j-1, -one,
215  $ a( j, 1 ), lda, one, a( j, j ), lda )
216  CALL dpotrf2( 'Lower', jb, a( j, j ), lda, info )
217  IF( info.NE.0 )
218  $ GO TO 30
219  IF( j+jb.LE.n ) THEN
220 *
221 * Compute the current block column.
222 *
223  CALL dgemm( 'No transpose', 'Transpose', n-j-jb+1, jb,
224  $ j-1, -one, a( j+jb, 1 ), lda, a( j, 1 ),
225  $ lda, one, a( j+jb, j ), lda )
226  CALL dtrsm( 'Right', 'Lower', 'Transpose', 'Non-unit',
227  $ n-j-jb+1, jb, one, a( j, j ), lda,
228  $ a( j+jb, j ), lda )
229  END IF
230  20 CONTINUE
231  END IF
232  END IF
233  GO TO 40
234 *
235  30 CONTINUE
236  info = info + j - 1
237 *
238  40 CONTINUE
239  RETURN
240 *
241 * End of DPOTRF
242 *
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: ilaenv.f:162
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine dtrsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
DTRSM
Definition: dtrsm.f:181
subroutine dsyrk(UPLO, TRANS, N, K, ALPHA, A, LDA, BETA, C, LDC)
DSYRK
Definition: dsyrk.f:169
subroutine dgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
DGEMM
Definition: dgemm.f:187
recursive subroutine dpotrf2(UPLO, N, A, LDA, INFO)
DPOTRF2
Definition: dpotrf2.f:106
Here is the call graph for this function:
Here is the caller graph for this function: