potrf#

Computes the Cholesky factorization of a symmetric (Hermitian) positive-definite matrix.

Description

potrf supports the following precisions.

T

float

double

std::complex<float>

std::complex<double>

The routine forms the Cholesky factorization of a symmetric positive-definite or, for complex data, Hermitian positive-definite matrix \(A\):

\(A\) = \(U^{T}U\) for real data, \(A = U^{H}U\) for complex data

if upper_lower=oneapi::mkl::uplo::upper

\(A\) = \(LL^{T}\) for real data, \(A = LL^{H}\) for complex data

if upper_lower=oneapi::mkl::uplo::lower

where \(L\) is a lower triangular matrix and \(U\) is upper triangular.

potrf (Buffer Version)#

Syntax

namespace oneapi::mkl::lapack {
  void potrf(sycl::queue &queue, oneapi::mkl::uplo upper_lower, std::int64_t n, sycl::buffer<T,1> &a, std::int64_t lda, sycl::buffer<T,1> &scratchpad, std::int64_t scratchpad_size)
}

Input Parameters

queue

The queue where the routine should be executed.

upper_lower

Indicates whether the upper or lower triangular part of \(A\) is stored and how \(A\) is factored:

If upper_lower=oneapi::mkl::uplo::upper, the array a stores the upper triangular part of the matrix \(A\), and the strictly lower triangular part of the matrix is not referenced.

If upper_lower=oneapi::mkl::uplo::lower, the array a stores the lower triangular part of the matrix \(A\), and the strictly upper triangular part of the matrix is not referenced.

n

Specifies the order of the matrix \(A\) (\(0 \le n\)).

a

Buffer holding input matrix \(A\). The buffer a contains either the upper or the lower triangular part of the matrix \(A\) (see upper_lower). The second dimension of a must be at least \(\max(1, n)\).

lda

The leading dimension of a.

scratchpad_size

Size of scratchpad memory as a number of floating point elements of type T. Size should not be less than the value returned by potrf_scratchpad_size function.

Output Parameters

a

The buffer a is overwritten by the Cholesky factor \(U\) or \(L\), as specified by upper_lower.

scratchpad

Buffer holding scratchpad memory to be used by routine for storing intermediate results.

potrf (USM Version)#

Syntax

namespace oneapi::mkl::lapack {
  sycl::event potrf(sycl::queue &queue, oneapi::mkl::uplo upper_lower, std::int64_t n, T *a, std::int64_t lda, T *scratchpad, std::int64_t scratchpad_size, const std::vector<sycl::event> &events = {})
}

Input Parameters

queue

The queue where the routine should be executed.

upper_lower

Indicates whether the upper or lower triangular part of \(A\) is stored and how \(A\) is factored:

If upper_lower=oneapi::mkl::uplo::upper, the array a stores the upper triangular part of the matrix \(A\), and the strictly lower triangular part of the matrix is not referenced.

If upper_lower=oneapi::mkl::uplo::lower, the array a stores the lower triangular part of the matrix \(A\), and the strictly upper triangular part of the matrix is not referenced.

n

Specifies the order of the matrix \(A\) (\(0 \le n\)).

a

Pointer to input matrix \(A\). The array a contains either the upper or the lower triangular part of the matrix \(A\) (see upper_lower). The second dimension of a must be at least \(\max(1, n)\).

lda

The leading dimension of a.

scratchpad_size

Size of scratchpad memory as a number of floating point elements of type T. Size should not be less than the value returned by potrf_scratchpad_size function.

events

List of events to wait for before starting computation. Defaults to empty list.

Output Parameters

a

The memory pointer to by pointer a is overwritten by the Cholesky factor \(U\) or \(L\), as specified by upper_lower.

scratchpad

Pointer to scratchpad memory to be used by routine for storing intermediate results.

Return Values

Output event to wait on to ensure computation is complete.

Parent topic: LAPACK Linear Equation Routines