Matrix Decomposition

geqrf

geqrf overload (1)

#include "MatrixDecomposition/geqrf.hpp"
template <
    typename T,
    int NRMAX,
    int NCMAX,
    int NCU
    >
int geqrf (
    int m,
    int n,
    T* A,
    int lda,
    T* tau
    )

This function computes QR decomposition of matrix \(A\)

\[\begin{equation*} {A = Q R}\end{equation*}\]

where \(A\) is a dense matrix of size \(m \times n\) , \(Q\) is a \(m \times n\) matrix with orthonormal columns, and \(R\) is an upper triangular matrix.

The maximum matrix size supported in FPGA is templated by NRMAX and NCMAX.

Parameters:

T data type (support float and double)
NRMAX maximum number of rows of input matrix
NCMAX maximum number of columns of input matrix
NCU number of computation unit, NCU = n will be n times faster than NCU = 1, and cost n times resource.
m number of rows of matrix A
n number of cols of matrix A
A input matrix of size \(m \times lda\) , and overwritten by the output triangular R matrix and min(m,n) elementary reflectors
lda leading dimension of matrix A
tau scalar factors for elementary reflectors

gesvdj

#include "MatrixDecomposition/gesvdj.hpp"
template <
    typename T,
    int NMAX,
    int NCU
    >
void gesvdj (
    int m,
    T* A,
    int lda,
    T* S,
    T* U,
    int ldu,
    T* V,
    int ldv,
    int& info
    )

Symmetric Matrix Jacobi based Singular Value Decomposition (GESVDJ) .

\[\begin{equation*} {A = U \Sigma {V}^T}\end{equation*}\]

where \(A\) is a dense symmetric matrix of size \(m \times m\) , \(U\) and \(V\) are \(m \times m\) matrix with orthonormal columns, and \(\Sigma\) is diagonal matrix.

The maximum matrix size supported in FPGA is templated by NMAX.

Parameters:

T data type (support float and double).
NMAX maximum number of rows/columns of input matrix
NCU number of computation unit
m number of rows/cols of matrix A
A input matrix of size \(m \times m\)
S decomposed diagonal singular matrix of size \(m \times m\)
U left U matrix of SVD
V right V matrix of SVD
lda leading dimension of matrix A
ldu leading dimension of matrix U
ldv leading dimension of matrix V
info output info (unused)

gesvj

#include "MatrixDecomposition/gesvj.hpp"
template <
    typename T,
    int NRMAX,
    int NCMAX,
    int MCU,
    int NCU
    >
void gesvj (
    int m,
    int n,
    T* A,
    T* U,
    T* S,
    T* V
    )

This function implements singular value decomposition of matrix A using one-sided Jacobi algorihtm.

\[\begin{equation*} {A = U \Sigma {V}^T}\end{equation*}\]

where \(A\) is a dense matrix of size \(m \times n\) , \(U\) is \(m \times m\) matrix with orthonormal columns, \(V\) is \(n \times n\) matrix with orthonormal columns, and \(\Sigma\) is diagonal matrix.

The maximum matrix size supported in FPGA is templated by NCMAX, NRMAX.

Parameters:

T  
: the data type of gesvj
NRMAX maximum number of rows of input matrix
NCMAX maximum number of columns of input matrix
MCU number of computation unit of M
NCU number of computation unit of N
m number of rows of matrix A
n number of cols of matrix A
A input matrix of size \(m \times n\)
S decomposed diagonal singular matrix of size n
U left U matrix of SVD of size \(m \times m\)
V right V matrix of SVD \(n \times n\)

getrf

#include "MatrixDecomposition/getrf.hpp"
template <
    typename T,
    int NMAX,
    int NCU
    >
void getrf (
    int n,
    T* A,
    int lda,
    int* ipiv,
    int& info
    )

This function computes the LU decomposition (with partial pivoting) of matrix \(A\)

\[\begin{equation*} {P A = L U}\end{equation*}\]

where \(P\) is a permutation matrix, \(A\) is a dense matrix of size \(n \times n\) , \(L\) is a lower triangular matrix with unit diagonal, and \(U\) is an upper triangular matrix. This function does not implement pivoting.

The maximum matrix size supported in FPGA is templated by NMAX.

