![]() |
Eigen-unsupported
3.3.4
|
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). | |
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*.
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:
mat | matrix of linear system of equations |
Rhs | right hand side vector of linear system of equations |
x | on input: initial guess, on output: solution |
precond | preconditioner used |
iters | on input: maximum number of iterations to perform on output: number of iterations performed |
restart | number of iterations for a restart |
tol_error | on input: relative residual tolerance on output: residuum achieved |
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.
void Eigen::internal::matrix_exp_pade13 | ( | const MatA & | A, |
MatU & | U, | ||
MatV & | V | ||
) |
Compute the (13,13)-Padé approximant to the exponential.
After exit, is the Padé approximant of
around
.
void Eigen::internal::matrix_exp_pade3 | ( | const MatA & | A, |
MatU & | U, | ||
MatV & | V | ||
) |
Compute the (3,3)-Padé approximant to the exponential.
After exit, is the Padé approximant of
around
.
void Eigen::internal::matrix_exp_pade5 | ( | const MatA & | A, |
MatU & | U, | ||
MatV & | V | ||
) |
Compute the (5,5)-Padé approximant to the exponential.
After exit, is the Padé approximant of
around
.
void Eigen::internal::matrix_exp_pade7 | ( | const MatA & | A, |
MatU & | U, | ||
MatV & | V | ||
) |
Compute the (7,7)-Padé approximant to the exponential.
After exit, is the Padé approximant of
around
.
void Eigen::internal::matrix_exp_pade9 | ( | const MatA & | A, |
MatU & | U, | ||
MatV & | V | ||
) |
Compute the (9,9)-Padé approximant to the exponential.
After exit, is the Padé approximant of
around
.
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.
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.
ListOfClusters::iterator Eigen::internal::matrix_function_find_cluster | ( | Index | key, |
ListOfClusters & | clusters | ||
) |
Find cluster in clusters
containing some value.
[in] | key | Value to find |
key
, or clusters.end()
if no cluster in m_clusters
contains key
. 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.
[in] | eivals | Eigenvalues |
[out] | clusters | Resulting partition of eigenvalues |
The partition satisfies the following two properties:
in the same cluster.
The implementation follows Algorithm 4.1 in the paper of Davies and Higham.
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.
[in] | A | the matrix A; should be square and upper triangular |
[in] | B | the matrix B; should be square and upper triangular |
[in] | C | the matrix C; should have correct size. |
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
This can be re-arranged to yield:
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 and
.
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.
void Eigen::internal::sortWithPermutation | ( | VectorType & | vec, |
IndexType & | perm, | ||
typename IndexType::Scalar & | ncut | ||
) |
Computes a permutation vector to have a sorted sequence.
vec | The vector to reorder. |
perm | gives the sorted sequence on output. Must be initialized with 0..n-1 |
ncut | Put 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 |