potrs_batch (Buffer Strided Version)#
Solves a system of linear equations with a Cholesky-factored symmetric
(Hermitian) positive-definite coefficient matrices. This routine belongs
to the oneapi::mkl::lapack namespace.
Description#
The routine solves for Xi the system of linear equations
Ai*Xi = Bi with a
symmetric positive-definite or, for complex data, Hermitian
positive-definite matrices Ai, given the Cholesky
factorization of Ai, i ϵ{1...batch_size} :
Ai =UiT*Ui for real data,Ai =UiH*Ui for complex data if uplo=mkl::uplo::upperAi =Li*LiT for real data,Ai =Li*LiH for complex data if uplo=mkl::uplo::lower
where Li is a lower triangular matrix and
Ui is upper triangular. The system is solved with
multiple right-hand sides stored in the columns of the matrix
Bi.
Before calling this routine, matrices Ai should be
factorized by a call to potrf_batch (Buffer Strided Version).
API#
Syntax#
namespace oneapi::mkl::lapack {
void potrs_batch(sycl::queue &queue,
mkl::uplo uplo,
int64_t n,
int64_t nrhs,
sycl::buffer<T> &a,
int64_t lda,
int64_t stride_a,
sycl::buffer<T> &b,
int64_t ldb,
int64_t stride_b,
int64_t batch_size,
sycl::buffer<T> &scratchpad,
int64_t scratchpad_size)
}
This function supports the following precisions and devices:
T |
Devices supported |
|---|---|
|
CPU and GPU |
|
CPU and GPU |
|
CPU and GPU |
|
CPU and GPU |
Input Parameters#
- queue
Device queue where calculations will be performed.
- uplo
Indicates how the input matrix has been factored:
If uplo=
mkl::uplo::upper, the upper triangleUi ofAi is stored, whereAi =UiT*Ui for real data,Ai =UiH*Ui for complex data.If uplo=
mkl::uplo::lower, the upper triangleLi ofAi is stored, whereAi =Li*LiT for real data,Ai =Li*LiH for complex data.- n
The order of the matrices
Ai (n ≥ 0).- nrhs
The number of right hand sides
(nrhs ≥ 0).- a
Array containing the batch of factorizations of the matrices
Ai, as returned by potrf_batch (Buffer Strided Version).- lda
The leading dimension of
Ai (lda≥max(1, n)).- stride_a
The stride between the beginnings of matrices inside the batch array
a(stride_a≥max(1, lda * n)).- b
The array containing the batch of matrices
Bi whose columns are the right-hand sides for the systems of equations.- ldb
The leading dimensions of
Bi (ldb≥max(1, n)).- stride_b
The stride between the beginnings of matrices
Bi inside the batch arrayb(stride_b≥max(1, ldb * nrhs)).- batch_size
The number of problems in a batch (
batch_size≥ 0).- scratchpad
Scratchpad memory to be used by routine for storing intermediate results.
- scratchpad_size
Size of scratchpad memory as a number of floating point elements of type T. Size should not be less then the value returned by stride version of potrs_batch_scratchpad_size (Strided Version) function.
Output Parameters#
- b
The batch array b is overwritten by the solution matrix
Xi.
Exceptions#
Exception |
Description |
|---|---|
|
This exception is thrown when problems occur during calculations. You can obtain the info code of the problem using the info() method of the exception object: If If If |