Classes | Functions
Eigen::internal Namespace Reference

Classes

class  BandMatrix
 Represents a rectangular matrix with a banded storage. More...
struct  HessenbergDecompositionMatrixHReturnType
 Expression type for return value of HessenbergDecomposition::matrixH() More...
class  TridiagonalMatrix
 Represents a tridiagonal matrix with a compact banded storage. More...

Functions

template<typename LhsScalar , typename RhsScalar , int KcFactor>
void computeProductBlockingSizes (std::ptrdiff_t &k, std::ptrdiff_t &m, std::ptrdiff_t &n)
 Computes the blocking parameters for a m x k times k x n matrix product.
int nbThreads ()
void setNbThreads (int v)
template<typename RealScalar >
std::complex< RealScalar > sqrt (const std::complex< RealScalar > &z)
template<typename MatrixType , typename DiagonalType , typename SubDiagonalType >
void tridiagonalization_inplace (MatrixType &mat, DiagonalType &diag, SubDiagonalType &subdiag, bool extractQ)
 Performs a full tridiagonalization in place.

Detailed Description

Defines the maximal loop size to enable meta unrolling of loops. Note that the value here is expressed in Eigen's own notion of "number of FLOPS", it does not correspond to the number of iterations or the number of instructions Defines the threshold between a "small" and a "large" matrix. This threshold is mainly used to select the proper product implementation. Defines the maximal width of the blocks used in the triangular product and solver for vectors (level 2 blas xTRMV and xTRSV). The default is 8. Defines the default number of registers available for that architecture. Currently it must be 8 or 16. Other values will fail.

This is defined in the Jacobi module.

 #include <Eigen/Jacobi> 

Applies the clock wise 2D rotation j to the set of 2D vectors of cordinates x and y: $ \left ( \begin{array}{cc} x \\ y \end{array} \right ) = J \left ( \begin{array}{cc} x \\ y \end{array} \right ) $

See also:
MatrixBase::applyOnTheLeft(), MatrixBase::applyOnTheRight()

Function Documentation

void Eigen::internal::computeProductBlockingSizes ( std::ptrdiff_t &  k,
std::ptrdiff_t &  m,
std::ptrdiff_t &  n 
)

Computes the blocking parameters for a m x k times k x n matrix product.

Parameters:
[in,out]kInput: the third dimension of the product. Output: the blocking size along the same dimension.
[in,out]mInput: the number of rows of the left hand side. Output: the blocking size along the same dimension.
[in,out]nInput: the number of columns of the right hand side. Output: the blocking size along the same dimension.

Given a m x k times k x n matrix product of scalar types LhsScalar and RhsScalar, this function computes the blocking size parameters along the respective dimensions for matrix products and related algorithms. The blocking sizes depends on various parameters:

  • the L1 and L2 cache sizes,
  • the register level blocking sizes defined by gebp_traits,
  • the number of scalars that fit into a packet (when vectorization is enabled).
See also:
setCpuCacheSizes
int Eigen::internal::nbThreads ( ) [inline]
Returns:
the max number of threads reserved for Eigen
See also:
setNbThreads
void Eigen::internal::setNbThreads ( int  v) [inline]

Sets the max number of threads reserved for Eigen

See also:
nbThreads
std::complex< RealScalar > sqrt ( const std::complex< RealScalar > &  z)

Computes the principal value of the square root of the complex z.

void tridiagonalization_inplace ( MatrixType &  mat,
DiagonalType &  diag,
SubDiagonalType &  subdiag,
bool  extractQ 
)

Performs a full tridiagonalization in place.

Parameters:
[in,out]matOn input, the selfadjoint matrix whose tridiagonal decomposition is to be computed. Only the lower triangular part referenced. The rest is left unchanged. On output, the orthogonal matrix Q in the decomposition if extractQ is true.
[out]diagThe diagonal of the tridiagonal matrix T in the decomposition.
[out]subdiagThe subdiagonal of the tridiagonal matrix T in the decomposition.
[in]extractQIf true, the orthogonal matrix Q in the decomposition is computed and stored in mat.

Computes the tridiagonal decomposition of the selfadjoint matrix mat in place such that $ mat = Q T Q^* $ where $ Q $ is unitary and $ T $ a real symmetric tridiagonal matrix.

The tridiagonal matrix T is passed to the output parameters diag and subdiag. If extractQ is true, then the orthogonal matrix Q is passed to mat. Otherwise the lower part of the matrix mat is destroyed.

The vectors diag and subdiag are not resized. The function assumes that they are already of the correct size. The length of the vector diag should equal the number of rows in mat, and the length of the vector subdiag should be one left.

This implementation contains an optimized path for 3-by-3 matrices which is especially useful for plane fitting.

Note:
Currently, it requires two temporary vectors to hold the intermediate Householder coefficients, and to reconstruct the matrix Q from the Householder reflectors.

Example (this uses the same matrix as the example in Tridiagonalization::Tridiagonalization(const MatrixType&)):

MatrixXd X = MatrixXd::Random(5,5);
MatrixXd A = X + X.transpose();
cout << "Here is a random symmetric 5x5 matrix:" << endl << A << endl << endl;

VectorXd diag(5);
VectorXd subdiag(4);
internal::tridiagonalization_inplace(A, diag, subdiag, true);
cout << "The orthogonal matrix Q is:" << endl << A << endl;
cout << "The diagonal of the tridiagonal matrix T is:" << endl << diag << endl;
cout << "The subdiagonal of the tridiagonal matrix T is:" << endl << subdiag << endl;

Output:

Here is a random symmetric 5x5 matrix:
  1.36 -0.816  0.521   1.43 -0.144
-0.816 -0.659  0.794 -0.173 -0.406
 0.521  0.794 -0.541  0.461  0.179
  1.43 -0.173  0.461  -1.43  0.822
-0.144 -0.406  0.179  0.822  -1.37

The orthogonal matrix Q is:
       1        0        0        0        0
       0   -0.471    0.127   -0.671   -0.558
       0    0.301   -0.195    0.437   -0.825
       0    0.825   0.0459   -0.563 -0.00872
       0  -0.0832   -0.971   -0.202   0.0922
The diagonal of the tridiagonal matrix T is:
1.36
-1.2
-1.28
-1.69
0.164
The subdiagonal of the tridiagonal matrix T is:
1.73
-0.966
0.214
0.345
See also:
class Tridiagonalization