?potrf
?potrf
Computes the Cholesky factorization of a symmetric (Hermitian) positive-definite matrix.
Syntax
lapack_int
LAPACKE_spotrf
(
int
matrix_layout
,
char
uplo
,
lapack_int
n
,
float
*
a
,
lapack_int
lda
);
lapack_int
LAPACKE_dpotrf
(
int
matrix_layout
,
char
uplo
,
lapack_int
n
,
double
*
a
,
lapack_int
lda
);
lapack_int
LAPACKE_cpotrf
(
int
matrix_layout
,
char
uplo
,
lapack_int
n
,
lapack_complex_float
*
a
,
lapack_int
lda
);
lapack_int
LAPACKE_zpotrf
(
int
matrix_layout
,
char
uplo
,
lapack_int
n
,
lapack_complex_double
*
a
,
lapack_int
lda
);
Include Files
- mkl.h
Description
The routine forms the Cholesky factorization of a symmetric positive-definite or, for complex data, Hermitian positive-definite matrix
A
:A = U T U A = U H U | if uplo ='U' |
A = L *L T A = L *L H | if uplo ='L' |
where
L
is a lower triangular matrix and U
is upper triangular. This routine supports the Progress Routine feature. See Progress Function
for details.
Input Parameters
- matrix_layout
- Specifies whether matrix storage layout is row major (LAPACK_ROW_MAJOR) or column major (LAPACK_COL_MAJOR).
- uplo
- Must be'U'or'L'.Indicates whether the upper or lower triangular part ofAis stored and howAis factored:Ifuplo='U', the arrayastores the upper triangular part of the matrixA, and the strictly lower triangular part of the matrix is not referenced.Ifuplo='L', the arrayastores the lower triangular part of the matrixA, and the strictly upper triangular part of the matrix is not referenced.
- n
- Specifies the order of the matrixA. The value ofnmust be at least zero.
- a
- Array, size max(1,lda*n. The arrayacontains either the upper or the lower triangular part of the matrixA(seeuplo).
- lda
- The leading dimension ofa. Must be at least max(1,n).
Output Parameters
- a
- The upper or lower triangular part ofais overwritten by the Cholesky factorUorL, as specified byuplo.
Return Values
This function returns a value
info
.If , the execution is successful.
info
=0If , parameter
info
= -i
i
had an illegal value. If , the leading minor of order
info
= i
i
(and therefore the matrix A
itself) is not positive-definite, and the factorization could not be completed. This may indicate an error in forming the matrix A
.Application Notes
If , where
uplo
= 'U'
, the computed factor U
is the exact factor of a perturbed matrix A
+ E

c
(n
)n
, and ε
is the machine precision.A similar estimate holds for .
uplo
= 'L'
The total number of floating-point operations is approximately
(1/3)
for real flavors or n
3
(4/3)
for complex
flavors.n
3
After calling this routine, you can call the following routines: