SOFA API  c76874b7
Open source framework for multi-physics simuation
sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector > Class Template Reference

#include <CGLinearSolver.h>

Linear system solver using the conjugate gradient iterative algorithm. More...

Inheritance diagram for sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >:

Detailed Description

template<class TMatrix, class TVector>
class sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >

Linear system solver using the conjugate gradient iterative algorithm.

Public Attributes

Data< unsigned > f_maxIter
 maximum number of iterations of the Conjugate Gradient solution More...
 
Data< SReal > f_tolerance
 desired precision of the Conjugate Gradient Solution (ratio of current residual norm over initial residual norm) More...
 
Data< SReal > f_smallDenominatorThreshold
 minimum value of the denominator in the conjugate Gradient solution More...
 
Data< boolf_warmStart
 Use previous solution as initial solution. More...
 
Data< boolf_verbose
 Dump system state at each iteration. More...
 
Data< std::map< std::string, sofa::helper::vector< SReal > > > f_graph
 Graph of residuals at each iteration. More...
 

Protected Attributes

int timeStepCount
 
bool equilibriumReached
 

Public Member Functions

 SOFA_CLASS (SOFA_TEMPLATE2(CGLinearSolver, TMatrix, TVector), SOFA_TEMPLATE2(sofa::component::linearsolver::MatrixLinearSolver, TMatrix, TVector))
 
void init () override
 
void reinit () override
 
void resetSystem () override
 
void setSystemMBKMatrix (const sofa::core::MechanicalParams *mparams) override
 
void solve (Matrix &M, Vector &x, Vector &b) override
 Solve Mx=b. More...
 

Protected Member Functions

 CGLinearSolver ()
 Linear system solver using the conjugate gradient iterative algorithm. More...
 
void cgstep_beta (const core::ExecParams *params, Vector &p, Vector &r, SReal beta)
 
void cgstep_alpha (const core::ExecParams *params, Vector &x, Vector &r, Vector &p, Vector &q, SReal alpha)
 
template<>
SOFA_BASE_LINEAR_SOLVER_API void cgstep_beta (const core::ExecParams *, Vector &p, Vector &r, SReal beta)
 
template<>
SOFA_BASE_LINEAR_SOLVER_API void cgstep_alpha (const core::ExecParams *params, Vector &x, Vector &r, Vector &p, Vector &q, SReal alpha)
 
template<>
void cgstep_beta (const core::ExecParams *, Vector &p, Vector &r, SReal beta)
 
template<>
void cgstep_alpha (const core::ExecParams *params, Vector &x, Vector &r, Vector &p, Vector &q, SReal alpha)
 

Attribute details

template<class TMatrix , class TVector >
bool sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::equilibriumReached
protected
template<class TMatrix , class TVector >
Data<std::map < std::string, sofa::helper::vector<SReal> > > sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::f_graph

Graph of residuals at each iteration.

template<class TMatrix , class TVector >
Data<unsigned> sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::f_maxIter

maximum number of iterations of the Conjugate Gradient solution

template<class TMatrix , class TVector >
Data<SReal> sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::f_smallDenominatorThreshold

minimum value of the denominator in the conjugate Gradient solution

template<class TMatrix , class TVector >
Data<SReal> sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::f_tolerance

desired precision of the Conjugate Gradient Solution (ratio of current residual norm over initial residual norm)

template<class TMatrix , class TVector >
Data<bool> sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::f_verbose

Dump system state at each iteration.

template<class TMatrix , class TVector >
Data<bool> sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::f_warmStart

Use previous solution as initial solution.

template<class TMatrix , class TVector >
int sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::timeStepCount
protected

Constructor details

template<class TMatrix , class TVector >
sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::CGLinearSolver ( )
protected

Linear system solver using the conjugate gradient iterative algorithm.

Function details

template<>
SOFA_BASE_LINEAR_SOLVER_API void sofa::component::linearsolver::CGLinearSolver< component::linearsolver::GraphScatteredMatrix, component::linearsolver::GraphScatteredVector >::cgstep_alpha ( const core::ExecParams params,
Vector x,
Vector r,
Vector p,
Vector q,
SReal  alpha 
)
inlineprotected
template<class TMatrix , class TVector >
void sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::cgstep_alpha ( const core::ExecParams params,
Vector x,
Vector r,
Vector p,
Vector q,
SReal  alpha 
)
inlineprotected

This method is separated from the rest to be able to use custom/optimized versions depending on the types of vectors. It computes: x += p*alpha, r -= q*alpha

template<>
void sofa::component::linearsolver::CGLinearSolver< component::linearsolver::GraphScatteredMatrix, component::linearsolver::GraphScatteredVector >::cgstep_alpha ( const core::ExecParams params,
Vector x,
Vector r,
Vector p,
Vector q,
SReal  alpha 
)
inlineprotected
template<>
SOFA_BASE_LINEAR_SOLVER_API void sofa::component::linearsolver::CGLinearSolver< component::linearsolver::GraphScatteredMatrix, component::linearsolver::GraphScatteredVector >::cgstep_beta ( const core::ExecParams ,
Vector p,
Vector r,
SReal  beta 
)
inlineprotected
template<class TMatrix , class TVector >
void sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::cgstep_beta ( const core::ExecParams params,
Vector p,
Vector r,
SReal  beta 
)
inlineprotected

This method is separated from the rest to be able to use custom/optimized versions depending on the types of vectors. It computes: p = p*beta + r

template<class TMatrix , class TVector >
void sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::init ( )
override
template<class TMatrix , class TVector >
void sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::reinit ( )
override
template<class TMatrix , class TVector >
void sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::resetSystem ( )
override
template<class TMatrix , class TVector >
void sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::setSystemMBKMatrix ( const sofa::core::MechanicalParams mparams)
override
template<class TMatrix , class TVector >
sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::SOFA_CLASS ( SOFA_TEMPLATE2(CGLinearSolver< TMatrix, TVector >, TMatrix, TVector)  ,
SOFA_TEMPLATE2(sofa::component::linearsolver::MatrixLinearSolver, TMatrix, TVector)   
)
template<class TMatrix , class TVector >
void sofa::component::linearsolver::CGLinearSolver< TMatrix, TVector >::solve ( Matrix M,
Vector x,
Vector b 
)
override

Solve Mx=b.