Caribou
Public Types | Public Member Functions | Static Public Member Functions | List of all members
SofaCaribou::solver::EigenSparseSolver< EigenSolver_t > Class Template Reference

Detailed Description

template<class EigenSolver_t>
class SofaCaribou::solver::EigenSparseSolver< EigenSolver_t >

Base class for sparse direct solvers using Eigen as a backend solver.

Note that the complete system matrix has to be assemble for any of the derivated solvers. This means that the mass, damping and forcefield components in the scene graph need to implement the addMtoMatrix, addBtoMatrix and addMtoMatrix methods respectively.

Template Parameters
EigenSolver_tEigen solver type (eg.: SimplicialLLT, SparseLU, PardisoLLT, etc.)

#include <EigenSparseSolver.h>

Inheritance diagram for SofaCaribou::solver::EigenSparseSolver< EigenSolver_t >:
SofaCaribou::solver::LinearSolver

Public Types

using EigenSolver = EigenSolver_t
 
using SparseMatrix = std::remove_cv_t< typename EigenSolver::MatrixType >
 
using Vector = Eigen::Matrix< FLOATING_POINT_TYPE, Eigen::Dynamic, 1 >
 

Public Member Functions

 SOFA_CLASS (SOFA_TEMPLATE(EigenSparseSolver, EigenSolver_t), sofa::core::behavior::LinearSolver)
 
virtual void assemble (const sofa::core::MechanicalParams *mparams)
 Assemble the system matrix A = (mM + bB + kK) inside the SparseMatrix p_A. More...
 
void resetSystem () final
 Reset the complete system (A, x and b are cleared). More...
 
void setSystemMBKMatrix (const sofa::core::MechanicalParams *mparams) final
 Set the linear system matrix A = (mM + bB + kK), storing the coefficients m, b and k of the mechanical M,B,K matrices. More...
 
void setSystemRHVector (sofa::core::MultiVecDerivId b_id) final
 Gives the identifier of the right-hand side vector b. More...
 
void setSystemLHVector (sofa::core::MultiVecDerivId x_id) final
 Gives the identifier of the left-hand side vector x. More...
 
void solveSystem () override
 Solves the system using the Eigen solver.
 
virtual auto symmetric () const -> bool
 States if the system matrix is symmetric. More...
 
virtual void set_symmetric (bool is_symmetric)
 Explicitly states if this matrix is symmetric.
 
auto solver () const -> const EigenSolver &
 Get a readonly reference to the backend solver.
 
auto solver () -> EigenSolver &
 Get a reference to the backend solver.
 
auto mechanical_params () const -> const sofa::core::MechanicalParams &
 Get a readonly reference to the mechanical parameters.
 
auto matrix_accessor () const -> const sofa::component::linearsolver::DefaultMultiMatrixAccessor &
 Get a readonly reference to the multi-matrix accessor. More...
 
auto A () const -> const SparseMatrix &
 Get a readonly reference to the global assembled system matrix.
 
auto b_id () const -> const sofa::core::MultiVecDerivId &
 Get a readonly reference to the right-hand side vector identifier.
 
auto x_id () const -> const sofa::core::MultiVecDerivId &
 Get a readonly reference to the left-hand side unknown vector identifier.
 
auto A_is_factorized () const -> bool
 True if the solver has successfully factorize the system matrix.
 
- Public Member Functions inherited from SofaCaribou::solver::LinearSolver
virtual sofa::defaulttype::BaseMatrix * create_new_matrix (unsigned int rows, unsigned int cols) const =0
 Creates a new BaseMatrix of size rows x cols. More...
 
virtual sofa::defaulttype::BaseVector * create_new_vector (unsigned int n) const =0
 Creates a new BaseVector of size n. More...
 

Static Public Member Functions

static auto GetCustomTemplateName () -> std::string
 
template<class Derived >
static auto canCreate (Derived *o, sofa::core::objectmodel::BaseContext *context, sofa::core::objectmodel::BaseObjectDescription *arg) -> bool
 

Member Function Documentation

◆ assemble()

template<class EigenSolver >
void SofaCaribou::solver::EigenSparseSolver< EigenSolver >::assemble ( const sofa::core::MechanicalParams *  mparams)
virtual

Assemble the system matrix A = (mM + bB + kK) inside the SparseMatrix p_A.

Parameters
mparamsMechanical parameters containing the m, b and k factors.

Reimplemented in SofaCaribou::solver::LUSolver< EigenSolver >.

◆ matrix_accessor()

template<class EigenSolver_t >
auto SofaCaribou::solver::EigenSparseSolver< EigenSolver_t >::matrix_accessor ( ) const -> const sofa::component::linearsolver::DefaultMultiMatrixAccessor &
inline

Get a readonly reference to the multi-matrix accessor.

This accessor doesn't actually hold the system matrix, but will be the one responsible for accumulating it. It can also be use to compute the stiffness or the mass matrix of a given mechanical node or mapped node.

Use the EigenSparseSolver::matrix() method to get a reference to the global system matrix.

◆ resetSystem()

template<class EigenSolver >
void SofaCaribou::solver::EigenSparseSolver< EigenSolver >::resetSystem
final

Reset the complete system (A, x and b are cleared).

When using no preconditioner (None), this does absolutely nothing here since the complete system is never built.

This method is called by the MultiMatrix::clear() and MultiMatrix::reset() methods.

◆ setSystemLHVector()

template<class EigenSolver >
void SofaCaribou::solver::EigenSparseSolver< EigenSolver >::setSystemLHVector ( sofa::core::MultiVecDerivId  x_id)
final

Gives the identifier of the left-hand side vector x.

This identifier will be used to find the actual vector in the mechanical objects of the system. The complete dense vector is accumulated from the mechanical objects found in the graph subtree of the current context.

◆ setSystemMBKMatrix()

template<class EigenSolver >
void SofaCaribou::solver::EigenSparseSolver< EigenSolver >::setSystemMBKMatrix ( const sofa::core::MechanicalParams *  mparams)
final

Set the linear system matrix A = (mM + bB + kK), storing the coefficients m, b and k of the mechanical M,B,K matrices.

Parameters
mparamsContains the coefficients m, b and k of the matrices M, B and K

◆ setSystemRHVector()

template<class EigenSolver >
void SofaCaribou::solver::EigenSparseSolver< EigenSolver >::setSystemRHVector ( sofa::core::MultiVecDerivId  b_id)
final

Gives the identifier of the right-hand side vector b.

This identifier will be used to find the actual vector in the mechanical objects of the system. The complete dense vector is accumulated from the mechanical objects found in the graph subtree of the current context.

◆ symmetric()

template<class EigenSolver_t >
virtual auto SofaCaribou::solver::EigenSparseSolver< EigenSolver_t >::symmetric ( ) const -> bool
inlinevirtual

States if the system matrix is symmetric.

Note that this value isn't set automatically, the user must explicitly specify it using set_symmetric(true). When it is true, some optimizations will be enabled.


The documentation for this class was generated from the following files: