Linear System Solver

LAPACK general linear system functions (gesv)

void blaze::gesv (int n, int nrhs, float *A, int lda, int *ipiv, float *B, int ldb, int *info)
 LAPACK kernel for solving a general single precision linear system of equations ( $ A*X=B $). More...
 
void blaze::gesv (int n, int nrhs, double *A, int lda, int *ipiv, double *B, int ldb, int *info)
 LAPACK kernel for solving a general double precision linear system of equations ( $ A*X=B $). More...
 
void blaze::gesv (int n, int nrhs, complex< float > *A, int lda, int *ipiv, complex< float > *B, int ldb, int *info)
 LAPACK kernel for solving a general single precision complex linear system of equations ( $ A*X=B $). More...
 
void blaze::gesv (int n, int nrhs, complex< double > *A, int lda, int *ipiv, complex< double > *B, int ldb, int *info)
 LAPACK kernel for solving a general double precision complex linear system of equations ( $ A*X=B $). More...
 
template<typename MT , bool SO, typename VT , bool TF>
void blaze::gesv (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &b, int *ipiv)
 LAPACK kernel for solving a general linear system of equations ( $ A*x=b $). More...
 
template<typename MT1 , bool SO1, typename MT2 , bool SO2>
void blaze::gesv (DenseMatrix< MT1, SO1 > &A, DenseMatrix< MT2, SO2 > &B, int *ipiv)
 LAPACK kernel for solving a general linear system of equations ( $ A*X=B $). More...
 

LAPACK Hermitian indefinite linear system functions (hesv)

void blaze::hesv (char uplo, int n, int nrhs, complex< float > *A, int lda, int *ipiv, complex< float > *B, int ldb, complex< float > *work, int lwork, int *info)
 LAPACK kernel for solving a Hermitian indefinite single precision complex linear system of equations ( $ A*X=B $). More...
 
void blaze::hesv (char uplo, int n, int nrhs, complex< double > *A, int lda, int *ipiv, complex< double > *B, int ldb, complex< double > *work, int lwork, int *info)
 LAPACK kernel for solving a Hermitian indefinite double precision complex linear system of equations ( $ A*X=B $). More...
 
template<typename MT , bool SO, typename VT , bool TF>
void blaze::hesv (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &b, char uplo, int *ipiv)
 LAPACK kernel for solving a Hermitian indefinite linear system of equations ( $ A*x=b $). More...
 
template<typename MT1 , bool SO1, typename MT2 , bool SO2>
void blaze::hesv (DenseMatrix< MT1, SO1 > &A, DenseMatrix< MT2, SO2 > &B, char uplo, int *ipiv)
 LAPACK kernel for solving a Hermitian indefinite linear system of equations ( $ A*X=B $). More...
 

LAPACK positive definite linear system functions (posv)

void blaze::posv (char uplo, int n, int nrhs, float *A, int lda, float *B, int ldb, int *info)
 LAPACK kernel for solving a positive definite single precision linear system of equations ( $ A*X=B $). More...
 
void blaze::posv (char uplo, int n, int nrhs, double *A, int lda, double *B, int ldb, int *info)
 LAPACK kernel for solving a positive definite double precision linear system of equations ( $ A*X=B $). More...
 
void blaze::posv (char uplo, int n, int nrhs, complex< float > *A, int lda, complex< float > *B, int ldb, int *info)
 LAPACK kernel for solving a positive definite single precision complex linear system of equations ( $ A*X=B $). More...
 
void blaze::posv (char uplo, int n, int nrhs, complex< double > *A, int lda, complex< double > *B, int ldb, int *info)
 LAPACK kernel for solving a positive definite double precision complex linear system of equations ( $ A*X=B $). More...
 
template<typename MT , bool SO, typename VT , bool TF>
void blaze::posv (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &b, char uplo)
 LAPACK kernel for solving a positive definite linear system of equations ( $ A*x=b $). More...
 
template<typename MT1 , bool SO1, typename MT2 , bool SO2>
void blaze::posv (DenseMatrix< MT1, SO1 > &A, DenseMatrix< MT2, SO2 > &B, char uplo)
 LAPACK kernel for solving a positive definite linear system of equations ( $ A*X=B $). More...
 

LAPACK symmetric indefinite linear system functions (sysv)

void blaze::sysv (char uplo, int n, int nrhs, float *A, int lda, int *ipiv, float *B, int ldb, float *work, int lwork, int *info)
 LAPACK kernel for solving a symmetric indefinite single precision linear system of equations ( $ A*X=B $). More...
 
void blaze::sysv (char uplo, int n, int nrhs, double *A, int lda, int *ipiv, double *B, int ldb, double *work, int lwork, int *info)
 LAPACK kernel for solving a symmetric indefinite double precision linear system of equations ( $ A*X=B $). More...
 
void blaze::sysv (char uplo, int n, int nrhs, complex< float > *A, int lda, int *ipiv, complex< float > *B, int ldb, complex< float > *work, int lwork, int *info)
 LAPACK kernel for solving a symmetric indefinite single precision complex linear system of equations ( $ A*X=B $). More...
 
void blaze::sysv (char uplo, int n, int nrhs, complex< double > *A, int lda, int *ipiv, complex< double > *B, int ldb, complex< double > *work, int lwork, int *info)
 LAPACK kernel for solving a symmetric indefinite double precision complex linear system of equations ( $ A*X=B $). More...
 