Parameters:

T data type (support float and double)
NMAX maximum number of rows/columns of input matrix
NCU number of computation unit
n number of rows/cols of matrix A
A input matrix, and overwritten by the output upper and lower triangular matrix
lda leading dimention of input matrix A
pivot indices, row i of matrix A is stored in row[i]
info output info (unused)

getrf_nopivot

#include "MatrixDecomposition/getrf_nopivot.hpp"
template <
    typename T,
    int NMAX,
    int NCU
    >
void getrf_nopivot (
    int n,
    T* A,
    int lda,
    int& info
    )

This function computes the LU decomposition (without pivoting) of matrix \(A\)

\[\begin{equation*} {A = L U}\end{equation*}\]

where \(A\) is a dense matrix of size \(n \times n\) , \(L\) is a lower triangular matrix with unit diagonal, and \(U\) is an upper triangular matrix. This function does not implement pivoting.

The maximum matrix size supported in FPGA is templated by NMAX.

Parameters:

T data type (support float and double)
NMAX maximum number of rows/cols of input matrix
NCU number of computation unit
n number of rows/cols of matrix A
A input matrix
lda leading dimention of input matrix A
info output info (unused)

potrf

#include "MatrixDecomposition/potrf.hpp"
template <
    typename T,
    int NMAX,
    int NCU
    >
void potrf (
    int m,
    T* A,
    int lda,
    int& info
    )

This function computes the Cholesky decomposition of matrix \(A\)

\[\begin{equation*} {A = L {L}^T}\end{equation*}\]

where \(A\) is a dense symmetric positive-definite matrix of size \(m \times m\) , \(L\) is a lower triangular matrix, and \({L}^T\) is the transposed matrix of \(L\) .

The maximum matrix size supported in FPGA is templated by NMAX.

Parameters:

T data type (support float and double)
NMAX maximum number of rows/columns of input matrix
NCU number of computation unit
m number of rows/cols of matrix A
A input matrix of size \(m \times m\)
lda leading dimention of input matrix A
info output info (unused)

Linear Solver

gelinearsolver

#include "LinearSolver/gelinearsolver.hpp"
template <
    typename T,
    int NMAX,
    int NCU
    >
void gelinearsolver (
    int n,
    T* A,
    int b,
    T* B,
    int lda,
    int ldb,
    int& info
    )

This function solves a system of linear equation with general matrix along with multiple right-hand side vector

\[\begin{equation*} {Ax=B}\end{equation*}\]

where \(A\) is a dense general matrix of size \(n \times n\) , \(x\) is a vector need to be computed, and \(B\) is input vector.

The maximum matrix size supported in FPGA is templated by NMAX.

Parameters:

T data type (support float and double)
NMAX maximum number of rows/columns of input matrix
NCU number of computation unit
n number of rows/cols of matrix A
A input matrix of size \(n \times n\)
b number of columns of matrix B
B input matrix of size \(b \times n\) , and overwritten by the output matrix x
lda leading dimention of input matrix A
ldb leading dimention of input matrix B
info output info (unused)

gematrixinverse

#include "LinearSolver/gematrixinverse.hpp"
template <
    typename T,
    int NMAX,
    int NCU
    >
void gematrixinverse (
    int m,
    T* A,
    int lda,
    int& info
    )

This function computes the inverse matrix of \(A\)

\[\begin{equation*} {A}^{-1}\end{equation*}\]

where \(A\) is a dense general matrix of size \(m \times m\) . The maximum matrix size supported in FPGA is templated by NMAX.

Parameters:

T data type (support float and double)
NMAX maximum number of rows/columns of input matrix
NCU number of computation unit
m number of rows/cols of matrix A
A input matrix of size \(n \times n\)
lda leading dimention of input matrix A
info output info (unused)

gtsv

#include "LinearSolver/gtsv_pcr.hpp"
template <
    typename T,
    unsigned int NMAX,
    unsigned int NCU
    >
int gtsv (
    unsigned int n,
    T* matDiagLow,
    T* matDiag,
    T* matDiagUp,
    T* rhs
    )

Tri-diagonal linear solver. Compute solution to linear system with a tridiagonal matrix. Parallel Cyclic Reduction method.

Parameters:

