Blaze  3.6
Eigenvalue

LAPACK general matrix eigenvalue functions (geev)

void blaze::geev (char jobvl, char jobvr, int n, float *A, int lda, float *wr, float *wi, float *VL, int ldvl, float *VR, int ldvr, float *work, int lwork, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense general single precision column-major matrix. More...
 
void blaze::geev (char jobvl, char jobvr, int n, double *A, int lda, double *wr, double *wi, double *VL, int ldvl, double *VR, int ldvr, double *work, int lwork, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense general double precision column-major matrix. More...
 
void blaze::geev (char jobvl, char jobvr, int n, complex< float > *A, int lda, complex< float > *w, complex< float > *VL, int ldvl, complex< float > *VR, int ldvr, complex< float > *work, int lwork, float *rwork, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense general single precision complex column-major matrix. More...
 
void blaze::geev (char jobvl, char jobvr, int n, complex< double > *A, int lda, complex< double > *w, complex< double > *VL, int ldvl, complex< double > *VR, int ldvr, complex< double > *work, int lwork, double *rwork, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense general double precision complex column-major matrix. More...
 
template<typename MT , bool SO, typename VT , bool TF>
void blaze::geev (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &w)
 LAPACK kernel for computing the eigenvalues of the given dense general matrix. More...
 
template<typename MT1 , bool SO1, typename MT2 , bool SO2, typename VT , bool TF>
void blaze::geev (DenseMatrix< MT1, SO1 > &A, DenseMatrix< MT2, SO2 > &VL, DenseVector< VT, TF > &w)
 LAPACK kernel for computing the eigenvalues of the given dense general matrix. More...
 
template<typename MT1 , bool SO1, typename VT , bool TF, typename MT2 , bool SO2>
void blaze::geev (DenseMatrix< MT1, SO1 > &A, DenseVector< VT, TF > &w, DenseMatrix< MT2, SO2 > &VR)
 LAPACK kernel for computing the eigenvalues of the given dense general matrix. More...
 
template<typename MT1 , bool SO1, typename MT2 , bool SO2, typename VT , bool TF, typename MT3 , bool SO3>
void blaze::geev (DenseMatrix< MT1, SO1 > &A, DenseMatrix< MT2, SO2 > &VL, DenseVector< VT, TF > &w, DenseMatrix< MT3, SO3 > &VR)
 LAPACK kernel for computing the eigenvalues of the given dense general matrix. More...
 

LAPACK Hermitian matrix eigenvalue functions (heev)

void blaze::heev (char jobz, char uplo, int n, complex< float > *A, int lda, float *w, complex< float > *work, int lwork, float *rwork, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense Hermitian single precision column-major matrix. More...
 
void blaze::heev (char jobz, char uplo, int n, complex< double > *A, int lda, double *w, complex< double > *work, int lwork, double *rwork, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense Hermitian double precision column-major matrix. More...
 
template<typename MT , bool SO, typename VT , bool TF>
void blaze::heev (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &w, char jobz, char uplo)
 LAPACK kernel for computing the eigenvalues of the given dense Hermitian matrix. More...
 

LAPACK Hermitian matrix eigenvalue functions (heevd)

void blaze::heevd (char jobz, char uplo, int n, complex< float > *A, int lda, float *w, complex< float > *work, int lwork, float *rwork, int lrwork, int *iwork, int liwork, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense Hermitian single precision column-major matrix. More...
 
void blaze::heevd (char jobz, char uplo, int n, complex< double > *A, int lda, double *w, complex< double > *work, int lwork, double *rwork, int lrwork, int *iwork, int liwork, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense Hermitian double precision column-major matrix. More...
 
template<typename MT , bool SO, typename VT , bool TF>
void blaze::heevd (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &w, char jobz, char uplo)
 LAPACK kernel for computing the eigenvalues of the given dense Hermitian matrix. More...
 

LAPACK Hermitian matrix eigenvalue functions (heevx)

void blaze::heevx (char jobz, char range, char uplo, int n, complex< float > *A, int lda, float vl, float vu, int il, int iu, float abstol, int *m, float *w, complex< float > *Z, int ldz, complex< float > *work, int lwork, float *rwork, int *iwork, int *ifail, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense Hermitian single precision complex column-major matrix. More...
 
void blaze::heevx (char jobz, char range, char uplo, int n, complex< double > *A, int lda, double vl, double vu, int il, int iu, double abstol, int *m, double *w, complex< double > *Z, int ldz, complex< double > *work, int lwork, double *rwork, int *iwork, int *ifail, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense Hermitian double precision complex column-major matrix. More...
 
template<typename MT , bool SO, typename VT , bool TF>
size_t blaze::heevx (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &w, char uplo)
 LAPACK kernel for the eigenvalue computation of the given dense Hermitian matrix. More...
 
template<typename MT , bool SO, typename VT , bool TF, typename ST >
size_t blaze::heevx (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &w, char uplo, ST low, ST upp)
 LAPACK kernel for computing the eigenvalues of the given dense Hermitian matrix. More...
 
template<typename MT1 , bool SO1, typename VT , bool TF, typename MT2 , bool SO2>
size_t blaze::heevx (DenseMatrix< MT1, SO1 > &A, DenseVector< VT, TF > &w, DenseMatrix< MT2, SO2 > &Z, char uplo)
 LAPACK kernel for computing the eigenvalues of the given dense Hermitian matrix. More...
 
template<typename MT1 , bool SO1, typename VT , bool TF, typename MT2 , bool SO2, typename ST >
size_t blaze::heevx (DenseMatrix< MT1, SO1 > &A, DenseVector< VT, TF > &w, DenseMatrix< MT2, SO2 > &Z, char uplo, ST low, ST upp)
 LAPACK kernel for computing the eigenvalues of the given dense Hermitian matrix. More...
 

LAPACK symmetric matrix eigenvalue functions (syev)

void blaze::syev (char jobz, char uplo, int n, float *A, int lda, float *w, float *work, int lwork, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense symmetric single precision column-major matrix. More...
 
void blaze::syev (char jobz, char uplo, int n, double *A, int lda, double *w, double *work, int lwork, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense symmetric double precision column-major matrix. More...
 
template<typename MT , bool SO, typename VT , bool TF>
void blaze::syev (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &w, char jobz, char uplo)
 LAPACK kernel for computing the eigenvalues of the given dense symmetric matrix. More...
 

LAPACK symmetric matrix eigenvalue functions (syevd)

void blaze::syevd (char jobz, char uplo, int n, float *A, int lda, float *w, float *work, int lwork, int *iwork, int liwork, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense symmetric single precision column-major matrix. More...
 
void blaze::syevd (char jobz, char uplo, int n, double *A, int lda, double *w, double *work, int lwork, int *iwork, int liwork, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense symmetric double precision column-major matrix. More...
 
template<typename MT , bool SO, typename VT , bool TF>
void blaze::syevd (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &w, char jobz, char uplo)
 LAPACK kernel for computing the eigenvalues of the given dense symmetric matrix. More...
 

LAPACK symmetric matrix eigenvalue functions (syevx)

void blaze::syevx (char jobz, char range, char uplo, int n, float *A, int lda, float vl, float vu, int il, int iu, float abstol, int *m, float *w, float *Z, int ldz, float *work, int lwork, int *iwork, int *ifail, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense symmetric single precision column-major matrix. More...
 
void blaze::syevx (char jobz, char range, char uplo, int n, double *A, int lda, double vl, double vu, int il, int iu, double abstol, int *m, double *w, double *Z, int ldz, double *work, int lwork, int *iwork, int *ifail, int *info)
 LAPACK kernel for computing the eigenvalues of the given dense symmetric double precision column-major matrix. More...
 
template<typename MT , bool SO, typename VT , bool TF>
size_t blaze::syevx (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &w, char uplo)
 LAPACK kernel for computing the eigenvalues of the given dense symmetric matrix. More...
 
template<typename MT , bool SO, typename VT , bool TF, typename ST >
size_t blaze::syevx (DenseMatrix< MT, SO > &A, DenseVector< VT, TF > &w, char uplo, ST low, ST upp)
 LAPACK kernel for computing the eigenvalues of the given dense symmetric matrix. More...
 
template<typename MT1 , bool SO1, typename VT , bool TF, typename MT2 , bool SO2>
size_t blaze::syevx (DenseMatrix< MT1, SO1 > &A, DenseVector< VT, TF > &w, DenseMatrix< MT2, SO2 > &Z, char uplo)
 LAPACK kernel for computing the eigenvalues of the given dense symmetric matrix. More...
 
template<typename MT1 , bool SO1, typename VT , bool TF, typename MT2 , bool SO2, typename ST >
size_t blaze::syevx (DenseMatrix< MT1, SO1 > &A, DenseVector< VT, TF > &w, DenseMatrix< MT2, SO2 > &Z, char uplo, ST low, ST upp)
 LAPACK kernel for computing the eigenvalues of the given dense symmetric matrix. More...
 

Detailed Description

Function Documentation

◆ geev() [1/8]

template<typename MT , bool SO, typename VT , bool TF>
void blaze::geev ( DenseMatrix< MT, SO > &  A,
DenseVector< VT, TF > &  w 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense general matrix.

Parameters
AThe given general matrix.
wThe resulting vector of complex eigenvalues.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::runtime_errorEigenvalue computation failed.

This function computes the eigenvalues of a non-symmetric, non-Hermitian n-by-n matrix based on the LAPACK geev() functions. The complex eigenvalues are returned in the given vector w, which is resized to the correct size (if possible and necessary). Please note that no order of eigenvalues can be assumed, except that complex conjugate pairs of eigenvalues appear consecutively with the eigenvalue having the positive imaginary part first.

Note that this function can only be used for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with any adapted matrix or matrices of any other element type results in a compile time error!

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the eigenvalue computation fails.

In all failure cases an exception is thrown.

Examples:

DynamicMatrix<double,rowMajor> A( 5UL, 5UL ); // The general matrix A
// ... Initialization
DynamicVector<complex<double>,columnVector> w( 5UL ); // The vector for the complex eigenvalues
geev( A, w );

For more information on the geev() functions (i.e. sgeev(), dgeev(), cgeev(), and zgeev()) 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.

◆ geev() [2/8]

template<typename MT1 , bool SO1, typename MT2 , bool SO2, typename VT , bool TF>
void blaze::geev ( DenseMatrix< MT1, SO1 > &  A,
DenseMatrix< MT2, SO2 > &  VL,
DenseVector< VT, TF > &  w 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense general matrix.

Parameters
AThe given general matrix.
VLThe resulting matrix of left eigenvectors.
wThe resulting vector of complex eigenvalues.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::invalid_argumentMatrix cannot be resized.
std::runtime_errorEigenvalue computation failed.

This function computes the eigenvalues of a non-symmetric, non-Hermitian n-by-n matrix based on the LAPACK geev() functions. Additionally, it computes the left eigenvectors. The left eigenvector $u[j]$ of A satisfies

\[ u[j]^{H} * A = lambda[j] * u[j]^{H}, \]

where $u[j]^{H}$ denotes the conjugate transpose of $u[j]$.

The complex eigenvalues are returned in the given vector w. The left eigenvectors are returned in the rows of VL in case VL is a row-major matrix and in the columns of VL in case VL is a column-major matrix. Both w and VL are resized to the correct dimensions (if possible and necessary). Please note that no order of eigenvalues can be assumed, except that complex conjugate pairs of eigenvalues appear consecutively with the eigenvalue having the positive imaginary part first. The computed eigenvectors are normalized to have Euclidean norm equal to 1 and largest component real.

Note that this function can only be used for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with any adapted matrix or matrices of any other element type results in a compile time error!

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given matrix VL is a fixed size matrix and the dimensions don't match;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the eigenvalue computation fails.

In all failure cases an exception is thrown.

Examples:

DynamicMatrix<double,rowMajor> A( 5UL, 5UL ); // The general matrix A
// ... Initialization
DynamicMatrix<complex<double>,rowMajor> VL( 5UL, 5UL ); // The matrix for the left eigenvectors
DynamicVector<complex<double>,columnVector> w( 5UL ); // The vector for the complex eigenvalues
geev( A, VL, w );

For more information on the geev() functions (i.e. sgeev(), dgeev(), cgeev(), and zgeev()) 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.

◆ geev() [3/8]

template<typename MT1 , bool SO1, typename VT , bool TF, typename MT2 , bool SO2>
void blaze::geev ( DenseMatrix< MT1, SO1 > &  A,
DenseVector< VT, TF > &  w,
DenseMatrix< MT2, SO2 > &  VR 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense general matrix.

Parameters
AThe given general matrix.
wThe resulting vector of complex eigenvalues.
VRThe resulting matrix of right eigenvectors.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::invalid_argumentMatrix cannot be resized.
std::runtime_errorEigenvalue computation failed.

This function computes the eigenvalues of a non-symmetric, non-Hermitian n-by-n matrix based on the LAPACK geev() functions. Additionally, it computes the right eigenvectors. The right eigenvector $v[j]$ of A satisfies

\[ A * v[j] = lambda[j] * v[j], \]

where $lambda[j]$ is its eigenvalue.

The complex eigenvalues are returned in the given vector w. The right eigenvectors are returned in the rows of VR in case VR is a row-major matrix and in the columns of VR in case VR is a column-major matrix. Both w and VL are resized to the correct dimensions (if possible and necessary). Please note that no order of eigenvalues can be assumed, except that complex conjugate pairs of eigenvalues appear consecutively with the eigenvalue having the positive imaginary part first. The computed eigenvectors are normalized to have Euclidean norm equal to 1 and largest component real.

Note that this function can only be used for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with any adapted matrix or matrices of any other element type results in a compile time error!

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given matrix VR is a fixed size matrix and the dimensions don't match;
  • ... the eigenvalue computation fails.

In all failure cases an exception is thrown.

Examples:

DynamicMatrix<double,rowMajor> A( 5UL, 5UL ); // The general matrix A
// ... Initialization
DynamicVector<complex<double>,columnVector> w( 5UL ); // The vector for the complex eigenvalues
DynamicMatrix<complex<double>,rowMajor> VR( 5UL, 5UL ); // The matrix for the right eigenvectors
geev( A, w, VR );

For more information on the geev() functions (i.e. sgeev(), dgeev(), cgeev(), and zgeev()) 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.

◆ geev() [4/8]

template<typename MT1 , bool SO1, typename MT2 , bool SO2, typename VT , bool TF, typename MT3 , bool SO3>
void blaze::geev ( DenseMatrix< MT1, SO1 > &  A,
DenseMatrix< MT2, SO2 > &  VL,
DenseVector< VT, TF > &  w,
DenseMatrix< MT3, SO3 > &  VR 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense general matrix.

Parameters
AThe given general matrix.
VLThe resulting matrix of left eigenvectors.
wThe resulting vector of complex eigenvalues.
VRThe resulting matrix of right eigenvectors.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::invalid_argumentMatrix cannot be resized.
std::runtime_errorEigenvalue computation failed.

This function computes the eigenvalues of a non-symmetric, non-Hermitian n-by-n matrix based on the LAPACK geev() functions. Additionally, it computes the left and right eigenvectors. The right eigenvector $v[j]$ of A satisfies

\[ A * v[j] = lambda[j] * v[j], \]

where $lambda[j]$ is its eigenvalue. The left eigenvector $u[j]$ of A satisfies

\[ u[j]^{H} * A = lambda[j] * u[j]^{H}, \]

where $u[j]^{H}$ denotes the conjugate transpose of $u[j]$.

The complex eigenvalues are returned in the given vector w. The left eigenvectors are returned in the rows of VL in case VL is a row-major matrix and in the columns of VL in case VL is a column-major matrix. The right eigenvectors are returned in the rows of VR in case VR is a row-major matrix and in the columns of VR in case VR is a column-major matrix. w, VL, and VR are resized to the correct dimensions (if possible and necessary). Please note that no order of eigenvalues can be assumed, except that complex conjugate pairs of eigenvalues appear consecutively with the eigenvalue having the positive imaginary part first. The computed eigenvectors are normalized to have Euclidean norm equal to 1 and largest component real.

Note that this function can only be used for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with any adapted matrix or matrices of any other element type results in a compile time error!

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given matrix VL is a fixed size matrix and the dimensions don't match;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given matrix VR is a fixed size matrix and the dimensions don't match;
  • ... the eigenvalue computation fails.

In all failure cases an exception is thrown.

Examples:

DynamicMatrix<double,rowMajor> A( 5UL, 5UL ); // The general matrix A
// ... Initialization
DynamicMatrix<complex<double>,rowMajor> VL( 5UL, 5UL ); // The matrix for the left eigenvectors
DynamicVector<complex<double>,columnVector> w( 5UL ); // The vector for the complex eigenvalues
DynamicMatrix<complex<double>,rowMajor> VR( 5UL, 5UL ); // The matrix for the right eigenvectors
geev( A, VL, w, VR );

For more information on the geev() functions (i.e. sgeev(), dgeev(), cgeev(), and zgeev()) 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.

◆ geev() [5/8]

void blaze::geev ( char  jobvl,
char  jobvr,
int  n,
float *  A,
int  lda,
float *  wr,
float *  wi,
float *  VL,
int  ldvl,
float *  VR,
int  ldvr,
float *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense general single precision column-major matrix.

Parameters
jobvl'V' to compute the left eigenvectors of A, 'N' to not compute them.
jobvr'V' to compute the right eigenvectors of A, 'N' to not compute them.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
wrPointer to the first element of the vector for the real part of the eigenvalues.
wiPointer to the first element of the vector for the imaginary part of the eigenvalues.
VLPointer to the first element of the column-major matrix for the left eigenvectors.
ldvlThe total number of elements between two columns of the matrix VL $[0..\infty)$.
VRPointer to the first element of the column-major matrix for the right eigenvectors.
ldvrThe total number of elements between two columns of the matrix VR $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
infoReturn code of the function call.
Returns
void

This function computes the eigenvalues of a non-symmetric n-by-n single precision column-major matrix based on the LAPACK sgeev() function. Optionally, it computes the left and/or right eigenvectors. The right eigenvector $v[j]$ of satisfies

\[ A * v[j] = lambda[j] * v[j], \]

where $lambda[j]$ is its eigenvalue. The left eigenvector $u[j]$ of A satisfies

\[ u[j]^{H} * A = lambda[j] * u[j]^{H}, \]

where $u[j]^{H}$ denotes the conjugate transpose of $u[j]$.

Complex conjugate pairs of eigenvalues appear consecutively with the eigenvalue having the positive imaginary part first. The computed eigenvectors are normalized to have Euclidean norm equal to 1 and largest component real.

The parameter jobvl specifies the computation of the left eigenvectors:

  • 'V': The left eigenvectors of A are computed and returned in VL.
  • 'N': The left eigenvectors of A are not computed.

The parameter jobvr specifies the computation of the right eigenvectors:

  • 'V': The right eigenvectors of A are computed and returned in VR.
  • 'N': The right eigenvectors of A are not computed.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the QR algorithm failed to compute all the eigenvalues and no eigenvectors have been computed; the elements with index larger than i have converged.

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

◆ geev() [6/8]

void blaze::geev ( char  jobvl,
char  jobvr,
int  n,
double *  A,
int  lda,
double *  wr,
double *  wi,
double *  VL,
int  ldvl,
double *  VR,
int  ldvr,
double *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense general double precision column-major matrix.

Parameters
jobvl'V' to compute the left eigenvectors of A, 'N' to not compute them.
jobvr'V' to compute the right eigenvectors of A, 'N' to not compute them.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
wrPointer to the first element of the vector for the real part of the eigenvalues.
wiPointer to the first element of the vector for the imaginary part of the eigenvalues.
VLPointer to the first element of the column-major matrix for the left eigenvectors.
ldvlThe total number of elements between two columns of the matrix VL $[0..\infty)$.
VRPointer to the first element of the column-major matrix for the right eigenvectors.
ldvrThe total number of elements between two columns of the matrix VR $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
infoReturn code of the function call.
Returns
void

This function computes the eigenvalues of a non-symmetric n-by-n double precision column-major matrix based on the LAPACK dgeev() function. Optionally, it computes the left and/or right eigenvectors. The right eigenvector $v[j]$ of satisfies

\[ A * v[j] = lambda[j] * v[j], \]

where $lambda[j]$ is its eigenvalue. The left eigenvector $u[j]$ of A satisfies

\[ u[j]^{H} * A = lambda[j] * u[j]^{H}, \]

where $u[j]^{H}$ denotes the conjugate transpose of $u[j]$.

Complex conjugate pairs of eigenvalues appear consecutively with the eigenvalue having the positive imaginary part first. The computed eigenvectors are normalized to have Euclidean norm equal to 1 and largest component real.

The parameter jobvl specifies the computation of the left eigenvectors:

  • 'V': The left eigenvectors of A are computed and returned in VL.
  • 'N': The left eigenvectors of A are not computed.

The parameter jobvr specifies the computation of the right eigenvectors:

  • 'V': The right eigenvectors of A are computed and returned in VR.
  • 'N': The right eigenvectors of A are not computed.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the QR algorithm failed to compute all the eigenvalues and no eigenvectors have been computed; the elements with index larger than i have converged.

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

◆ geev() [7/8]

void blaze::geev ( char  jobvl,
char  jobvr,
int  n,
complex< float > *  A,
int  lda,
complex< float > *  w,
complex< float > *  VL,
int  ldvl,
complex< float > *  VR,
int  ldvr,
complex< float > *  work,
int  lwork,
float *  rwork,
int *  info 
)
inline

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

Parameters
jobvl'V' to compute the left eigenvectors of A, 'N' to not compute them.
jobvr'V' to compute the right eigenvectors of A, 'N' to not compute them.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
wPointer to the first element of the vector for the eigenvalues.
VLPointer to the first element of the column-major matrix for the left eigenvectors.
ldvlThe total number of elements between two columns of the matrix VL $[0..\infty)$.
VRPointer to the first element of the column-major matrix for the right eigenvectors.
ldvrThe total number of elements between two columns of the matrix VR $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
rworkAuxiliary array; size >= 2*n.
infoReturn code of the function call.
Returns
void

This function computes the eigenvalues of a non-symmetric n-by-n single precision complex column-major matrix based on the LAPACK cgeev() function. Optionally, it computes the left and/or right eigenvectors. The right eigenvector $v[j]$ of satisfies

\[ A * v[j] = lambda[j] * v[j], \]

where $lambda[j]$ is its eigenvalue. The left eigenvector $u[j]$ of A satisfies

\[ u[j]^{H} * A = lambda[j] * u[j]^{H}, \]

where $u[j]^{H}$ denotes the conjugate transpose of $u[j]$.

Complex conjugate pairs of eigenvalues appear consecutively with the eigenvalue having the positive imaginary part first. The computed eigenvectors are normalized to have Euclidean norm equal to 1 and largest component real.

The parameter jobvl specifies the computation of the left eigenvectors:

  • 'V': The left eigenvectors of A are computed and returned in VL.
  • 'N': The left eigenvectors of A are not computed.

The parameter jobvr specifies the computation of the right eigenvectors:

  • 'V': The right eigenvectors of A are computed and returned in VR.
  • 'N': The right eigenvectors of A are not computed.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the QR algorithm failed to compute all the eigenvalues and no eigenvectors have been computed; the elements with index larger than i have converged.

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

◆ geev() [8/8]

void blaze::geev ( char  jobvl,
char  jobvr,
int  n,
complex< double > *  A,
int  lda,
complex< double > *  w,
complex< double > *  VL,
int  ldvl,
complex< double > *  VR,
int  ldvr,
complex< double > *  work,
int  lwork,
double *  rwork,
int *  info 
)
inline

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

Parameters
jobvl'V' to compute the left eigenvectors of A, 'N' to not compute them.
jobvr'V' to compute the right eigenvectors of A, 'N' to not compute them.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
wPointer to the first element of the vector for the eigenvalues.
VLPointer to the first element of the column-major matrix for the left eigenvectors.
ldvlThe total number of elements between two columns of the matrix VL $[0..\infty)$.
VRPointer to the first element of the column-major matrix for the right eigenvectors.
ldvrThe total number of elements between two columns of the matrix VR $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
rworkAuxiliary array; size >= 2*n.
infoReturn code of the function call.
Returns
void

This function computes the eigenvalues of a non-symmetric n-by-n double precision complex column-major matrix based on the LAPACK zgeev() function. Optionally, it computes the left and/or right eigenvectors. The right eigenvector $v[j]$ of satisfies

\[ A * v[j] = lambda[j] * v[j], \]

where $lambda[j]$ is its eigenvalue. The left eigenvector $u[j]$ of A satisfies

\[ u[j]^{H} * A = lambda[j] * u[j]^{H}, \]

where $u[j]^{H}$ denotes the conjugate transpose of $u[j]$.

Complex conjugate pairs of eigenvalues appear consecutively with the eigenvalue having the positive imaginary part first. The computed eigenvectors are normalized to have Euclidean norm equal to 1 and largest component real.

The parameter jobvl specifies the computation of the left eigenvectors:

  • 'V': The left eigenvectors of A are computed and returned in VL.
  • 'N': The left eigenvectors of A are not computed.

The parameter jobvr specifies the computation of the right eigenvectors:

  • 'V': The right eigenvectors of A are computed and returned in VR.
  • 'N': The right eigenvectors of A are not computed.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the QR algorithm failed to compute all the eigenvalues and no eigenvectors have been computed; the elements with index larger than i have converged.

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

◆ heev() [1/3]

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

LAPACK kernel for computing the eigenvalues of the given dense Hermitian matrix.

Parameters
AThe given Hermitian matrix.
wThe resulting vector of eigenvalues.
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
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_argumentVector cannot be resized.
std::invalid_argumentInvalid jobz argument provided.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorEigenvalue computation failed.

This function computes the eigenvalues of a Hermitian n-by-n matrix based on the LAPACK heev() functions. Optionally, it computes the left and right eigenvectors.

The real eigenvalues are returned in ascending order in the given vector w. w is resized to the correct size (if possible and necessary). In case A is a row-major matrix, the left eigenvectors are returned in the rows of A, in case A is a column-major matrix, the right eigenvectors are returned in the columns of A.

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given jobz argument is neither 'V' nor 'N';
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the eigenvalue computation fails.

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

Examples:

DynamicMatrix<complex<double>,rowMajor> A( 5UL, 5UL ); // The Hermitian matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 5UL ); // The vector for the real eigenvalues
heev( A, w, 'V', 'L' );

For more information on the heev() functions (i.e. cheev() and zheev()) 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.

◆ heev() [2/3]

void blaze::heev ( char  jobz,
char  uplo,
int  n,
complex< float > *  A,
int  lda,
float *  w,
complex< float > *  work,
int  lwork,
float *  rwork,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense Hermitian single precision column-major matrix.

Parameters
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
wPointer to the first element of the vector for the eigenvalues.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
rworkAuxiliary array; size >= max(1,3*n-2).
infoReturn code of the function call.
Returns
void

This function computes the eigenvalues of an Hermitian n-by-n single precision complex column-major matrix based on the LAPACK cheev() function. Optionally, it computes the left and right eigenvectors. The real eigenvalues are returned in ascending order in the given n-dimensional vector w.

The parameter jobz specifies the computation of the orthonormal eigenvectors of A:

  • 'V': The eigenvectors of A are computed and returned within the matrix A.
  • 'N': The eigenvectors of A are not computed.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the algorithm failed to converge; i values did not converge to zero.

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

◆ heev() [3/3]

void blaze::heev ( char  jobz,
char  uplo,
int  n,
complex< double > *  A,
int  lda,
double *  w,
complex< double > *  work,
int  lwork,
double *  rwork,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense Hermitian double precision column-major matrix.

Parameters
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
wPointer to the first element of the vector for the eigenvalues.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
rworkAuxiliary array; size >= max(1,3*n-2).
infoReturn code of the function call.
Returns
void

This function computes the eigenvalues of an Hermitian n-by-n double precision complex column-major matrix based on the LAPACK zheev() function. Optionally, it computes the left and right eigenvectors. The real eigenvalues are returned in ascending order in the given n-dimensional vector w.

The parameter jobz specifies the computation of the orthonormal eigenvectors of A:

  • 'V': The eigenvectors of A are computed and returned within the matrix A.
  • 'N': The eigenvectors of A are not computed.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the algorithm failed to converge; i values did not converge to zero.

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

◆ heevd() [1/3]

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

LAPACK kernel for computing the eigenvalues of the given dense Hermitian matrix.

Parameters
AThe given Hermitian matrix.
wThe resulting vector of eigenvalues.
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
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_argumentVector cannot be resized.
std::invalid_argumentInvalid jobz argument provided.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorEigenvalue computation failed.

This function computes the eigenvalues of a Hermitian n-by-n matrix based on the LAPACK heevd() functions, which use a divide-and-conquer strategy. Optionally, it computes the left and right eigenvectors.

The real eigenvalues are returned in ascending order in the given vector w. w is resized to the correct size (if possible and necessary). In case A is a row-major matrix, the left eigenvectors are returned in the rows of A, in case A is a column-major matrix, the right eigenvectors are returned in the columns of A.

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given jobz argument is neither 'V' nor 'N';
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the eigenvalue computation fails.

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

Examples:

DynamicMatrix<complex<double>,rowMajor> A( 5UL, 5UL ); // The Hermitian matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 5UL ); // The vector for the real eigenvalues
heevd( A, w, 'V', 'L' );

For more information on the heevd() functions (i.e. cheevd() and zheevd()) 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.

◆ heevd() [2/3]

void blaze::heevd ( char  jobz,
char  uplo,
int  n,
complex< float > *  A,
int  lda,
float *  w,
complex< float > *  work,
int  lwork,
float *  rwork,
int  lrwork,
int *  iwork,
int  liwork,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense Hermitian single precision column-major matrix.

Parameters
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
wPointer to the first element of the vector for the eigenvalues.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
rworkAuxiliary array; size >= max( 1, lrwork ).
lrworkThe dimension of the array rwork; see online reference for details.
iworkAuxiliary array; size >= max( 1, liwork ).
liworkThe dimension of the array iwork; see online reference for details.
infoReturn code of the function call.
Returns
void

This function computes the eigenvalues of an Hermitian n-by-n single precision complex column-major matrix based on the LAPACK cheevd() function. Optionally, it computes the left and right eigenvectors using a divide-and-conquer strategy. The real eigenvalues are returned in ascending order in the given n-dimensional vector w.

The parameter jobz specifies the computation of the orthonormal eigenvectors of A:

  • 'V': The eigenvectors of A are computed and returned within the matrix A.
  • 'N': The eigenvectors of A are not computed.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the algorithm failed to converge; i values did not converge to zero.

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

◆ heevd() [3/3]

void blaze::heevd ( char  jobz,
char  uplo,
int  n,
complex< double > *  A,
int  lda,
double *  w,
complex< double > *  work,
int  lwork,
double *  rwork,
int  lrwork,
int *  iwork,
int  liwork,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense Hermitian double precision column-major matrix.

Parameters
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
wPointer to the first element of the vector for the eigenvalues.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
rworkAuxiliary array; size >= max( 1, lrwork ).
lrworkThe dimension of the array rwork; see online reference for details.
iworkAuxiliary array; size >= max( 1, liwork ).
liworkThe dimension of the array iwork; see online reference for details.
infoReturn code of the function call.
Returns
void

This function computes the eigenvalues of an Hermitian n-by-n double precision complex column-major matrix based on the LAPACK cheevd() function. Optionally, it computes the left and right eigenvectors using a divide-and-conquer strategy. The real eigenvalues are returned in ascending order in the given n-dimensional vector w.

The parameter jobz specifies the computation of the orthonormal eigenvectors of A:

  • 'V': The eigenvectors of A are computed and returned within the matrix A.
  • 'N': The eigenvectors of A are not computed.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the algorithm failed to converge; i values did not converge to zero.

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

◆ heevx() [1/6]

template<typename MT , bool SO, typename VT , bool TF>
size_t blaze::heevx ( DenseMatrix< MT, SO > &  A,
DenseVector< VT, TF > &  w,
char  uplo 
)
inline

LAPACK kernel for the eigenvalue computation of the given dense Hermitian matrix.

Parameters
AThe given Hermitian matrix.
wThe resulting vector of eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
Returns
The total number of eigenvalues found.
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorEigenvalue computation failed.

This function computes all eigenvalues of a Hermitian n-by-n matrix based on the LAPACK heevx() functions. The real eigenvalues are returned in ascending order in the given vector w. w is resized to the correct size (if possible and necessary).

Note that this function can only be used for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with any adapted matrix or matrices of any other element type results in a compile time error!

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the eigenvalue computation fails.

In all failure cases an exception is thrown.

Examples:

DynamicMatrix<complex<double>,rowMajor> A( 5UL, 5UL ); // The Hermitian matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 5UL ); // The vector for the real eigenvalues
heevx( A, w, 'L' );

For more information on the heevx() functions (i.e. cheevx() and zheevx()) 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.

◆ heevx() [2/6]

template<typename MT , bool SO, typename VT , bool TF, typename ST >
size_t blaze::heevx ( DenseMatrix< MT, SO > &  A,
DenseVector< VT, TF > &  w,
char  uplo,
ST  low,
ST  upp 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense Hermitian matrix.

Parameters
AThe given Hermitian matrix.
wThe resulting vector of eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
lowThe lower bound of the interval to be searched for eigenvalues.
uppThe upper bound of the interval to be searched for eigenvalues.
Returns
The total number of eigenvalues found.
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::invalid_argumentInvalid uplo argument provided.
std::invalid_argumentInvalid value range provided.
std::invalid_argumentInvalid index range provided.
std::runtime_errorEigenvalue computation failed.

This function computes a specific number of eigenvalues of a Hermitian n-by-n matrix based on the LAPACK heevx() functions. The number of eigenvalues to be computed is specified by the lower bound low and the upper bound upp, which either form an integral or a floating point range.

In case low and upp are of integral type, the function computes all eigenvalues in the index range $[low..upp]$. The num resulting real eigenvalues are stored in ascending order in the given vector w, which is either resized (if possible) or expected to be a num-dimensional vector.

In case low and upp are of floating point type, the function computes all eigenvalues in the half-open interval $(low..upp]$. The resulting real eigenvalues are stored in ascending order in the given vector w, which is either resized (if possible) or expected to be an n-dimensional vector.

Note that this function can only be used for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with any adapted matrix or matrices of any other element type results in a compile time error!

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given scalar values don't form a proper range;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the eigenvalue computation fails.

In all failure cases an exception is thrown.

Examples:

DynamicMatrix<complex<double>,rowMajor> A( 5UL, 5UL ); // The Hermitian matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 5UL ); // The vector for the real eigenvalues
heevx( A, w, 'L', 1.0, 2.0 ); // Computes all eigenvalues in the interval (1..2]
DynamicMatrix<complex<double>,rowMajor> A( 5UL, 5UL ); // The Hermitian matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 3UL ); // The vector for the real eigenvalues
heevx( A, w, 'L', 0, 2 ); // Computes the first three eigenvalues

For more information on the heevx() functions (i.e. cheevx() and zheevx()) 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.

◆ heevx() [3/6]

template<typename MT1 , bool SO1, typename VT , bool TF, typename MT2 , bool SO2>
size_t blaze::heevx ( DenseMatrix< MT1, SO1 > &  A,
DenseVector< VT, TF > &  w,
DenseMatrix< MT2, SO2 > &  Z,
char  uplo 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense Hermitian matrix.

Parameters
AThe given Hermitian matrix.
wThe resulting vector of eigenvalues.
ZThe resulting matrix of eigenvectors.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
Returns
The total number of eigenvalues found.
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::invalid_argumentMatrix cannot be resized.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorEigenvalue computation failed.

This function computes all eigenvalues of a Hermitian n-by-n matrix based on the LAPACK heevx() functions. Additionally, it computes all eigenvectors of A. The real eigenvalues are returned in ascending order in the given vector w. The eigenvectors are returned in the rows of Z in case Z is a row-major matrix and in the columns of Z in case Z is a column-major matrix. Both w and Z are resized to the correct dimensions (if possible and necessary).

Note that this function can only be used for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with any adapted matrix or matrices of any other element type results in a compile time error!

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given matrix Z is a fixed size matrix and the dimensions don't match;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the eigenvalue computation fails.

In all failure cases an exception is thrown.

Examples:

DynamicMatrix<complex<double>,rowMajor> A( 5UL, 5UL ); // The Hermitian matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 5UL ); // The vector for the real eigenvalues
DynamicMatrix<complex<double>,rowMajor> Z( 5UL, 5UL ); // The matrix for the complex eigenvectors
heevx( A, w, 'L' );

For more information on the heevx() functions (i.e. cheevx() and zheevx()) 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.

◆ heevx() [4/6]

template<typename MT1 , bool SO1, typename VT , bool TF, typename MT2 , bool SO2, typename ST >
size_t blaze::heevx ( DenseMatrix< MT1, SO1 > &  A,
DenseVector< VT, TF > &  w,
DenseMatrix< MT2, SO2 > &  Z,
char  uplo,
ST  low,
ST  upp 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense Hermitian matrix.

Parameters
AThe given Hermitian matrix.
wThe resulting vector of eigenvalues.
ZThe resulting matrix of eigenvectors.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
lowThe lower bound of the interval to be searched for eigenvalues.
uppThe upper bound of the interval to be searched for eigenvalues.
Returns
The total number of eigenvalues found.
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::invalid_argumentMatrix cannot be resized.
std::invalid_argumentInvalid uplo argument provided.
std::invalid_argumentInvalid value range provided.
std::invalid_argumentInvalid index range provided.
std::runtime_errorEigenvalue computation failed.

This function computes a specific number of eigenvalues of a Hermitian n-by-n matrix based on the LAPACK heevx() functions. Additionally, it computes a specific number of eigenvectors of A. The number of eigenvalues to be computed is specified by the lower bound low and the upper bound upp, which either form an integral or a floating point range.

In case low and upp are of integral type, the function computes all eigenvalues in the index range $[low..upp]$. The num resulting real eigenvalues are stored in ascending order in the given vector w, which is either resized (if possible) or expected to be a num-dimensional vector. The eigenvectors are returned in the rows of Z in case Z is row-major matrix and in the columns of Z in case Z is a column-major matrix. Z is resized (if possible) or expected to be a num-by-n row-major matrix or a n-by-num column-major matrix.

In case low and upp are of floating point type, the function computes all eigenvalues in the half-open interval $(low..upp]$. The resulting real eigenvalues are stored in ascending order in the given vector w, which is either resized (if possible) or expected to be an n-dimensional vector. The eigenvectors are returned in the rows of Z in case Z is a row-major matrix and in the columns of Z in case Z is a column-major matrix. Z is resized (if possible) or expected to be a n-by-n matrix.

Note that this function can only be used for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with any adapted matrix or matrices of any other element type results in a compile time error!

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given matrix Z is a fixed size matrix and the dimensions don't match;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the eigenvalue computation fails.

In all failure cases an exception is thrown.

Examples:

DynamicMatrix<complex<double>,rowMajor> A( 5UL, 5UL ); // The Hermitian matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 5UL ); // The vector for the real eigenvalues
DynamicMatrix<complex<double>,rowMajor> Z( 5UL, 5UL ); // The matrix for the complex eigenvectors
heevx( A, w, Z, 'L', 1.0, 2.0 ); // Computes all eigenvalues in the interval (1..2]
DynamicMatrix<complex<double>,rowMajor> A( 5UL, 5UL ); // The Hermitian matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 3UL ); // The vector for the real eigenvalues
DynamicMatrix<complex<double>,rowMajor> Z( 3UL, 5UL ); // The matrix for the complex eigenvectors
heevx( A, w, Z, 'L', 0, 2 ); // Computes the first three eigenvalues

For more information on the heevx() functions (i.e. cheevx() and zheevx()) 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.

◆ heevx() [5/6]

void blaze::heevx ( char  jobz,
char  range,
char  uplo,
int  n,
complex< float > *  A,
int  lda,
float  vl,
float  vu,
int  il,
int  iu,
float  abstol,
int *  m,
float *  w,
complex< float > *  Z,
int  ldz,
complex< float > *  work,
int  lwork,
float *  rwork,
int *  iwork,
int *  ifail,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense Hermitian single precision complex column-major matrix.

Parameters
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
rangeSpecifies the range of eigenvalues to find ('A', 'V', or 'I').
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
vlThe lower bound of the interval to be searched for eigenvalues (vl < vu).
vuThe upper bound of the interval to be searched for eigenvalues (vl < vu).
ilThe index of the smallest eigenvalue to be returned (0 <= il <= iu).
iuThe index of the largest eigenvalue to be returned (0 <= il <= iu).
abstolThe absolute error tolerance for the computation of the eigenvalues.
mThe total number of eigenvalues found (0 <= m <= n).
wPointer to the first element of the vector for the eigenvalues.
ZPointer to the first element of the column-major matrix for the eigenvectors.
ldzThe total number of elements between two columns of the matrix Z $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
rworkAuxiliary array; size >= 7*n.
iworkAuxiliary array; size >= 5*n.
ifailIndex array of eigenvectors that failed to converge.
infoReturn code of the function call.
Returns
void

This function computes a specified number of eigenvalues of a Hermitian n-by-n single precision complex column-major matrix based on the LAPACK cheevx() function. Optionally, it computes a specified number of eigenvectors. The selected real eigenvalues are returned in ascending order in the given n-dimensional vector w.

The parameter jobz specifies the computation of the orthonormal eigenvectors of A:

  • 'V': The eigenvectors of A are computed and returned within the matrix Z.
  • 'N': The eigenvectors of A are not computed.

The parameter range specifies the amount of eigenvalues/eigenvectors to be found:

  • 'A': All eigenvalues will be found.
  • 'V': All eigenvalues in the half-open interval $(vl..vu]$ will be found.
  • 'I': The il-th through iu-th eigenvalues will be found.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the algorithm failed to converge; i values did not converge to zero.

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

◆ heevx() [6/6]

void blaze::heevx ( char  jobz,
char  range,
char  uplo,
int  n,
complex< double > *  A,
int  lda,
double  vl,
double  vu,
int  il,
int  iu,
double  abstol,
int *  m,
double *  w,
complex< double > *  Z,
int  ldz,
complex< double > *  work,
int  lwork,
double *  rwork,
int *  iwork,
int *  ifail,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense Hermitian double precision complex column-major matrix.

Parameters
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
rangeSpecifies the range of eigenvalues to find ('A', 'V', or 'I').
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
vlThe lower bound of the interval to be searched for eigenvalues (vl < vu).
vuThe upper bound of the interval to be searched for eigenvalues (vl < vu).
ilThe index of the smallest eigenvalue to be returned (0 <= il <= iu).
iuThe index of the largest eigenvalue to be returned (0 <= il <= iu).
abstolThe absolute error tolerance for the computation of the eigenvalues.
mThe total number of eigenvalues found (0 <= m <= n).
wPointer to the first element of the vector for the eigenvalues.
ZPointer to the first element of the column-major matrix for the eigenvectors.
ldzThe total number of elements between two columns of the matrix Z $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
rworkAuxiliary array; size >= 7*n.
iworkAuxiliary array; size >= 5*n.
ifailIndex array of eigenvectors that failed to converge.
infoReturn code of the function call.
Returns
void

This function computes a specified number of eigenvalues of a Hermitian n-by-n double precision complex column-major matrix based on the LAPACK cheevx() function. Optionally, it computes a specified number of eigenvectors. The selected real eigenvalues are returned in ascending order in the given n-dimensional vector w.

The parameter jobz specifies the computation of the orthonormal eigenvectors of A:

  • 'V': The eigenvectors of A are computed and returned within the matrix Z.
  • 'N': The eigenvectors of A are not computed.

The parameter range specifies the amount of eigenvalues/eigenvectors to be found:

  • 'A': All eigenvalues will be found.
  • 'V': All eigenvalues in the half-open interval $(vl..vu]$ will be found.
  • 'I': The il-th through iu-th eigenvalues will be found.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the algorithm failed to converge; i values did not converge to zero.

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

◆ syev() [1/3]

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

LAPACK kernel for computing the eigenvalues of the given dense symmetric matrix.

Parameters
AThe given symmetric matrix.
wThe resulting vector of eigenvalues.
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
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_argumentVector cannot be resized.
std::invalid_argumentInvalid jobz argument provided.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorEigenvalue computation failed.

This function computes the eigenvalues of a symmetric n-by-n matrix based on the LAPACK syev() functions. Optionally, it computes the left or right eigenvectors.

The real eigenvalues are returned in ascending order in the given vector w. w is resized to the correct size (if possible and necessary). In case A is a row-major matrix, the left eigenvectors are returned in the rows of A, in case A is a column-major matrix, the right eigenvectors are returned in the columns of A.

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given jobz argument is neither 'V' nor 'N';
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the eigenvalue computation fails.

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

Examples:

DynamicMatrix<double,rowMajor> A( 5UL, 5UL ); // The symmetric matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 5UL ); // The vector for the real eigenvalues
syev( A, w, 'V', 'L' );

For more information on the syev() functions (i.e. ssyev() and dsyev()) 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.

◆ syev() [2/3]

void blaze::syev ( char  jobz,
char  uplo,
int  n,
float *  A,
int  lda,
float *  w,
float *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense symmetric single precision column-major matrix.

Parameters
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
wPointer to the first element of the vector for the eigenvalues.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
infoReturn code of the function call.
Returns
void

This function computes the eigenvalues of a symmetric n-by-n single precision column-major matrix based on the LAPACK ssyev() function. Optionally, it computes the left and right eigenvectors. The real eigenvalues are returned in ascending order in the given n-dimensional vector w.

The parameter jobz specifies the computation of the orthonormal eigenvectors of A:

  • 'V': The eigenvectors of A are computed and returned within the matrix A.
  • 'N': The eigenvectors of A are not computed.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the algorithm failed to converge; i values did not converge to zero.

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

◆ syev() [3/3]

void blaze::syev ( char  jobz,
char  uplo,
int  n,
double *  A,
int  lda,
double *  w,
double *  work,
int  lwork,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense symmetric double precision column-major matrix.

Parameters
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
wPointer to the first element of the vector for the eigenvalues.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
infoReturn code of the function call.
Returns
void

This function computes the eigenvalues of a symmetric n-by-n double precision column-major matrix based on the LAPACK dsyev() function. Optionally, it computes the left and right eigenvectors. The real eigenvalues are returned in ascending order in the given n-dimensional vector w.

The parameter jobz specifies the computation of the orthonormal eigenvectors of A:

  • 'V': The eigenvectors of A are computed and returned within the matrix A.
  • 'N': The eigenvectors of A are not computed.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the algorithm failed to converge; i values did not converge to zero.

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

◆ syevd() [1/3]

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

LAPACK kernel for computing the eigenvalues of the given dense symmetric matrix.

Parameters
AThe given symmetric matrix.
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
wThe resulting vector of eigenvalues.
Returns
void
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::invalid_argumentInvalid jobz argument provided.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorEigenvalue computation failed.

This function computes the eigenvalues of a symmetric n-by-n matrix based on the LAPACK syevd() functions, which use a divide-and-conquer strategy. Optionally, it computes the left or right eigenvectors.

The real eigenvalues are returned in ascending order in the given vector w. w is resized to the correct size (if possible and necessary). In case A is a row-major matrix, the left eigenvectors are returned in the rows of A, in case A is a column-major matrix, the right eigenvectors are returned in the columns of A.

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given jobz argument is neither 'V' nor 'N';
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the eigenvalue computation fails.

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

Examples:

DynamicMatrix<double,rowMajor> A( 5UL, 5UL ); // The symmetric matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 5UL ); // The vector for the real eigenvalues
syevd( A, w, 'V', 'L' );

For more information on the syevd() functions (i.e. ssyevd() and dsyevd()) 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.

◆ syevd() [2/3]

void blaze::syevd ( char  jobz,
char  uplo,
int  n,
float *  A,
int  lda,
float *  w,
float *  work,
int  lwork,
int *  iwork,
int  liwork,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense symmetric single precision column-major matrix.

Parameters
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
wPointer to the first element of the vector for the eigenvalues.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
iworkAuxiliary array; size >= max( 1, liwork ).
liworkThe dimension of the array iwork; see online reference for details.
infoReturn code of the function call.
Returns
void

This function computes the eigenvalues of a symmetric n-by-n single precision column-major matrix based on the LAPACK ssyevd() function. Optionally, it computes the left and right eigenvectors using a divide-and-conquer strategy. The real eigenvalues are returned in ascending order in the given n-dimensional vector w.

The parameter jobz specifies the computation of the orthonormal eigenvectors of A:

  • 'V': The eigenvectors of A are computed and returned within the matrix A.
  • 'N': The eigenvectors of A are not computed.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the algorithm failed to converge; i values did not converge to zero.

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

◆ syevd() [3/3]

void blaze::syevd ( char  jobz,
char  uplo,
int  n,
double *  A,
int  lda,
double *  w,
double *  work,
int  lwork,
int *  iwork,
int  liwork,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense symmetric double precision column-major matrix.

Parameters
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
wPointer to the first element of the vector for the eigenvalues.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
iworkAuxiliary array; size >= max( 1, liwork ).
liworkThe dimension of the array iwork; see online reference for details.
infoReturn code of the function call.
Returns
void

This function computes the eigenvalues of a symmetric n-by-n double precision column-major matrix based on the LAPACK ssyevd() function. Optionally, it computes the left and right eigenvectors using a divide-and-conquer strategy. The real eigenvalues are returned in ascending order in the given n-dimensional vector w.

The parameter jobz specifies the computation of the orthonormal eigenvectors of A:

  • 'V': The eigenvectors of A are computed and returned within the matrix A.
  • 'N': The eigenvectors of A are not computed.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the algorithm failed to converge; i values did not converge to zero.

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

◆ syevx() [1/6]

template<typename MT , bool SO, typename VT , bool TF>
size_t blaze::syevx ( DenseMatrix< MT, SO > &  A,
DenseVector< VT, TF > &  w,
char  uplo 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense symmetric matrix.

Parameters
AThe given symmetric matrix.
wThe resulting vector of eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
Returns
The total number of eigenvalues found.
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorEigenvalue computation failed.

This function computes all eigenvalues of a symmetric n-by-n matrix based on the LAPACK syevx() functions. The real eigenvalues are returned in ascending order in the given vector w. w is resized to the correct size (if possible and necessary).

Note that this function can only be used for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with any adapted matrix or matrices of any other element type results in a compile time error!

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the eigenvalue computation fails.

In all failure cases an exception is thrown.

Examples:

DynamicMatrix<double,rowMajor> A( 5UL, 5UL ); // The symmetric matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 5UL ); // The vector for the real eigenvalues
syevx( A, w, 'L' );

For more information on the syevx() functions (i.e. ssyevx() and dsyevx()) 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.

◆ syevx() [2/6]

template<typename MT , bool SO, typename VT , bool TF, typename ST >
size_t blaze::syevx ( DenseMatrix< MT, SO > &  A,
DenseVector< VT, TF > &  w,
char  uplo,
ST  low,
ST  upp 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense symmetric matrix.

Parameters
AThe given symmetric matrix.
wThe resulting vector of eigenvalues.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
lowThe lower bound of the interval to be searched for eigenvalues.
uppThe upper bound of the interval to be searched for eigenvalues.
Returns
The total number of eigenvalues found.
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::invalid_argumentInvalid uplo argument provided.
std::invalid_argumentInvalid value range provided.
std::invalid_argumentInvalid index range provided.
std::runtime_errorEigenvalue computation failed.

This function computes a specific number of eigenvalues of a symmetric n-by-n matrix based on the LAPACK syevx() functions. The number of eigenvalues to be computed is specified by the lower bound low and the upper bound upp, which either form an integral or a floating point range.

In case low and upp are of integral type, the function computes all eigenvalues in the index range $[low..upp]$. The num resulting real eigenvalues are stored in ascending order in the given vector w, which is either resized (if possible) or expected to be a num-dimensional vector.

In case low and upp are of floating point type, the function computes all eigenvalues in the half-open interval $(low..upp]$. The resulting real eigenvalues are stored in ascending order in the given vector w, which is either resized (if possible) or expected to be an n-dimensional vector.

Note that this function can only be used for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with any adapted matrix or matrices of any other element type results in a compile time error!

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given scalar values don't form a proper range;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the eigenvalue computation fails.

In all failure cases an exception is thrown.

Examples:

DynamicMatrix<double,rowMajor> A( 5UL, 5UL ); // The symmetric matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 5UL ); // The vector for the real eigenvalues
syevx( A, w, 'L', 1.0, 2.0 ); // Computes all eigenvalues in the interval (1..2]
DynamicMatrix<double,rowMajor> A( 5UL, 5UL ); // The symmetric matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 3UL ); // The vector for the real eigenvalues
syevx( A, w, 'L', 0, 2 ); // Computes the first three eigenvalues

For more information on the syevx() functions (i.e. ssyevx() and dsyevx()) 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.

◆ syevx() [3/6]

template<typename MT1 , bool SO1, typename VT , bool TF, typename MT2 , bool SO2>
size_t blaze::syevx ( DenseMatrix< MT1, SO1 > &  A,
DenseVector< VT, TF > &  w,
DenseMatrix< MT2, SO2 > &  Z,
char  uplo 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense symmetric matrix.

Parameters
AThe given symmetric matrix.
wThe resulting vector of eigenvalues.
ZThe resulting matrix of eigenvectors.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
Returns
The total number of eigenvalues found.
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::invalid_argumentMatrix cannot be resized.
std::invalid_argumentInvalid uplo argument provided.
std::runtime_errorEigenvalue computation failed.

This function computes all eigenvalues of a symmetric n-by-n matrix based on the LAPACK syevx() functions. Additionally, it computes all eigenvectors of A. The real eigenvalues are returned in ascending order in the given vector w. The eigenvectors are returned in the rows of Z in case Z is a row-major matrix and in the columns of Z in case Z is a column-major matrix. Both w and Z are resized to the correct dimensions (if possible and necessary).

Note that this function can only be used for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with any adapted matrix or matrices of any other element type results in a compile time error!

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given matrix Z is a fixed size matrix and the dimensions don't match;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the eigenvalue computation fails.

In all failure cases an exception is thrown.

Examples:

DynamicMatrix<double,rowMajor> A( 5UL, 5UL ); // The symmetric matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 5UL ); // The vector for the real eigenvalues
DynamicMatrix<double,rowMajor> Z( 5UL, 5UL ); // The matrix for the real eigenvectors
syevx( A, w, 'L' );

For more information on the syevx() functions (i.e. ssyevx() and dsyevx()) 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.

◆ syevx() [4/6]

template<typename MT1 , bool SO1, typename VT , bool TF, typename MT2 , bool SO2, typename ST >
size_t blaze::syevx ( DenseMatrix< MT1, SO1 > &  A,
DenseVector< VT, TF > &  w,
DenseMatrix< MT2, SO2 > &  Z,
char  uplo,
ST  low,
ST  upp 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense symmetric matrix.

Parameters
AThe given symmetric matrix.
wThe resulting vector of eigenvalues.
ZThe resulting matrix of eigenvectors.
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
lowThe lower bound of the interval to be searched for eigenvalues.
uppThe upper bound of the interval to be searched for eigenvalues.
Returns
The total number of eigenvalues found.
Exceptions
std::invalid_argumentInvalid non-square matrix provided.
std::invalid_argumentVector cannot be resized.
std::invalid_argumentMatrix cannot be resized.
std::invalid_argumentInvalid uplo argument provided.
std::invalid_argumentInvalid value range provided.
std::invalid_argumentInvalid index range provided.
std::runtime_errorEigenvalue computation failed.

This function computes a specific number of eigenvalues of a symmetric n-by-n matrix based on the LAPACK syevx() functions. Additionally, it computes a specific number of eigenvectors of A. The number of eigenvalues to be computed is specified by the lower bound low and the upper bound upp, which either form an integral or a floating point range.

In case low and upp are of integral type, the function computes all eigenvalues in the index range $[low..upp]$. The num resulting real eigenvalues are stored in ascending order in the given vector w, which is either resized (if possible) or expected to be a num-dimensional vector. The eigenvectors are returned in the rows of Z in case Z is row-major matrix and in the columns of Z in case Z is a column-major matrix. Z is resized (if possible) or expected to be a num-by-n row-major matrix or a n-by-num column-major matrix.

In case low and upp are of floating point type, the function computes all eigenvalues in the half-open interval $(low..upp]$. The resulting real eigenvalues are stored in ascending order in the given vector w, which is either resized (if possible) or expected to be an n-dimensional vector. The eigenvectors are returned in the rows of Z in case Z is a row-major matrix and in the columns of Z in case Z is a column-major matrix. Z is resized (if possible) or expected to be a n-by-n matrix.

Note that this function can only be used for general, non-adapted matrices with float, double, complex<float>, or complex<double> element type. The attempt to call the function with any adapted matrix or matrices of any other element type results in a compile time error!

The function fails if ...

  • ... the given matrix A is not a square matrix;
  • ... the given vector w is a fixed size vector and the size doesn't match;
  • ... the given matrix Z is a fixed size matrix and the dimensions don't match;
  • ... the given uplo argument is neither 'L' nor 'U';
  • ... the eigenvalue computation fails.

In all failure cases an exception is thrown.

Examples:

DynamicMatrix<double,rowMajor> A( 5UL, 5UL ); // The Hermitian matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 5UL ); // The vector for the real eigenvalues
DynamicMatrix<double,rowMajor> Z( 5UL, 5UL ); // The matrix for the real eigenvectors
syevx( A, w, Z, 'L', 1.0, 2.0 ); // Computes all eigenvalues in the interval (1..2]
DynamicMatrix<double,rowMajor> A( 5UL, 5UL ); // The Hermitian matrix A
// ... Initialization
DynamicVector<double,columnVector> w( 3UL ); // The vector for the real eigenvalues
DynamicMatrix<double,rowMajor> Z( 3UL, 5UL ); // The matrix for the real eigenvectors
syevx( A, w, Z, 'L', 0, 2 ); // Computes the first three eigenvalues

For more information on the syevx() functions (i.e. ssyevx() and dsyevx()) 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.

◆ syevx() [5/6]

void blaze::syevx ( char  jobz,
char  range,
char  uplo,
int  n,
float *  A,
int  lda,
float  vl,
float  vu,
int  il,
int  iu,
float  abstol,
int *  m,
float *  w,
float *  Z,
int  ldz,
float *  work,
int  lwork,
int *  iwork,
int *  ifail,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense symmetric single precision column-major matrix.

Parameters
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
rangeSpecifies the range of eigenvalues to find ('A', 'V', or 'I').
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
vlThe lower bound of the interval to be searched for eigenvalues (vl < vu).
vuThe upper bound of the interval to be searched for eigenvalues (vl < vu).
ilThe index of the smallest eigenvalue to be returned (0 <= il <= iu).
iuThe index of the largest eigenvalue to be returned (0 <= il <= iu).
abstolThe absolute error tolerance for the computation of the eigenvalues.
mThe total number of eigenvalues found (0 <= m <= n).
wPointer to the first element of the vector for the eigenvalues.
ZPointer to the first element of the column-major matrix for the eigenvectors.
ldzThe total number of elements between two columns of the matrix Z $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
iworkAuxiliary array; size >= 5*n.
ifailIndex array of eigenvectors that failed to converge.
infoReturn code of the function call.
Returns
void

This function computes a specified number of eigenvalues of a symmetric n-by-n single precision column-major matrix based on the LAPACK ssyevx() function. Optionally, it computes a specified number of eigenvectors. The selected real eigenvalues are returned in ascending order in the given n-dimensional vector w.

The parameter jobz specifies the computation of the orthonormal eigenvectors of A:

  • 'V': The eigenvectors of A are computed and returned within the matrix Z.
  • 'N': The eigenvectors of A are not computed.

The parameter range specifies the amount of eigenvalues/eigenvectors to be found:

  • 'A': All eigenvalues will be found.
  • 'V': All eigenvalues in the half-open interval $(vl..vu]$ will be found.
  • 'I': The il-th through iu-th eigenvalues will be found.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the algorithm failed to converge; i values did not converge to zero.

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

◆ syevx() [6/6]

void blaze::syevx ( char  jobz,
char  range,
char  uplo,
int  n,
double *  A,
int  lda,
double  vl,
double  vu,
int  il,
int  iu,
double  abstol,
int *  m,
double *  w,
double *  Z,
int  ldz,
double *  work,
int  lwork,
int *  iwork,
int *  ifail,
int *  info 
)
inline

LAPACK kernel for computing the eigenvalues of the given dense symmetric double precision column-major matrix.

Parameters
jobz'V' to compute the eigenvectors of A, 'N' to only compute the eigenvalues.
rangeSpecifies the range of eigenvalues to find ('A', 'V', or 'I').
uplo'L' to use the lower part of the matrix, 'U' to use the upper part.
nThe number of rows and columns of the given 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 A $[0..\infty)$.
vlThe lower bound of the interval to be searched for eigenvalues (vl < vu).
vuThe upper bound of the interval to be searched for eigenvalues (vl < vu).
ilThe index of the smallest eigenvalue to be returned (0 <= il <= iu).
iuThe index of the largest eigenvalue to be returned (0 <= il <= iu).
abstolThe absolute error tolerance for the computation of the eigenvalues.
mThe total number of eigenvalues found (0 <= m <= n).
wPointer to the first element of the vector for the eigenvalues.
ZPointer to the first element of the column-major matrix for the eigenvectors.
ldzThe total number of elements between two columns of the matrix Z $[0..\infty)$.
workAuxiliary array; size >= max( 1, lwork ).
lworkThe dimension of the array work; see online reference for details.
iworkAuxiliary array; size >= 5*n.
ifailIndex array of eigenvectors that failed to converge.
infoReturn code of the function call.
Returns
void

This function computes a specified number of eigenvalues of a symmetric n-by-n double precision column-major matrix based on the LAPACK dsyevx() function. Optionally, it computes a specified number of eigenvectors. The selected real eigenvalues are returned in ascending order in the given n-dimensional vector w.

The parameter jobz specifies the computation of the orthonormal eigenvectors of A:

  • 'V': The eigenvectors of A are computed and returned within the matrix Z.
  • 'N': The eigenvectors of A are not computed.

The parameter range specifies the amount of eigenvalues/eigenvectors to be found:

  • 'A': All eigenvalues will be found.
  • 'V': All eigenvalues in the half-open interval $(vl..vu]$ will be found.
  • 'I': The il-th through iu-th eigenvalues will be found.

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

  • = 0: The computation finished successfully.
  • < 0: If info = -i, the i-th argument had an illegal value.
  • > 0: If info = i, the algorithm failed to converge; i values did not converge to zero.

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