template<typename MT , bool SO, typename VT , bool TF>
void blaze::sysv (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &b, char uplo, int *ipiv)
 LAPACK kernel for solving a symmetric indefinite linear system of equations ( $ A*x=b $). More...
 
template<typename MT1 , bool SO1, typename MT2 , bool SO2>
void blaze::sysv (DenseMatrix< MT1, SO1 > &A, DenseMatrix< MT2, SO2 > &B, char uplo, int *ipiv)
 LAPACK kernel for solving a symmetric indefinite linear system of equations ( $ A*X=B $). More...
 

LAPACK triangular linear system functions (trsv)

void blaze::trsv (char uplo, char trans, char diag, int n, const float *A, int lda, float *x, int incX)
 LAPACK kernel for solving a triangular single precision linear system of equations ( $ A*X=B $). More...
 
void blaze::trsv (char uplo, char trans, char diag, int n, const double *A, int lda, double *x, int incX)
 LAPACK kernel for solving a triangular double precision linear system of equations ( $ A*x=b $). More...
 
void blaze::trsv (char uplo, char trans, char diag, int n, const complex< float > *A, int lda, complex< float > *x, int incX)
 LAPACK kernel for solving a triangular single precision complex linear system of equations ( $ A*x=b $). More...
 
void blaze::trsv (char uplo, char trans, char diag, int n, const complex< double > *A, int lda, complex< double > *x, int incX)
 LAPACK kernel for solving a triangular double precision complex linear system of equations ( $ A*x=b $). More...
 
template<typename MT , bool SO, typename VT , bool TF>
void blaze::trsv (const DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &b, char uplo, char trans, char diag)
 LAPACK kernel for solving a triangular linear system of equations ( $ A*x=b $). ( $ A*x=b $). More...
 

Detailed Description

Function Documentation

◆ gesv() [1/6]

template<typename MT , bool SO, typename VT , bool TF>
void blaze::gesv ( DenseMatrix< MT, SO > &  A,
DenseVector< VT, TF > &  b,
int *  ipiv 
)
inline

LAPACK kernel for solving a general linear system of equations ( $ A*x=b $).

Parameters
AThe column-major system matrix.
bThe right-hand side vector.
ipivAuxiliary array of size n for the pivot indices.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentInvalid right-hand side vector provided.
std::runtime_errorInversion of singular matrix failed.

This function uses the LAPACK gesv() functions to compute the solution to the system of general linear equations:

  • $ A *x=b $ if A is column-major
  • $ A^T*x=b $ if A is row-major

In this context the general system matrix A n-by-n matrix and x and b are n-dimensional vectors. Note that the function only works for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with adaptors or matrices of any other element type results in a compile time error!

If the function exits successfully, the vector b contains the solution of the linear system of equations and A has been decomposed by means of an LU decomposition with partial pivoting and row interchanges. The decomposition has the form

\[ A = P \cdot L \cdot U, \]

where P is a permutation matrix, L is a lower unitriangular matrix, and U is an upper triangular matrix. L is stored in the lower part of A and U is stored in the upper part. The unit diagonal elements of L are not stored. The factored form of A is then used to solve the system of equations.

The function fails if ...

  • ... the given system matrix is not a square matrix;
  • ... the given system matrix is singular and not invertible.

In all failure cases an exception is thrown.

For more information on the gesv() functions (i.e. sgesv(), dgesv(), cgesv(), and zgesv()), see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
This function does only provide the basic exception safety guarantee, i.e. in case of an exception A may already have been modified.

◆ gesv() [2/6]

template<typename MT1 , bool SO1, typename MT2 , bool SO2>
void blaze::gesv ( DenseMatrix< MT1, SO1 > &  A,
DenseMatrix< MT2, SO2 > &  B,
int *  ipiv 
)
inline

LAPACK kernel for solving a general linear system of equations ( $ A*X=B $).

Parameters
AThe column-major system matrix.
BThe matrix of right-hand sides.
ipivAuxiliary array of size n for the pivot indices.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentInvalid right-hand side matrix provided.
std::runtime_errorInversion of singular matrix failed.

This function uses the LAPACK gesv() functions to compute the solution to the general system of linera equations:

  • $ A *X =B $ if both A and B are column-major
  • $ A^T*X =B $ if A is row-major and B is column-major
  • $ A *X^T=B^T $ if A is column-major and B is row-major
  • $ A^T*X^T=B^T $ if both A and B are row-major

In this context the general system matrix A is a n-by-n matrix and X and B are either row-major m-by-n matrices or column-major n-by-m matrices. Note that the function only works for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with adaptors or matrices of any other element type results in a compile time error!

If the function exits successfully, the matrix B contains the solutions of the linear system of equations and A has been decomposed by means of an LU decomposition with partial pivoting and row interchanges. The decomposition has the form

\[ A = P \cdot L \cdot U, \]

where P is a permutation matrix, L is a lower unitriangular matrix, and U is an upper triangular matrix. L is stored in the lower part of A and U is stored in the upper part. The unit diagonal elements of L are not stored. The factored form of A is then used to solve the system of equations.

The function fails if ...

  • ... the given system matrix is not a square matrix;
  • ... the sizes of the two given matrices do not match;
  • ... the given system matrix is singular and not invertible.

In all failure cases an exception is thrown.

For more information on the gesv() functions (i.e. sgesv(), dgesv(), cgesv(), and zgesv()), see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
This function does only provide the basic exception safety guarantee, i.e. in case of an exception A may already have been modified.

