Matrix Inversion

LAPACK LU-based inversion functions (getri)

void blaze::getri (int n, float *A, int lda, const int *ipiv, float *work, int lwork, int *info)
 LAPACK kernel for the inversion of the given dense general single precision column-major square matrix. More...
 
void blaze::getri (int n, double *A, int lda, const int *ipiv, double *work, int lwork, int *info)
 LAPACK kernel for the inversion of the given dense general double precision column-major square matrix. More...
 
void blaze::getri (int n, complex< float > *A, int lda, const int *ipiv, complex< float > *work, int lwork, int *info)
 LAPACK kernel for the inversion of the given dense general single precision complex column-major square matrix. More...
 
void blaze::getri (int n, complex< double > *A, int lda, const int *ipiv, complex< double > *work, int lwork, int *info)
 LAPACK kernel for the inversion of the given dense general double precision complex column-major square matrix. More...
 
template<typename MT , bool SO>
void blaze::getri (DenseMatrix< MT, SO > &A, const int *ipiv)
 LAPACK kernel for the inversion of the given dense general matrix. More...
 

LAPACK LDLH-based inversion functions (hetri)

void blaze::hetri (char uplo, int n, complex< float > *A, int lda, const int *ipiv, complex< float > *work, int *info)
 LAPACK kernel for the inversion of the given dense Hermitian indefinite single precision complex column-major square matrix. More...
 
void blaze::hetri (char uplo, int n, complex< double > *A, int lda, const int *ipiv, complex< double > *work, int *info)
 LAPACK kernel for the inversion of the given dense Hermitian indefinite double precision complex column-major square matrix. More...
 
template<typename MT , bool SO>
void blaze::hetri (DenseMatrix< MT, SO > &A, char uplo, const int *ipiv)
 LAPACK kernel for the inversion of the given dense Hermitian indefinite matrix. More...
 

LAPACK LLH-based inversion functions (potri)

void blaze::potri (char uplo, int n, float *A, int lda, int *info)
 LAPACK kernel for the inversion of the given dense positive definite single precision column-major square matrix. More...
 
void blaze::potri (char uplo, int n, double *A, int lda, int *info)
 LAPACK kernel for the inversion of the given dense positive definite double precision column-major square matrix. More...
 
void blaze::potri (char uplo, int n, complex< float > *A, int lda, int *info)
 LAPACK kernel for the inversion of the given dense positive definite single precision complex column-major square matrix. More...
 
void blaze::potri (char uplo, int n, complex< double > *A, int lda, int *info)
 LAPACK kernel for the inversion of the given dense positive definite double precision complex column-major square matrix. More...
 
template<typename MT , bool SO>
void blaze::potri (DenseMatrix< MT, SO > &A, char uplo)
 LAPACK kernel for the inversion of the given dense positive definite matrix. More...
 

LAPACK LDLT-based inversion functions (sytri)

void blaze::sytri (char uplo, int n, float *A, int lda, const int *ipiv, float *work, int *info)
 LAPACK kernel for the inversion of the given dense symmetric indefinite single precision column-major square matrix. More...
 
void blaze::sytri (char uplo, int n, double *A, int lda, const int *ipiv, double *work, int *info)
 LAPACK kernel for the inversion of the given dense symmetric indefinite double precision column-major square matrix. More...
 
void blaze::sytri (char uplo, int n, complex< float > *A, int lda, const int *ipiv, complex< float > *work, int *info)
 LAPACK kernel for the inversion of the given dense symmetric indefinite single precision complex column-major square matrix. More...
 
void blaze::sytri (char uplo, int n, complex< double > *A, int lda, const int *ipiv, complex< double > *work, int *info)
 LAPACK kernel for the inversion of the given dense symmetric indefinite double precision complex column-major square matrix. More...
 
template<typename MT , bool SO>
void blaze::sytri (DenseMatrix< MT, SO > &A, char uplo, const int *ipiv)
 LAPACK kernel for the inversion of the given dense symmetric indefinite matrix. More...
 

LAPACK triangular matrix inversion functions (trtri)

void blaze::trtri (char uplo, char diag, int n, float *A, int lda, int *info)
 LAPACK kernel for the inversion of the given dense triangular single precision column-major matrix. More...
 
void blaze::trtri (char uplo, char diag, int n, double *A, int lda, int *info)
 LAPACK kernel for the inversion of the given dense triangular double precision column-major matrix. More...
 
