Linear Solvers¶
The linear solver interface provides methods for users to run a linear solver using either cuBLAS or cuSolver as the backend.
Cached API¶
-
template<typename OutputTensor, typename ATensor>
void matx::chol(OutputTensor &out, const ATensor &a, cudaStream_t stream = 0, cublasFillMode_t uplo = CUBLAS_FILL_MODE_UPPER)¶ Perform a Cholesky decomposition using a cached plan
See documentation of matxDnCholSolverPlan_t for a description of how the algorithm works. This function provides a simple interface to the cuSolver library by deducing all parameters needed to perform a Cholesky decomposition from only the matrix A. The input and output parameters may be the same tensor. In that case, the input is destroyed and the output is stored in-place.
- Template Parameters
T1 – Data type of matrix A
RANK – Rank of matrix A
- Parameters
out – Output tensor
a – Input tensor
stream – CUDA stream
uplo – Part of matrix to fill
-
template<typename OutputTensor, typename PivotTensor, typename ATensor>
void matx::lu(OutputTensor &out, PivotTensor &piv, const ATensor &a, const cudaStream_t stream = 0)¶ Perform a LU decomposition using a cached plan
See documentation of matxDnLUSolverPlan_t for a description of how the algorithm works. This function provides a simple interface to the cuSolver library by deducing all parameters needed to perform an LU decomposition from only the matrix A. The input and output parameters may be the same tensor. In that case, the input is destroyed and the output is stored in-place.
- Template Parameters
T1 – Data type of matrix A
RANK – Rank of matrix A
- Parameters
out – Output tensor view
piv – Output of pivot indices
a – Input matrix A
stream – CUDA stream
-
template<typename OutTensor, typename TauTensor, typename ATensor>
void matx::qr(OutTensor &out, TauTensor &tau, const ATensor &a, cudaStream_t stream = 0)¶ Perform a QR decomposition using a cached plan
See documentation of matxDnQRSolverPlan_t for a description of how the algorithm works. This function provides a simple interface to the cuSolver library by deducing all parameters needed to perform a QR decomposition from only the matrix A. The input and output parameters may be the same tensor. In that case, the input is destroyed and the output is stored in-place.
- Template Parameters
T1 – Data type of matrix A
RANK – Rank of matrix A
- Parameters
out – Output tensor view
tau – Output of reflection scalar values
a – Input tensor A
stream – CUDA stream
-
template<typename UTensor, typename STensor, typename VTensor, typename ATensor>
void matx::svd(UTensor &u, STensor &s, VTensor &v, const ATensor &a, cudaStream_t stream = 0, const char jobu = 'A', const char jobvt = 'A')¶ Perform a SVD decomposition using a cached plan
See documentation of matxDnSVDSolverPlan_t for a description of how the algorithm works. This function provides a simple interface to the cuSolver library by deducing all parameters needed to perform a SVD decomposition from only the matrix A.
- Template Parameters
T1 – Data type of matrix A
RANK – Rank of matrix A
- Parameters
u – U matrix output
s – Sigma matrix output
v – V matrix output
a – Input matrix A
stream – CUDA stream
jobu – Specifies options for computing all or part of the matrix U: = ‘A’. See cuSolver documentation for more info
jobvt – specifies options for computing all or part of the matrix V**T. See cuSolver documentation for more info
-
template<typename OutputTensor, typename WTensor, typename ATensor>
void matx::eig([[maybe_unused]] OutputTensor &out, WTensor &w, const ATensor &a, cudaStream_t stream = 0, cusolverEigMode_t jobz = CUSOLVER_EIG_MODE_VECTOR, cublasFillMode_t uplo = CUBLAS_FILL_MODE_UPPER)¶ Perform a Eig decomposition using a cached plan
See documentation of matxDnEigSolverPlan_t for a description of how the algorithm works. This function provides a simple interface to the cuSolver library by deducing all parameters needed to perform a eigen decomposition from only the matrix A. The input and output parameters may be the same tensor. In that case, the input is destroyed and the output is stored in-place.
- Template Parameters
T1 – Data type of matrix A
RANK – Rank of matrix A
- Parameters
out – Output tensor view
w – Eigenvalues output
a – Input matrix A
stream – CUDA stream
jobz – CUSOLVER_EIG_MODE_VECTOR to compute eigenvectors or CUSOLVER_EIG_MODE_NOVECTOR to not compute
uplo – Where to store data in A
Non-Cached API¶
-
template<typename OutputTensor, typename ATensor>
class matx::matxDnCholSolverPlan_t : public matx::matxDnSolver_t¶ Public Functions
-
inline matxDnCholSolverPlan_t(const ATensor &a, cublasFillMode_t uplo = CUBLAS_FILL_MODE_UPPER)¶
Plan for solving \(\textbf{A} = \textbf{L} * \textbf{L^{H}}\) or \(\textbf{A} = \textbf{U} * \textbf{U^{H}}\) using the Cholesky method
Creates a handle for solving the factorization of A = M * M^H of a dense matrix using the Cholesky method, where M is either the upper or lower triangular portion of A. Input matrix A must be a square Hermitian matrix positive-definite where only the upper or lower triangle is used.
- Template Parameters
T1 – Data type of A matrix
RANK – Rank of A matrix
- Parameters
a – Input tensor view
uplo – Use upper or lower triangle for computation
-
inline ~matxDnCholSolverPlan_t()¶
Cholesky solver handle destructor
Destroys any helper data used for provider type and any workspace memory created
-
inline matxDnCholSolverPlan_t(const ATensor &a, cublasFillMode_t uplo = CUBLAS_FILL_MODE_UPPER)¶
-
template<typename OutputTensor, typename PivotTensor, typename ATensor>
class matx::matxDnLUSolverPlan_t : public matx::matxDnSolver_t¶ Public Functions
-
inline matxDnLUSolverPlan_t(PivotTensor &piv, const ATensor &a)¶
Plan for factoring A such that \(\textbf{P} * \textbf{A} = \textbf{L} * \textbf{U}\)
Creates a handle for factoring matrix A into the format above. Matrix must not be singular.
- Template Parameters
T1 – Data type of A matrix
RANK – Rank of A matrix
- Parameters
piv – Pivot indices
a – Input tensor view
-
inline ~matxDnLUSolverPlan_t()¶
LU solver handle destructor
Destroys any helper data used for provider type and any workspace memory created
-
inline matxDnLUSolverPlan_t(PivotTensor &piv, const ATensor &a)¶
-
template<typename OutTensor, typename TauTensor, typename ATensor>
class matx::matxDnQRSolverPlan_t : public matx::matxDnSolver_t¶ Public Functions
-
inline matxDnQRSolverPlan_t(TauTensor &tau, const ATensor &a)¶
Plan for factoring A such that \(\textbf{A} = \textbf{Q} * \textbf{R}\)
Creates a handle for factoring matrix A into the format above. QR decomposition in cuBLAS/cuSolver does not return the Q matrix directly, and it must be computed separately used the Householder reflections in the tau output, along with the overwritten A matrix input. The input and output parameters may be the same tensor. In that case, the input is destroyed and the output is stored in-place.
- Template Parameters
T1 – Data type of A matrix
RANK – Rank of A matrix
- Parameters
tau – Scaling factors for reflections
a – Input tensor view
-
inline ~matxDnQRSolverPlan_t()¶
QR solver handle destructor
Destroys any helper data used for provider type and any workspace memory created
-
inline matxDnQRSolverPlan_t(TauTensor &tau, const ATensor &a)¶
-
template<typename OutputTensor, typename WTensor, typename ATensor>
class matx::matxDnEigSolverPlan_t : public matx::matxDnSolver_t¶ Public Functions
-
inline matxDnEigSolverPlan_t(WTensor &w, const ATensor &a, cusolverEigMode_t jobz = CUSOLVER_EIG_MODE_VECTOR, cublasFillMode_t uplo = CUBLAS_FILL_MODE_UPPER)¶
Plan computing eigenvalues/vectors on A such that:
\(\textbf{A} * textbf{V} = \textbf{V} * \textbf{\Lambda}\)
- Template Parameters
T1 – Data type of A matrix
T2 – Data type of W matrix
RANK – Rank of A matrix
- Parameters
w – Eigenvalues of A
a – Input tensor view
jobz – CUSOLVER_EIG_MODE_VECTOR to compute eigenvectors or CUSOLVER_EIG_MODE_NOVECTOR to not compute
uplo – Where to store data in A
-
inline ~matxDnEigSolverPlan_t()¶
Eigen solver handle destructor
Destroys any helper data used for provider type and any workspace memory created
-
inline matxDnEigSolverPlan_t(WTensor &w, const ATensor &a, cusolverEigMode_t jobz = CUSOLVER_EIG_MODE_VECTOR, cublasFillMode_t uplo = CUBLAS_FILL_MODE_UPPER)¶