◆ gesv() [3/6]

void blaze::gesv ( int  n,
int  nrhs,
float *  A,
int  lda,
int *  ipiv,
float *  B,
int  ldb,
int *  info 
)
inline

LAPACK kernel for solving a general single precision linear system of equations ( $ A*X=B $).

Parameters
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the single precision column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
infoReturn code of the function call.
Returns
void

This function uses the LAPACK sgesv() function to compute the solution to the general system of linear equations $ A*X=B $, where A is a n-by-n matrix and X and B are n-by-nrhs matrices.

The LU decomposition with partial pivoting and row interchanges is used to factor A as

\[ A = P \cdot L \cdot U, \]

where P is a permutation matrix, L is a lower unitriangular matrix, and U is an upper triangular matrix. The resulting decomposition is stored within A: L is stored in the lower part of A and U is stored in the upper part. The unit diagonal elements of L are not stored. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the decomposition has been completed, but since factor U(i,i) is exactly singular the solution could not be computed.

For more information on the sgesv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ gesv() [4/6]

void blaze::gesv ( int  n,
int  nrhs,
double *  A,
int  lda,
int *  ipiv,
double *  B,
int  ldb,
int *  info 
)
inline

LAPACK kernel for solving a general double precision linear system of equations ( $ A*X=B $).

Parameters
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the double precision column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
infoReturn code of the function call.
Returns
void

This function uses the LAPACK dgesv() function to compute the solution to the general system of linear equations $ A*X=B $, where A is a n-by-n matrix and X and B are n-by-nrhs matrices.

The LU decomposition with partial pivoting and row interchanges is used to factor A as

\[ A = P \cdot L \cdot U, \]

where P is a permutation matrix, L is a lower unitriangular matrix, and U is an upper triangular matrix. The resulting decomposition is stored within A: L is stored in the lower part of A and U is stored in the upper part. The unit diagonal elements of L are not stored. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the decomposition has been completed, but since factor U(i,i) is exactly singular the solution could not be computed.

For more information on the dgesv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ gesv() [5/6]

void blaze::gesv ( int  n,
int  nrhs,
complex< float > *  A,
int  lda,
int *  ipiv,
complex< float > *  B,
int  ldb,
int *  info 
)
inline

LAPACK kernel for solving a general single precision complex linear system of equations ( $ A*X=B $).

Parameters
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the single precision complex column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
infoReturn code of the function call.
Returns
void

This function uses the LAPACK cgesv() function to compute the solution to the general system of linear equations $ A*X=B $, where A is a n-by-n matrix and X and B are n-by-nrhs matrices.

The LU decomposition with partial pivoting and row interchanges is used to factor A as

\[ A = P \cdot L \cdot U, \]

where P is a permutation matrix, L is a lower unitriangular matrix, and U is an upper triangular matrix. The resulting decomposition is stored within A: L is stored in the lower part of A and U is stored in the upper part. The unit diagonal elements of L are not stored. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the decomposition has been completed, but since factor U(i,i) is exactly singular the solution could not be computed.

For more information on the cgesv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ gesv() [6/6]

void blaze::gesv ( int  n,
int  nrhs,
complex< double > *  A,
int  lda,
int *  ipiv,
complex< double > *  B,
int  ldb,
int *  info 
)
inline

LAPACK kernel for solving a general double precision complex linear system of equations ( $ A*X=B $).

Parameters
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the double precision complex column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
infoReturn code of the function call.
Returns
void

This function uses the LAPACK zgesv() function to compute the solution to the general system of linear equations $ A*X=B $, where A is a n-by-n matrix and X and B are n-by-nrhs matrices.

The LU decomposition with partial pivoting and row interchanges is used to factor A as

\[ A = P \cdot L \cdot U, \]

where P is a permutation matrix, L is a lower unitriangular matrix, and U is an upper triangular matrix. The resulting decomposition is stored within A: L is stored in the lower part of A and U is stored in the upper part. The unit diagonal elements of L are not stored. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the decomposition has been completed, but since factor U(i,i) is exactly singular the solution could not be computed.

For more information on the zgesv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ hesv() [1/4]

template<typename MT , bool SO, typename VT , bool TF>
void blaze::hesv ( DenseMatrix< MT, SO > &  A,
DenseVector< VT, TF > &  b,
char  uplo,
int *  ipiv 
)
inline

LAPACK kernel for solving a Hermitian indefinite linear system of equations ( $ A*x=b $).

Parameters
AThe column-major system matrix.
bThe right-hand side vector.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
ipivAuxiliary array of size n for the pivot indices.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentInvalid right-hand side vector provided.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorInversion of singular matrix failed.

This function uses the LAPACK hesv() functions to compute the solution to the Hermitian indefinite system of linear equations;

  • $ A *x=b $ if A is column-major
  • $ A^T*x=b $ if A is row-major

In this context the Hermitian indefinite system matrix A is a n-by-n matrix and x and b are n-dimensional vectors. Note that the function only works for general, non-adapted matrices with complex<float>, or complex<double> element type. The attempt to call the function with adaptors or matrices of any other element type results in a compile time error!

If the function exits successfully, the vector b contains the solution of the linear system of equations and A has been decomposed by means of the Bunch-Kaufman decomposition. The decomposition has the form

