All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Groups Pages
Eigen::internal Namespace Reference

Functions

template<typename TMatrix , typename CMatrix , typename VectorX , typename VectorB , typename VectorF >
void constrained_cg (const TMatrix &A, const CMatrix &C, VectorX &x, const VectorB &b, const VectorF &f, IterationController &iter)
 
template<typename MatrixType , typename Rhs , typename Dest , typename Preconditioner >
bool gmres (const MatrixType &mat, const Rhs &rhs, Dest &x, const Preconditioner &precond, Index &iters, const Index &restart, typename Dest::RealScalar &tol_error)
 
template<typename MatA , typename MatU , typename MatV >
void matrix_exp_pade13 (const MatA &A, MatU &U, MatV &V)
 Compute the (13,13)-Padé approximant to the exponential. More...
 
template<typename MatA , typename MatU , typename MatV >
void matrix_exp_pade3 (const MatA &A, MatU &U, MatV &V)
 Compute the (3,3)-Padé approximant to the exponential. More...
 
template<typename MatA , typename MatU , typename MatV >
void matrix_exp_pade5 (const MatA &A, MatU &U, MatV &V)
 Compute the (5,5)-Padé approximant to the exponential. More...
 
template<typename MatA , typename MatU , typename MatV >
void matrix_exp_pade7 (const MatA &A, MatU &U, MatV &V)
 Compute the (7,7)-Padé approximant to the exponential. More...
 
template<typename MatA , typename MatU , typename MatV >
void matrix_exp_pade9 (const MatA &A, MatU &U, MatV &V)
 Compute the (9,9)-Padé approximant to the exponential. More...
 
template<typename MatrixType , typename VectorType >
void matrix_function_compute_above_diagonal (const MatrixType &T, const VectorType &blockStart, const VectorType &clusterSize, MatrixType &fT)
 Compute part of matrix function above block diagonal. More...
 
template<typename MatrixType , typename AtomicType , typename VectorType >
void matrix_function_compute_block_atomic (const MatrixType &T, AtomicType &atomic, const VectorType &blockStart, const VectorType &clusterSize, MatrixType &fT)
 Compute block diagonal part of matrix function. More...
 
template<typename VectorType >
void matrix_function_compute_block_start (const VectorType &clusterSize, VectorType &blockStart)
 Compute start of each block using clusterSize.
 
template<typename ListOfClusters , typename Index >
void matrix_function_compute_cluster_size (const ListOfClusters &clusters, Matrix< Index, Dynamic, 1 > &clusterSize)
 Compute size of each cluster given a partitioning.
 
template<typename EivalsType , typename ListOfClusters , typename VectorType >
void matrix_function_compute_map (const EivalsType &eivals, const ListOfClusters &clusters, VectorType &eivalToCluster)
 Compute mapping of eigenvalue indices to cluster indices.
 
template<typename DynVectorType , typename VectorType >
void matrix_function_compute_permutation (const DynVectorType &blockStart, const DynVectorType &eivalToCluster, VectorType &permutation)
 Compute permutation which groups ei'vals in same cluster together.
 
template<typename Index , typename ListOfClusters >
ListOfClusters::iterator matrix_function_find_cluster (Index key, ListOfClusters &clusters)
 Find cluster in clusters containing some value. More...
 
template<typename EivalsType , typename Cluster >
void matrix_function_partition_eigenvalues (const EivalsType &eivals, std::list< Cluster > &clusters)
 Partition eigenvalues in clusters of ei'vals close to each other. More...
 
template<typename VectorType , typename MatrixType >
void matrix_function_permute_schur (VectorType &permutation, MatrixType &U, MatrixType &T)
 Permute Schur decomposition in U and T according to permutation.
 
template<typename MatrixType >
MatrixType matrix_function_solve_triangular_sylvester (const MatrixType &A, const MatrixType &B, const MatrixType &C)
 Solve a triangular Sylvester equation AX + XB = C. More...
 
template<typename MatrixType >
void matrix_log_compute_2x2 (const MatrixType &A, MatrixType &result)
 Compute logarithm of 2x2 triangular matrix.
 
template<typename MatrixType >
void matrix_log_compute_big (const MatrixType &A, MatrixType &result)
 Compute logarithm of triangular matrices with size > 2. More...
 
template<typename CMatrix , typename CINVMatrix >
void pseudo_inverse (const CMatrix &C, CINVMatrix &CINV)
 