void blaze::trtri (char uplo, char diag, int n, complex< float > *A, int lda, int *info)
 LAPACK kernel for the inversion of the given dense triangular single precision complex column-major matrix. More...
 
void blaze::trtri (char uplo, char diag, int n, complex< double > *A, int lda, int *info)
 LAPACK kernel for the inversion of the given dense triangular double precision complex column-major matrix. More...
 
template<typename MT , bool SO>
void blaze::trtri (DenseMatrix< MT, SO > &A, char uplo, char diag)
 LAPACK kernel for the inversion of the given dense triangular matrix. More...
 

Detailed Description

Function Documentation

◆ getri() [1/5]

template<typename MT , bool SO>
void blaze::getri ( DenseMatrix< MT, SO > &  A,
const int *  ipiv 
)
inline

LAPACK kernel for the inversion of the given dense general matrix.

Parameters
AThe matrix to be inverted.
ipivAuxiliary array for the pivot indices; size >= min( m, n ).
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::runtime_errorInversion of singular matrix failed.

This function performs the dense matrix inversion based on the LAPACK getri() functions for matrices that have already been factorized by the getrf() functions. 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!

The function fails if ...

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

In all failure cases an exception is thrown.

For more information on the getri() functions (i.e. sgetri(), dgetri(), cgetri(), and zgetri()) 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.

◆ getri() [2/5]

void blaze::getri ( int  n,
float *  A,
int  lda,
const int *  ipiv,
float *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense general single precision column-major square matrix.

Parameters
nThe number of rows/columns of the 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)$.
ipivAuxiliary array for the pivot indices; size >= min( m, n ).
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 performs the dense matrix inversion based on the LAPACK sgetri() function for single precision column-major matrices that have already been factorized by the sgetrf() function. The info argument provides feedback on the success of the function call:

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the inversion could not be computed since U(i,i) is exactly zero.

If the function exits successfully (i.e. info = 0) then the first element of the work array returns the optimal lwork. For optimal performance lwork >= N*NB, where NB is the optimal blocksize returned by the LAPACK function ilaenv(). If lwork = -1 then a workspace query is assumed. The function only calculates the optimal size of the work array and returns this value as the first entry of the work array.

For more information on the sgetri() 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.

◆ getri() [3/5]

void blaze::getri ( int  n,
double *  A,
int  lda,
const int *  ipiv,
double *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense general double precision column-major square matrix.

Parameters
nThe number of rows/columns of the 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)$.
ipivAuxiliary array for the pivot indices; size >= min( m, n ).
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 performs the dense matrix inversion based on the LAPACK dgetri() function for double precision column-major matrices that have already been factorized by the dgetrf() function. The info argument provides feedback on the success of the function call:

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the inversion could not be computed since U(i,i) is exactly zero.

If the function exits successfully (i.e. info = 0) then the first element of the work array returns the optimal lwork. For optimal performance lwork >= N*NB, where NB is the optimal blocksize returned by the LAPACK function ilaenv(). If lwork = -1 then a workspace query is assumed. The function only calculates the optimal size of the work array and returns this value as the first entry of the work array.

For more information on the sgetri() 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.

◆ getri() [4/5]

void blaze::getri ( int  n,
complex< float > *  A,
int  lda,
const int *  ipiv,
complex< float > *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense general single precision complex column-major square matrix.

Parameters
nThe number of rows/columns of the 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)$.
ipivAuxiliary array for the pivot indices; size >= min( m, n ).
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 performs the dense matrix inversion based on the LAPACK cgetri() function for single precision complex column-major matrices that have already been factorized by the cgetrf() function. The info argument provides feedback on the success of the function call:

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the inversion could not be computed since U(i,i) is exactly zero.

If the function exits successfully (i.e. info = 0) then the first element of the work array returns the optimal lwork. For optimal performance lwork >= N*NB, where NB is the optimal blocksize returned by the LAPACK function ilaenv(). If lwork = -1 then a workspace query is assumed. The function only calculates the optimal size of the work array and returns this value as the first entry of the work array.

For more information on the sgetri() 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.

◆ getri() [5/5]

void blaze::getri ( int  n,
complex< double > *  A,
int  lda,
const int *  ipiv,
complex< double > *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense general double precision complex column-major square matrix.

Parameters
nThe number of rows/columns of the 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)$.
ipivAuxiliary array for the pivot indices; size >= min( m, n ).
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 performs the dense matrix inversion based on the LAPACK cgetri() function for double precision complex column-major matrices that have already been factorized by the zgetrf() function. The info argument provides feedback on the success of the function call:

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the inversion could not be computed since U(i,i) is exactly zero.