T data type (support float and double)
NMAX matrix size
NCU number of compute units
matDiagLow lower diagonal of matrix
matDiag diagonal of matrix
matDiagUp upper diagonal of matrix
rhs right-hand side

polinearsolver

#include "LinearSolver/polinearsolver.hpp"
template <
    typename T,
    int NMAX,
    int NCU
    >
void polinearsolver (
    int n,
    T* A,
    int b,
    T* B,
    int lda,
    int ldb,
    int& info
    )

This function solves a system of linear equation with symmetric positive definite (SPD) matrix along with multiple right-hand side vector

\[\begin{equation*} {Ax=B}\end{equation*}\]

where \(A\) is a dense SPD triangular matrix of size \(m \times m\) , \(x\) is a vector need to be computed, and \(B\) is input vector.

The maximum matrix size supported in FPGA is templated by NMAX.

Parameters:

T data type (support float and double)
NMAX maximum number of rows/columns of input matrix
NCU number of computation unit
n number of rows/cols of matrix A
A input matrix of size \(n \times n\)
b number of columns of matrix B
B input matrix of size \(b \times n\) , and overwritten by the output matrix x
lda leading dimention of input matrix A
ldb leading dimention of input matrix B
info output info (unused)

pomatrixinverse

#include "LinearSolver/pomatrixinverse.hpp"
template <
    typename T,
    int NMAX,
    int NCU
    >
void pomatrixinverse (
    int m,
    T* A,
    int lda,
    int& info
    )

This function computes the inverse matrix of \(A\)

\[\begin{equation*} {A}^{-1}\end{equation*}\]

where \(A\) is a dense symmetric positive-definite matrix of size \(m \times m\) . The maximum matrix size supported in FPGA is templated by NMAX.

Parameters:

T data type (support float and double)
NMAX maximum number of rows/columns of input matrix
NCU number of computation unit
m number of rows/cols of matrix A
A input matrix of size \(n \times n\)
lda leading dimention of input matrix A
info output info (unused)

trtrs

#include "LinearSolver/trtrs.hpp"
template <
    typename T,
    int NMAX,
    int NCU
    >
void trtrs (
    char uplo,
    int m,
    T* A,
    int b,
    T* B,
    int lda,
    int ldb,
    int& info
    )

This function solves a system of linear equation with triangular coefficient matrix along with multiple right-hand side vector

\[\begin{equation*} {Ax=B}\end{equation*}\]

where \(A\) is a dense lower or upper triangular matrix of size \(m \times m\) , \(x\) is a vector need to be computed, and \(B\) is input vector.

The maximum matrix size supported in FPGA is templated by NMAX.

Parameters:

T data type (support float and double)
NMAX maximum number of rows/columns of input matrix
NCU number of computation unit
m number of rows/cols of matrix A
A input matrix of size \(n \times n\)
b number of columns of matrix B
B input matrix of size \(b \times n\) , and overwritten by the output matrix x
lda leading dimention of input matrix A
ldb leading dimention of input matrix B
info output info (unused)

Eigenvalue Solver

syevj

#include "EigenSolver/syevj.hpp"
template <
    typename T,
    int NMAX,
    int NCU
    >
void syevj (
    int m,
    T* A,
    int lda,
    T* S,
    T* U,
    int ldu,
    int& info
    )

Symmetric Matrix Jacobi based Eigenvalue Decomposition (SYEVJ) .

\[\begin{equation*} {A U = U \Sigma, }\end{equation*}\]

where \(A\) is a dense symmetric matrix of size \(m \times m\) , \(U\) is a \(m \times m\) matrix with orthonormal columns, each column of U is the eigenvector \(v_{i}\) , and \(\Sigma\) is diagonal matrix, which contains the eigenvalues \(\lambda_{i}\) of matrix A.

The maximum matrix size supported in FPGA is templated by NMAX.

Parameters:

T data type (support float and double).
NMAX maximum number of rows/columns of input matrix
NCU number of computation unit
m number of rows/cols of matrix A
A input matrix of size \(m \times m\)
S decomposed diagonal singular matrix of size \(m \times m\)
U left U matrix of SVD
lda leading dimension of matrix A
ldu leading dimension of matrix U
info output info (unused)