template<typename VectorType , typename IndexType >
void sortWithPermutation (VectorType &vec, IndexType &perm, typename IndexType::Scalar &ncut)
 Computes a permutation vector to have a sorted sequence. More...
 
template<typename Scalar >
Scalar stem_function_cos (Scalar x, int n)
 Cosine (and its derivatives).
 
template<typename Scalar >
Scalar stem_function_cosh (Scalar x, int n)
 Hyperbolic cosine (and its derivatives).
 
template<typename Scalar >
Scalar stem_function_exp (Scalar x, int)
 The exponential function (and its derivatives).
 
template<typename Scalar >
Scalar stem_function_sin (Scalar x, int n)
 Sine (and its derivatives).
 
template<typename Scalar >
Scalar stem_function_sinh (Scalar x, int n)
 Hyperbolic sine (and its derivatives).
 

Detailed Description

template <class> class MakePointer_ is added to convert the host pointer to the device pointer. It is added due to the fact that for our device compiler T* is not allowed. If we wanted to use the same Evaluator functions we have to convert that type to our pointer T. This is done through our MakePointer_ class. By default the Type in the MakePointer_<T> is T* . Therefore, by adding the default value, we managed to convert the type and it does not break any existing code as its default value is T*.

Function Documentation

template<typename MatrixType , typename Rhs , typename Dest , typename Preconditioner >
bool Eigen::internal::gmres ( const MatrixType &  mat,
const Rhs &  rhs,
Dest &  x,
const Preconditioner &  precond,
Index &  iters,
const Index &  restart,
typename Dest::RealScalar &  tol_error 
)

Generalized Minimal Residual Algorithm based on the Arnoldi algorithm implemented with Householder reflections.

Parameters:

Parameters
matmatrix of linear system of equations
Rhsright hand side vector of linear system of equations
xon input: initial guess, on output: solution
precondpreconditioner used
iterson input: maximum number of iterations to perform on output: number of iterations performed
restartnumber of iterations for a restart
tol_erroron input: relative residual tolerance on output: residuum achieved
See Also
IterativeMethods::bicgstab()

For references, please see:

Saad, Y. and Schultz, M. H. GMRES: A Generalized Minimal Residual Algorithm for Solving Nonsymmetric Linear Systems. SIAM J.Sci.Stat.Comp. 7, 1986, pp. 856 - 869.

Saad, Y. Iterative Methods for Sparse Linear Systems. Society for Industrial and Applied Mathematics, Philadelphia, 2003.

Walker, H. F. Implementations of the GMRES method. Comput.Phys.Comm. 53, 1989, pp. 311 - 320.

Walker, H. F. Implementation of the GMRES Method using Householder Transformations. SIAM J.Sci.Stat.Comp. 9, 1988, pp. 152 - 163.

template<typename MatA , typename MatU , typename MatV >
void Eigen::internal::matrix_exp_pade13 ( const MatA &  A,
MatU &  U,
MatV &  V 
)

Compute the (13,13)-Padé approximant to the exponential.

After exit, $ (V+U)(V-U)^{-1} $ is the Padé approximant of $ \exp(A) $ around $ A = 0 $.

template<typename MatA , typename MatU , typename MatV >
void Eigen::internal::matrix_exp_pade3 ( const MatA &  A,
MatU &  U,
MatV &  V 
)

Compute the (3,3)-Padé approximant to the exponential.

After exit, $ (V+U)(V-U)^{-1} $ is the Padé approximant of $ \exp(A) $ around $ A = 0 $.

template<typename MatA , typename MatU , typename MatV >
void Eigen::internal::matrix_exp_pade5 ( const MatA &  A,
MatU &  U,
MatV &  V 
)

Compute the (5,5)-Padé approximant to the exponential.

After exit, $ (V+U)(V-U)^{-1} $ is the Padé approximant of $ \exp(A) $ around $ A = 0 $.

template<typename MatA , typename MatU , typename MatV >
void Eigen::internal::matrix_exp_pade7 ( const MatA &  A,
MatU &  U,
MatV &  V 
)

Compute the (7,7)-Padé approximant to the exponential.

After exit, $ (V+U)(V-U)^{-1} $ is the Padé approximant of $ \exp(A) $ around $ A = 0 $.