If the function exits successfully (i.e. info = 0) then the first element of the work array returns the optimal lwork. For optimal performance lwork >= N*NB, where NB is the optimal blocksize returned by the LAPACK function ilaenv(). If lwork = -1 then a workspace query is assumed. The function only calculates the optimal size of the work array and returns this value as the first entry of the work array.

For more information on the sgetri() 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.

◆ hetri() [1/3]

template<typename MT , bool SO>
void blaze::hetri ( DenseMatrix< MT, SO > &  A,
char  uplo,
const int *  ipiv 
)
inline

LAPACK kernel for the inversion of the given dense Hermitian indefinite matrix.

Parameters
AThe triangular matrix to be inverted.
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
ipivAuxiliary array of size n for the pivot indices.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorInversion of singular matrix failed.

This function performs the dense matrix inversion based on the LAPACK hetri() functions for Hermitian indefinite matrices that have already been factorized by the hetrf() functions. 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!

The function fails if ...

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

In all failure cases an exception is thrown.

For more information on the hetri() functions (i.e. chetri() and zhetri()) 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.

◆ hetri() [2/3]

void blaze::hetri ( char  uplo,
int  n,
complex< float > *  A,
int  lda,
const int *  ipiv,
complex< float > *  work,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense Hermitian indefinite single precision complex column-major square matrix.

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
nThe number of rows/columns of the Hermitian matrix $[0..\infty)$.
APointer to the first element of the single precision complex column-major matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
workAuxiliary array of size n.
infoReturn code of the function call.
Returns
void

This function performs the dense matrix inversion based on the LAPACK chetri() function for Hermitian indefinite single precision complex column-major matrices that have already been factorized by the chetrf() function.

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element D(i,i) is exactly zero and the inverse could not be computed.

For more information on the chetri() 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.

◆ hetri() [3/3]

void blaze::hetri ( char  uplo,
int  n,
complex< double > *  A,
int  lda,
const int *  ipiv,
complex< double > *  work,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense Hermitian indefinite double precision complex column-major square matrix.

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
nThe number of rows/columns of the Hermitian matrix $[0..\infty)$.
APointer to the first element of the double precision complex column-major matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
workAuxiliary array of size n.
infoReturn code of the function call.
Returns
void

This function performs the dense matrix inversion based on the LAPACK zhetri() function for Hermitian indefinite double precision complex column-major matrices that have already been factorized by the zhetrf() function.

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element D(i,i) is exactly zero and the inverse could not be computed.

For more information on the zhetri() 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.

◆ potri() [1/5]

template<typename MT , bool SO>
void blaze::potri ( DenseMatrix< MT, SO > &  A,
char  uplo 
)
inline

LAPACK kernel for the inversion of the given dense positive definite matrix.

Parameters
AThe positive-definite matrix to be inverted.
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 uplo argument provided.
std::runtime_errorInversion of singular matrix failed.

This function performs the dense matrix inversion based on the LAPACK potri() functions for positive-definite matrices that have already been factorized by the potrf() functions. The resulting symmetric inverse of the given matrix A is stored either in the lower part of A (uplo = 'L') or in the upper part (uplo = 'U'). 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!

The function fails if ...

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

In all failure cases an exception is thrown.

For more information on the potri() functions (i.e. spotri(), dpotri(), cpotri(), and zpotri()) 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.

◆ potri() [2/5]

void blaze::potri ( char  uplo,
int  n,
float *  A,
int  lda,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense positive definite single precision column-major square matrix.

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

This function performs the dense matrix inversion based on the LAPACK spotri() function for positive-definite single precision column-major matrices that have already been factorized by the spotrf() function. The resulting symmetric inverse of A is stored either in the lower part of A (uplo = 'L') or in the upper part (uplo = 'U').

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element (i,i) of U or L is zero and the inverse could not be computed.

For more information on the spotri() 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.

◆ potri() [3/5]

void blaze::potri ( char  uplo,
int  n,
double *  A,
int  lda,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense positive definite double precision column-major square matrix.

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

This function performs the dense matrix inversion based on the LAPACK dpotri() function for positive-definite double precision column-major matrices that have already been factorized by the dpotrf() function. The resulting symmetric inverse of A is stored either in the lower part of A (uplo = 'L') or in the upper part (uplo = 'U').

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element (i,i) of U or L is zero and the inverse could not be computed.

For more information on the spotri() 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.

◆ potri() [4/5]

void blaze::potri ( char  uplo,
int  n,
complex< float > *  A,
int  lda,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense positive definite single precision complex column-major square matrix.

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

This function performs the dense matrix inversion based on the LAPACK cpotri() function for positive-definite single precision complex column-major matrices that have already been factorized by the cpotrf() function. The resulting symmetric inverse of A is stored either in the lower part of A (uplo = 'L') or in the upper part (uplo = 'U').

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element (i,i) of U or L is zero and the inverse could not be computed.

For more information on the cpotri() 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.

◆ potri() [5/5]

void blaze::potri ( char  uplo,
int  n,
complex< double > *  A,
int  lda,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense positive definite double precision complex column-major square matrix.

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

This function performs the dense matrix inversion based on the LAPACK zpotri() function for positive-definite double precision complex column-major matrices that have already been factorized by the zpotrf() function. The resulting symmetric inverse of A is stored either in the lower part of A (uplo = 'L') or in the upper part (uplo = 'U').

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element (i,i) of U or L is zero and the inverse could not be computed.

For more information on the zpotri() 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.

◆ sytri() [1/5]

template<typename MT , bool SO>
void blaze::sytri ( DenseMatrix< MT, SO > &  A,
char  uplo,
const int *  ipiv 
)
inline

LAPACK kernel for the inversion of the given dense symmetric indefinite matrix.

Parameters
AThe triangular matrix to be inverted.
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
ipivAuxiliary array of size n for the pivot indices.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorInversion of singular matrix failed.

This function performs the dense matrix inversion based on the LAPACK sytri() functions for symmetric indefinite matrices that have already been factorized by the sytrf() functions. 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!

The function fails if ...

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

In all failure cases an exception is thrown.

For more information on the sytri() functions (i.e. ssytri(), dsytri(), csytri(), and zsytri()) 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.

◆ sytri() [2/5]

void blaze::sytri ( char  uplo,
int  n,
float *  A,
int  lda,
const int *  ipiv,
float *  work,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense symmetric indefinite single precision column-major square matrix.

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
nThe number of rows/columns of the symmetric matrix $[0..\infty)$.
APointer to the first element of the single precision column-major matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
workAuxiliary array of size n.
infoReturn code of the function call.
Returns
void

This function performs the dense matrix inversion based on the LAPACK ssytri() function for symmetric indefinite single precision column-major matrices that have already been factorized by the ssytrf() function.

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element D(i,i) is exactly zero and the inverse could not be computed.

For more information on the ssytri() 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.

◆ sytri() [3/5]

void blaze::sytri ( char  uplo,
int  n,
double *  A,
int  lda,
const int *  ipiv,
double *  work,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense symmetric indefinite double precision column-major square matrix.

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
nThe number of rows/columns of the symmetric matrix $[0..\infty)$.
APointer to the first element of the double precision column-major matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
workAuxiliary array of size n.
infoReturn code of the function call.
Returns
void

This function performs the dense matrix inversion based on the LAPACK dsytri() function for symmetric indefinite double precision column-major matrices that have already been factorized by the dsytrf() function.

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element D(i,i) is exactly zero and the inverse could not be computed.

For more information on the dsytri() 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.

◆ sytri() [4/5]

void blaze::sytri ( char  uplo,
int  n,
complex< float > *  A,
int  lda,
const int *  ipiv,
complex< float > *  work,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense symmetric indefinite single precision complex column-major square matrix.

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
nThe number of rows/columns of the symmetric matrix $[0..\infty)$.
APointer to the first element of the single precision complex column-major matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
workAuxiliary array of size n.
infoReturn code of the function call.
Returns
void

This function performs the dense matrix inversion based on the LAPACK csytri() function for symmetric indefinite single precision complex column-major matrices that have already been factorized by the csytrf() function.

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element D(i,i) is exactly zero and the inverse could not be computed.

For more information on the csytri() 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.

◆ sytri() [5/5]

void blaze::sytri ( char  uplo,
int  n,
complex< double > *  A,
int  lda,
const int *  ipiv,
complex< double > *  work,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense symmetric indefinite double precision complex column-major square matrix.

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
nThe number of rows/columns of the symmetric matrix $[0..\infty)$.
APointer to the first element of the double precision complex column-major matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
ipivAuxiliary array of size n for the pivot indices.
workAuxiliary array of size n.
infoReturn code of the function call.
Returns
void

This function performs the dense matrix inversion based on the LAPACK zsytri() function for symmetric indefinite double precision complex column-major matrices that have already been factorized by the zsytrf() function.

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element D(i,i) is exactly zero and the inverse could not be computed.

For more information on the zsytri() 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.

◆ trtri() [1/5]

template<typename MT , bool SO>
void blaze::trtri ( DenseMatrix< MT, SO > &  A,
char  uplo,
char  diag 
)
inline

LAPACK kernel for the inversion of the given dense triangular matrix.

Parameters
AThe triangular matrix to be inverted.
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
diag'U' in case of a unitriangular matrix, 'N' otherwise.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentInvalid uplo argument provided.
std::invalid_argumentInvalid diag argument provided.
std::runtime_errorInversion of singular matrix failed.

This function performs the dense matrix inversion based on the LAPACK trtri() functions for a lower triangular (uplo = 'L') or upper triangular (uplo = 'U') matrix. 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!

The function fails if ...

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

In all failure cases an exception is thrown.

For more information on the trtri() functions (i.e. strtri(), dtrtri(), ctrtri(), and ztrtri()) 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.

◆ trtri() [2/5]

void blaze::trtri ( char  uplo,
char  diag,
int  n,
float *  A,
int  lda,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense triangular single precision column-major matrix.

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
diag'U' in case of a unitriangular matrix, 'N' otherwise.
nThe number of rows/columns of the triangular matrix $[0..\infty)$.
APointer to the first element of the single precision column-major matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
infoReturn code of the function call.
Returns
void

This function performs the dense matrix inversion based on the LAPACK strtri() function for lower triangular (uplo = 'L') or upper triangular (uplo = 'U') single precision column-major matrices.

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element A(i,i) is exactly zero and the inverse could not be computed.

For more information on the strtri() 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.

◆ trtri() [3/5]

void blaze::trtri ( char  uplo,
char  diag,
int  n,
double *  A,
int  lda,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense triangular double precision column-major matrix.

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
diag'U' in case of a unitriangular matrix, 'N' otherwise.
nThe number of rows/columns of the triangular matrix $[0..\infty)$.
APointer to the first element of the double precision column-major matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
infoReturn code of the function call.
Returns
void

This function performs the dense matrix inversion based on the LAPACK dtrtri() function for lower triangular (uplo = 'L') or upper triangular (uplo = 'U') double precision column-major matrices.

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element A(i,i) is exactly zero and the inverse could not be computed.

For more information on the dtrtri() 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.

◆ trtri() [4/5]

void blaze::trtri ( char  uplo,
char  diag,
int  n,
complex< float > *  A,
int  lda,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense triangular single precision complex column-major matrix.

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
diag'U' in case of a unitriangular matrix, 'N' otherwise.
nThe number of rows/columns of the triangular matrix $[0..\infty)$.
APointer to the first element of the single precision complex column-major matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
infoReturn code of the function call.
Returns
void

This function performs the dense matrix inversion based on the LAPACK ctrtri() function for lower triangular (uplo = 'L') or upper triangular (uplo = 'U') single precision complex column-major matrices.

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element A(i,i) is exactly zero and the inverse could not be computed.

For more information on the ctrtri() 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.

◆ trtri() [5/5]

void blaze::trtri ( char  uplo,
char  diag,
int  n,
complex< double > *  A,
int  lda,
int *  info 
)
inline

LAPACK kernel for the inversion of the given dense triangular double precision complex column-major matrix.

Parameters
uplo'L' in case of a lower matrix, 'U' in case of an upper matrix.
diag'U' in case of a unitriangular matrix, 'N' otherwise.
nThe number of rows/columns of the triangular matrix $[0..\infty)$.
APointer to the first element of the double precision complex column-major matrix.
ldaThe total number of elements between two columns of the matrix $[0..\infty)$.
infoReturn code of the function call.
Returns
void

This function performs the dense matrix inversion based on the LAPACK ztrtri() function for lower triangular (uplo = 'L') or upper triangular (uplo = 'U') double precision complex column-major matrices.

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

  • = 0: The inversion finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, element A(i,i) is exactly zero and the inverse could not be computed.

For more information on the ztrtri() 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.