\[ A = U D U^{H} \texttt{ (if uplo = 'U'), or } A = L D L^{H} \texttt{ (if uplo = 'L'), } \]

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is Hermitian and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The function fails if ...

  • ... the given system matrix is not a square matrix;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the given system matrix is singular and not invertible.

In all failure cases an exception is thrown.

For more information on the hesv() functions (i.e. chesv() and zhesv()) see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
This function does only provide the basic exception safety guarantee, i.e. in case of an exception A may already have been modified.

◆ hesv() [2/4]

template<typename MT1 , bool SO1, typename MT2 , bool SO2>
void blaze::hesv ( DenseMatrix< MT1, SO1 > &  A,
DenseMatrix< MT2, SO2 > &  B,
char  uplo,
int *  ipiv 
)
inline

LAPACK kernel for solving a Hermitian indefinite linear system of equations ( $ A*X=B $).

Parameters
AThe system matrix.
BThe matrix of right-hand sides.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
ipivAuxiliary array of size n for the pivot indices.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentInvalid right-hand side matrix provided.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorInversion of singular matrix failed.

This function uses the LAPACK hesv() functions to compute the solution to the Hermitian indefinite system of linear equations:

  • $ A *X =B $ if both A and B are column-major
  • $ A^T*X =B $ if A is row-major and B is column-major
  • $ A *X^T=B^T $ if A is column-major and B is row-major
  • $ A^T*X^T=B^T $ if both A and B are row-major

In this context the Hermitian indefinite system matrix A is a n-by-n matrix and X and B are either row-major m-by-n matrices or column-major n-by-m matrices. Note that the function only works for general, non-adapted matrices with complex<float> or complex<double> element type. The attempt to call the function with adaptors or matrices of any other element type results in a compile time error!

If the function exits successfully, the matrix B contains the solution of the linear system of equations and A has been decomposed by means of a Bunch-Kaufman decomposition. The decomposition has the form

\[ A = U D U^{H} \texttt{ (if uplo = 'U'), or } A = L D L^{H} \texttt{ (if uplo = 'L'), } \]

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is Hermitian and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The function fails if ...

  • ... the given system matrix is not a square matrix;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the sizes of the two given matrices do not match;
  • ... the given system matrix is singular and not invertible.

In all failure cases an exception is thrown.

For more information on the hesv() functions (i.e. chesv() and zhesv()) see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
This function does only provide the basic exception safety guarantee, i.e. in case of an exception A may already have been modified.

◆ hesv() [3/4]

void blaze::hesv ( char  uplo,
int  n,
int  nrhs,
complex< float > *  A,
int  lda,
int *  ipiv,
complex< float > *  B,
int  ldb,
complex< float > *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for solving a Hermitian indefinite single precision complex linear system of equations ( $ A*X=B $).

Parameters
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the single precision complex column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; size >= max( 1, n ).
infoReturn code of the function call.
Returns
void

This function uses the LAPACK chesv() function to compute the solution to the Hermitian indefinite system of linear equations $ A*X=B $, where A is a n-by-n matrix and X and B are n-by-nrhs matrices.

The Bunch-Kaufman decomposition is used to factor A as

\[ A = U D U^{H} \texttt{ (if uplo = 'U'), or } A = L D L^{H} \texttt{ (if uplo = 'L'), } \]

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is Hermitian and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the decomposition has been completed, but since factor D(i,i) is exactly zero the solution could not be computed.

For more information on the chesv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ hesv() [4/4]

void blaze::hesv ( char  uplo,
int  n,
int  nrhs,
complex< double > *  A,
int  lda,
int *  ipiv,
complex< double > *  B,
int  ldb,
complex< double > *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for solving a Hermitian indefinite double precision complex linear system of equations ( $ A*X=B $).

Parameters
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the double precision complex column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; size >= max( 1, n ).
infoReturn code of the function call.
Returns
void

This function uses the LAPACK zhesv() function to compute the solution to the Hermitian indefinite system of linear equations $ A*X=B $, where A is a n-by-n matrix and X and B are n-by-nrhs matrices.

The Bunch-Kaufman decomposition is used to factor A as

\[ A = U D U^{H} \texttt{ (if uplo = 'U'), or } A = L D L^{H} \texttt{ (if uplo = 'L'), } \]

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is Hermitian and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the decomposition has been completed, but since factor D(i,i) is exactly zero the solution could not be computed.

For more information on the zhesv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ posv() [1/6]

template<typename MT , bool SO, typename VT , bool TF>
void blaze::posv ( DenseMatrix< MT, SO > &  A,
DenseVector< VT, TF > &  b,
char  uplo 
)
inline

LAPACK kernel for solving a positive definite linear system of equations ( $ A*x=b $).

Parameters
AThe column-major system matrix.
bThe right-hand side vector.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentInvalid right-hand side vector provided.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorInversion of singular matrix failed.

This function uses the LAPACK posv() functions to compute the solution to the positive definite system of linear equations:

  • $ A *x=b $ if A is column-major
  • $ A^T*x=b $ if A is row-major

In this context the positive definite system matrix A is a n-by-n matrix and x and b are n-dimensional vectors. Note that the function only works for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with adaptors or matrices of any other element type results in a compile time error!

If the function exits successfully, the vector b contains the solution of the linear system of equations and A has been decomposed by means of a Cholesky decomposition. The decomposition has the form

\[ A = U^{T} U \texttt{ (if uplo = 'U'), or } A = L L^{T} \texttt{ (if uplo = 'L'), } \]

where U is an upper triangular matrix and L is a lower triangular matrix. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The function fails if ...

  • ... the given system matrix is not a square matrix;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the given system matrix is singular and not invertible.

In all failure cases an exception is thrown.

For more information on the posv() functions (i.e. sposv(), dposv(), cposv(), and zposv()), see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
This function does only provide the basic exception safety guarantee, i.e. in case of an exception A may already have been modified.

◆ posv() [2/6]

template<typename MT1 , bool SO1, typename MT2 , bool SO2>
void blaze::posv ( DenseMatrix< MT1, SO1 > &  A,
DenseMatrix< MT2, SO2 > &  B,
char  uplo 
)
inline

LAPACK kernel for solving a positive definite linear system of equations ( $ A*X=B $).

Parameters
AThe system matrix.
BThe matrix of right-hand sides.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentInvalid right-hand side matrix provided.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorInversion of singular matrix failed.

This function uses the LAPACK posv() functions to compute the solution to the positive definite system of linear equations:

  • $ A *X =B $ if both A and B are column-major
  • $ A^T*X =B $ if A is row-major and B is column-major
  • $ A *X^T=B^T $ if A is column-major and B is row-major
  • $ A^T*X^T=B^T $ if both A and B are row-major

In this context the positive definite system matrix A is a n-by-n matrix and X and B are either row-major m-by-n matrices or column-major n-by-m matrices. Note that the function only works for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with adaptors or matrices of any other element type results in a compile time error!

If the function exits successfully, the matrix B contains the solution of the linear system of equations and A has been decomposed by means of a Cholesky decomposition. The decomposition has the form

\[ A = U^{T} U \texttt{ (if uplo = 'U'), or } A = L L^{T} \texttt{ (if uplo = 'L'), } \]

where U is an upper triangular matrix and L is a lower triangular matrix. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The function fails if ...

  • ... the given system matrix is not a square matrix;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the sizes of the two given matrices do not match;
  • ... the given system matrix is singular and not invertible.

In all failure cases an exception is thrown.

For more information on the posv() functions (i.e. sposv(), dposv(), cposv(), and zposv()), see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
This function does only provide the basic exception safety guarantee, i.e. in case of an exception A may already have been modified.

◆ posv() [3/6]

void blaze::posv ( char  uplo,
int  n,
int  nrhs,
float *  A,
int  lda,
float *  B,
int  ldb,
int *  info 
)
inline

LAPACK kernel for solving a positive definite single precision linear system of equations ( $ A*X=B $).

Parameters
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the single precision column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
infoReturn code of the function call.
Returns
void

This function uses the LAPACK sposv() function to compute the solution to the positive definite system of linear equations $ A*X=B $, where A is a n-by-n positive definite matrix and X and B are n-by-nrhs matrices.

The Cholesky decomposition is used to factor A as

\[ A = U^{T} U \texttt{ (if uplo = 'U'), or } A = L L^{T} \texttt{ (if uplo = 'L'), } \]

where U is an upper triangular matrix and L is a lower triangular matrix. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the leading minor of order i is not positive definite, so the decomposition could not be completed and the solution has not been computed.

For more information on the sposv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ posv() [4/6]

void blaze::posv ( char  uplo,
int  n,
int  nrhs,
double *  A,
int  lda,
double *  B,
int  ldb,
int *  info 
)
inline

LAPACK kernel for solving a positive definite double precision linear system of equations ( $ A*X=B $).

Parameters
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the double precision column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
infoReturn code of the function call.
Returns
void

This function uses the LAPACK dposv() function to compute the solution to the positive definite system of linear equations $ A*X=B $, where A is a n-by-n positive definite matrix and X and B are n-by-nrhs matrices.

The Cholesky decomposition is used to factor A as

\[ A = U^{T} U \texttt{ (if uplo = 'U'), or } A = L L^{T} \texttt{ (if uplo = 'L'), } \]

where U is an upper triangular matrix and L is a lower triangular matrix. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the leading minor of order i is not positive definite, so the decomposition could not be completed and the solution has not been computed.

For more information on the dposv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ posv() [5/6]

void blaze::posv ( char  uplo,
int  n,
int  nrhs,
complex< float > *  A,
int  lda,
complex< float > *  B,
int  ldb,
int *  info 
)
inline

LAPACK kernel for solving a positive definite single precision complex linear system of equations ( $ A*X=B $).

Parameters
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the single precision complex column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
infoReturn code of the function call.
Returns
void

This function uses the LAPACK cposv() function to compute the solution to the positive definite system of linear equations $ A*X=B $, where A is a n-by-n positive definite matrix and X and B are n-by-nrhs matrices.

The Cholesky decomposition is used to factor A as

\[ A = U^{T} U \texttt{ (if uplo = 'U'), or } A = L L^{T} \texttt{ (if uplo = 'L'), } \]

where U is an upper triangular matrix and L is a lower triangular matrix. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the leading minor of order i is not positive definite, so the decomposition could not be completed and the solution has not been computed.

For more information on the cposv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ posv() [6/6]

void blaze::posv ( char  uplo,
int  n,
int  nrhs,
complex< double > *  A,
int  lda,
complex< double > *  B,
int  ldb,
int *  info 
)
inline

LAPACK kernel for solving a positive definite double precision complex linear system of equations ( $ A*X=B $).

Parameters
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the double precision complex column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
infoReturn code of the function call.
Returns
void

This function uses the LAPACK zposv() function to compute the solution to the positive definite system of linear equations $ A*X=B $, where A is a n-by-n positive definite matrix and X and B are n-by-nrhs matrices.

The Cholesky decomposition is used to factor A as

\[ A = U^{T} U \texttt{ (if uplo = 'U'), or } A = L L^{T} \texttt{ (if uplo = 'L'), } \]

where U is an upper triangular matrix and L is a lower triangular matrix. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the leading minor of order i is not positive definite, so the decomposition could not be completed and the solution has not been computed.

For more information on the zposv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ sysv() [1/6]

template<typename MT , bool SO, typename VT , bool TF>
void blaze::sysv ( DenseMatrix< MT, SO > &  A,
DenseVector< VT, TF > &  b,
char  uplo,
int *  ipiv 
)
inline

LAPACK kernel for solving a symmetric indefinite linear system of equations ( $ A*x=b $).

Parameters
AThe column-major system matrix.
bThe right-hand side vector.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
ipivAuxiliary array of size n for the pivot indices.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentInvalid right-hand side vector provided.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorInversion of singular matrix failed.

This function uses the LAPACK sysv() functions to compute the solution to the symmetric indefinite system of linear equations;

  • $ A *x=b $ if A is column-major
  • $ A^T*x=b $ if A is row-major

In this context the symmetric indefinite system matrix A is a n-by-n matrix and x and b are n-dimensional vectors. Note that the function only works for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with adaptors or matrices of any other element type results in a compile time error!

If the function exits successfully, the vector b contains the solution of the linear system of equations and A has been decomposed by means of the Bunch-Kaufman decomposition. The decomposition has the form

\[ A = U D U^{T} \texttt{ (if uplo = 'U'), or } A = L D L^{T} \texttt{ (if uplo = 'L'), } \]

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is symmetric and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The function fails if ...

  • ... the given system matrix is not a square matrix;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the given system matrix is singular and not invertible.

In all failure cases an exception is thrown.

For more information on the sysv() functions (i.e. ssysv(), dsysv(), csysv(), and zsysv()), see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
This function does only provide the basic exception safety guarantee, i.e. in case of an exception A may already have been modified.

◆ sysv() [2/6]

template<typename MT1 , bool SO1, typename MT2 , bool SO2>
void blaze::sysv ( DenseMatrix< MT1, SO1 > &  A,
DenseMatrix< MT2, SO2 > &  B,
char  uplo,
int *  ipiv 
)
inline

LAPACK kernel for solving a symmetric indefinite linear system of equations ( $ A*X=B $).

Parameters
AThe system matrix.
BThe matrix of right-hand sides.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
ipivAuxiliary array of size n for the pivot indices.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentInvalid right-hand side matrix provided.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorInversion of singular matrix failed.

This function uses the LAPACK sysv() functions to compute the solution to the symmetric indefinite system of linear equations:

  • $ A *X =B $ if both A and B are column-major
  • $ A^T*X =B $ if A is row-major and B is column-major
  • $ A *X^T=B^T $ if A is column-major and B is row-major
  • $ A^T*X^T=B^T $ if both A and B are row-major

In this context the symmetric indefinite system matrix A is a n-by-n matrix and X and B are either row-major m-by-n matrices or column-major n-by-m matrices. Note that the function only works for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with adaptors or matrices of any other element type results in a compile time error!

If the function exits successfully, the matrix B contains the solution of the linear system of equations and A has been decomposed by means of a Bunch-Kaufman decomposition. The decomposition has the form

\[ A = U D U^{T} \texttt{ (if uplo = 'U'), or } A = L D L^{T} \texttt{ (if uplo = 'L'), } \]

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is symmetric and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The function fails if ...

  • ... the given system matrix is not a square matrix;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the sizes of the two given matrices do not match;
  • ... the given system matrix is singular and not invertible.

In all failure cases an exception is thrown.

For more information on the sysv() functions (i.e. ssysv(), dsysv(), csysv(), and zsysv()), see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
This function does only provide the basic exception safety guarantee, i.e. in case of an exception A may already have been modified.

◆ sysv() [3/6]

void blaze::sysv ( char  uplo,
int  n,
int  nrhs,
float *  A,
int  lda,
int *  ipiv,
float *  B,
int  ldb,
float *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for solving a symmetric indefinite single precision linear system of equations ( $ A*X=B $).

Parameters
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the single precision column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; size >= max( 1, n ).
infoReturn code of the function call.
Returns
void

This function uses the LAPACK ssysv() function to compute the solution to the symmetric indefinite system of linear equations $ A*X=B $, where A is a n-by-n matrix and X and B are n-by-nrhs matrices.

The Bunch-Kaufman decomposition is used to factor A as

\[ A = U D U^{T} \texttt{ (if uplo = 'U'), or } A = L D L^{T} \texttt{ (if uplo = 'L'), } \]

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is symmetric and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the decomposition has been completed, but since factor D(i,i) is exactly zero the solution could not be computed.

For more information on the ssysv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ sysv() [4/6]

void blaze::sysv ( char  uplo,
int  n,
int  nrhs,
double *  A,
int  lda,
int *  ipiv,
double *  B,
int  ldb,
double *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for solving a symmetric indefinite double precision linear system of equations ( $ A*X=B $).

Parameters
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the double precision column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; size >= max( 1, n ).
infoReturn code of the function call.
Returns
void

This function uses the LAPACK dsysv() function to compute the solution to the symmetric indefinite system of linear equations $ A*X=B $, where A is a n-by-n matrix and X and B are n-by-nrhs matrices.

The Bunch-Kaufman decomposition is used to factor A as

\[ A = U D U^{T} \texttt{ (if uplo = 'U'), or } A = L D L^{T} \texttt{ (if uplo = 'L'), } \]

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is symmetric and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the decomposition has been completed, but since factor D(i,i) is exactly zero the solution could not be computed.

For more information on the dsysv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ sysv() [5/6]

void blaze::sysv ( char  uplo,
int  n,
int  nrhs,
complex< float > *  A,
int  lda,
int *  ipiv,
complex< float > *  B,
int  ldb,
complex< float > *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for solving a symmetric indefinite single precision complex linear system of equations ( $ A*X=B $).

Parameters
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the single precision complex column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; size >= max( 1, n ).
infoReturn code of the function call.
Returns
void

This function uses the LAPACK csysv() function to compute the solution to the symmetric indefinite system of linear equations $ A*X=B $, where A is a n-by-n matrix and X and B are n-by-nrhs matrices.

The Bunch-Kaufman decomposition is used to factor A as

\[ A = U D U^{T} \texttt{ (if uplo = 'U'), or } A = L D L^{T} \texttt{ (if uplo = 'L'), } \]

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is symmetric and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the decomposition has been completed, but since factor D(i,i) is exactly zero the solution could not be computed.

For more information on the csysv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ sysv() [6/6]

void blaze::sysv ( char  uplo,
int  n,
int  nrhs,
complex< double > *  A,
int  lda,
int *  ipiv,
complex< double > *  B,
int  ldb,
complex< double > *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for solving a symmetric indefinite double precision complex linear system of equations ( $ A*X=B $).

Parameters
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows/columns of matrix A $[0..\infty)$.
nrhsThe number of right-hand side vectors $[0..\infty)$.
APointer to the first element of the double precision complex column-major square matrix.
ldaThe total number of elements between two columns of matrix A $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
BPointer to the first element of the column-major matrix.
ldbThe total number of elements between two columns of matrix B $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; size >= max( 1, n ).
infoReturn code of the function call.
Returns
void

This function uses the LAPACK zsysv() function to compute the solution to the symmetric indefinite system of linear equations $ A*X=B $, where A is a n-by-n matrix and X and B are n-by-nrhs matrices.

The Bunch-Kaufman decomposition is used to factor A as

\[ A = U D U^{T} \texttt{ (if uplo = 'U'), or } A = L D L^{T} \texttt{ (if uplo = 'L'), } \]

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is symmetric and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. The resulting decomposition is stored within A: In case uplo is set to 'L' the result is stored in the lower part of the matrix and the upper part remains untouched, in case uplo is set to 'U' the result is stored in the upper part and the lower part remains untouched. The factored form of A is then used to solve the system of equations.

The info argument provides feedback on the success of the function call:

  • = 0: The function finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the decomposition has been completed, but since factor D(i,i) is exactly zero the solution could not be computed.

For more information on the zsysv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.

◆ trsv() [1/5]

template<typename MT , bool SO, typename VT , bool TF>
void blaze::trsv ( const DenseMatrix< MT, SO > &  A,
DenseVector< VT, TF > &  b,
char  uplo,
char  trans,
char  diag 
)
inline

LAPACK kernel for solving a triangular linear system of equations ( $ A*x=b $). ( $ A*x=b $).

Parameters
AThe system matrix.
bThe right-hand side vector.
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
trans'N' for $ A*x=b $, 'T' for $ A^T*x=b $, and C for $ A^H*x=b $.
diag'U' in case of a unitriangular matrix, 'N' otherwise.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentInvalid right-hand side vector provided.
std::invalid_argumentInvalid uplo argument provided.
std::invalid_argumentInvalid trans argument provided.
std::invalid_argumentInvalid diag argument provided.

This function uses the LAPACK trsv() functions to compute the solution to the triangular system of linear equations:

  • $ A *x=b $ if A is column-major
  • $ A^T*x=b $ if A is row-major

In this context the positive definite system matrix A is a n-by-n matrix and x and b are n-dimensional vectors. Note that the function only works for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with adaptors or matrices of any other element type results in a compile time error!

If the function exits successfully, the vector x contains the solution of the linear system of equations. The function fails if ...

  • ... the given system matrix is not a square matrix;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the given trans argument is neither 'N' nor 'T' nor 'C';
  • ... the given diag argument is neither 'U' nor 'N'.

In all failure cases a std::invalid_argument exception is thrown.

Examples:

DynamicMatrix<double,columnMajor> A( 2UL, 2UL ); // The system matrix A
DynamicVector<double,columnVector> b( 2UL ); // The right-hand side vector b
// ... Initialization
DynamicMatrix<double,columnMajor> D( A ); // Temporary matrix to be decomposed
DynamicVector<double,columnVector> x( b ); // Temporary vector for the solution
trsv( D, x, 'L', 'N', 'N' );
assert( A * x == b );
DynamicMatrix<double,rowMajor> A( 2UL, 2UL ); // The system matrix A
DynamicVector<double,columnVector> b( 2UL ); // The right-hand side vector b
// ... Initialization
DynamicMatrix<double,rowMajor> D( A ); // Temporary matrix to be decomposed
DynamicVector<double,columnVector> x( b ); // Temporary vector for the solution
trsv( D, x, 'L', 'N', 'N' );
assert( trans( A ) * x == b );

For more information on the trsv() functions (i.e. strsv(), dtrsv(), ctrsv(), and ztrsv()), see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
The function does not perform any test for singularity or near-singularity. Such tests must be performed prior to calling this function!

◆ trsv() [2/5]

void blaze::trsv ( char  uplo,
char  trans,
char  diag,
int  n,
const float *  A,
int  lda,
float *  x,
int  incX 
)
inline

LAPACK kernel for solving a triangular single precision linear system of equations ( $ A*X=B $).

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
trans'N' for $ A*X=B $, 'T' for $ A^T*X=B $, and C for $ A^H*X=B $.
diag'U' in case of a unitriangular matrix, 'N' otherwise.
nThe number of rows/columns of the column-major triangular matrix $[0..\infty)$.
APointer to the first element of the single precision column-major square matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
xPointer to the first element of vector x.
incXThe stride within vector x.
Returns
void

This function uses the LAPACK strsv() function to compute the solution to the triangular system of linear equations $ A*x=b $, where A is a n-by-n triangular matrix and x and b are n-dimensional vectors. The trans argument specifies the form of the linear system of equations:

  • 'N': $ A*x=b $ (no transpose)
  • 'T': $ A^{T}*x=b $ (transpose)
  • 'C': $ A^{H}*x=b $ (conjugate transpose)

For more information on the strsv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
The function does not perform any test for singularity or near-singularity. Such tests must be performed prior to calling this function!

◆ trsv() [3/5]

void blaze::trsv ( char  uplo,
char  trans,
char  diag,
int  n,
const double *  A,
int  lda,
double *  x,
int  incX 
)
inline

LAPACK kernel for solving a triangular double precision linear system of equations ( $ A*x=b $).

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
trans'N' for $ A*X=B $, 'T' for $ A^T*X=B $, and C for $ A^H*X=B $.
diag'U' in case of a unitriangular matrix, 'N' otherwise.
nThe number of rows/columns of the column-major triangular matrix $[0..\infty)$.
APointer to the first element of the double precision column-major square matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
xPointer to the first element of vector x.
incXThe stride within vector x.
Returns
void

This function uses the LAPACK dtrsv() function to compute the solution to the triangular system of linear equations $ A*x=b $, where A is a n-by-n triangular matrix and x and b are n-dimensional vectors. The trans argument specifies the form of the linear system of equations:

  • 'N': $ A*x=b $ (no transpose)
  • 'T': $ A^{T}*x=b $ (transpose)
  • 'C': $ A^{H}*x=b $ (conjugate transpose)

For more information on the dtrsv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
The function does not perform any test for singularity or near-singularity. Such tests must be performed prior to calling this function!

◆ trsv() [4/5]

void blaze::trsv ( char  uplo,
char  trans,
char  diag,
int  n,
const complex< float > *  A,
int  lda,
complex< float > *  x,
int  incX 
)
inline

LAPACK kernel for solving a triangular single precision complex linear system of equations ( $ A*x=b $).

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
trans'N' for $ A*X=B $, 'T' for $ A^T*X=B $, and C for $ A^H*X=B $.
diag'U' in case of a unitriangular matrix, 'N' otherwise.
nThe number of rows/columns of the column-major triangular matrix $[0..\infty)$.
APointer to the first element of the single precision complex column-major square matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
xPointer to the first element of vector x.
incXThe stride within vector x.
Returns
void

This function uses the LAPACK ctrsv() function to compute the solution to the triangular system of linear equations $ A*x=b $, where A is a n-by-n triangular matrix and x and b are n-dimensional vectors. The trans argument specifies the form of the linear system of equations:

  • 'N': $ A*x=b $ (no transpose)
  • 'T': $ A^{T}*x=b $ (transpose)
  • 'C': $ A^{H}*x=b $ (conjugate transpose)

For more information on the ctrsv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
The function does not perform any test for singularity or near-singularity. Such tests must be performed prior to calling this function!

◆ trsv() [5/5]

void blaze::trsv ( char  uplo,
char  trans,
char  diag,
int  n,
const complex< double > *  A,
int  lda,
complex< double > *  x,
int  incX 
)
inline

LAPACK kernel for solving a triangular double precision complex linear system of equations ( $ A*x=b $).

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
trans'N' for $ A*X=B $, 'T' for $ A^T*X=B $, and C for $ A^H*X=B $.
diag'U' in case of a unitriangular matrix, 'N' otherwise.
nThe number of rows/columns of the column-major triangular matrix $[0..\infty)$.
APointer to the first element of the double precision complex column-major square matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
xPointer to the first element of vector x.
incXThe stride within vector x.
Returns
void

This function uses the LAPACK ztrsv() function to compute the solution to the triangular system of linear equations $ A*x=b $, where A is a n-by-n triangular matrix and x and b are n-dimensional vectors. The trans argument specifies the form of the linear system of equations:

  • 'N': $ A*x=b $ (no transpose)
  • 'T': $ A^{T}*x=b $ (transpose)
  • 'C': $ A^{H}*x=b $ (conjugate transpose)

For more information on the ztrsv() function, see the LAPACK online documentation browser:

http://www.netlib.org/lapack/explore-html/

Note
This function can only be used if a fitting LAPACK library, which supports this function, is available and linked to the executable. Otherwise a call to this function will result in a linker error.
The function does not perform any test for singularity or near-singularity. Such tests must be performed prior to calling this function!