template<typename MatA , typename MatU , typename MatV >
void Eigen::internal::matrix_exp_pade9 ( const MatA &  A,
MatU &  U,
MatV &  V 
)

Compute the (9,9)-Padé approximant to the exponential.

After exit, $ (V+U)(V-U)^{-1} $ is the Padé approximant of $ \exp(A) $ around $ A = 0 $.

template<typename MatrixType , typename VectorType >
void Eigen::internal::matrix_function_compute_above_diagonal ( const MatrixType &  T,
const VectorType &  blockStart,
const VectorType &  clusterSize,
MatrixType &  fT 
)

Compute part of matrix function above block diagonal.

This routine completes the computation of fT, denoting a matrix function applied to the triangular matrix T. It assumes that the block diagonal part of fT has already been computed. The part below the diagonal is zero, because T is upper triangular.

template<typename MatrixType , typename AtomicType , typename VectorType >
void Eigen::internal::matrix_function_compute_block_atomic ( const MatrixType &  T,
AtomicType &  atomic,
const VectorType &  blockStart,
const VectorType &  clusterSize,
MatrixType &  fT 
)

Compute block diagonal part of matrix function.

This routine computes the matrix function applied to the block diagonal part of T (which should be upper triangular), with the blocking given by blockStart and clusterSize. The matrix function of each diagonal block is computed by atomic. The off-diagonal parts of fT are set to zero.

template<typename Index , typename ListOfClusters >
ListOfClusters::iterator Eigen::internal::matrix_function_find_cluster ( Index  key,
ListOfClusters &  clusters 
)

Find cluster in clusters containing some value.

Parameters
[in]keyValue to find
Returns
Iterator to cluster containing key, or clusters.end() if no cluster in m_clusters contains key.
template<typename EivalsType , typename Cluster >
void Eigen::internal::matrix_function_partition_eigenvalues ( const EivalsType &  eivals,
std::list< Cluster > &  clusters 
)

Partition eigenvalues in clusters of ei'vals close to each other.

Parameters
[in]eivalsEigenvalues
[out]clustersResulting partition of eigenvalues

The partition satisfies the following two properties:

Any eigenvalue in a certain cluster is at most matrix_function_separation() away from another eigenvalue

in the same cluster.

The distance between two eigenvalues in different clusters is more than matrix_function_separation().

The implementation follows Algorithm 4.1 in the paper of Davies and Higham.

template<typename MatrixType >
MatrixType Eigen::internal::matrix_function_solve_triangular_sylvester ( const MatrixType &  A,
const MatrixType &  B,
const MatrixType &  C 
)

Solve a triangular Sylvester equation AX + XB = C.

Parameters
[in]Athe matrix A; should be square and upper triangular
[in]Bthe matrix B; should be square and upper triangular
[in]Cthe matrix C; should have correct size.
Returns
the solution X.

If A is m-by-m and B is n-by-n, then both C and X are m-by-n. The (i,j)-th component of the Sylvester equation is

\[ \sum_{k=i}^m A_{ik} X_{kj} + \sum_{k=1}^j X_{ik} B_{kj} = C_{ij}. \]

This can be re-arranged to yield:

\[ X_{ij} = \frac{1}{A_{ii} + B_{jj}} \Bigl( C_{ij} - \sum_{k=i+1}^m A_{ik} X_{kj} - \sum_{k=1}^{j-1} X_{ik} B_{kj} \Bigr). \]

It is assumed that A and B are such that the numerator is never zero (otherwise the Sylvester equation does not have a unique solution). In that case, these equations can be evaluated in the order $ i=m,\ldots,1 $ and $ j=1,\ldots,n $.

template<typename MatrixType >
void Eigen::internal::matrix_log_compute_big ( const MatrixType &  A,
MatrixType &  result 
)

Compute logarithm of triangular matrices with size > 2.

This uses a inverse scale-and-square algorithm.

template<typename VectorType , typename IndexType >
void Eigen::internal::sortWithPermutation ( VectorType &  vec,
IndexType &  perm,
typename IndexType::Scalar &  ncut 
)

Computes a permutation vector to have a sorted sequence.

Parameters
vecThe vector to reorder.
permgives the sorted sequence on output. Must be initialized with 0..n-1
ncutPut the ncut smallest elements at the end of the vector WARNING This is an expensive sort, so should be used only for small size vectors TODO Use modified QuickSplit or std::nth_element to get the smallest values