CGLinearSolver
This component belongs to the category of LinearSolver. The role of the CGLinearSolver is to solve the linear system \(\mathbf{A}x=b\) without any a priori on this system.
In SOFA, the CGLinearSolver follows the well-known conjugate gradient method, which consists in iteratively solving \(r=b-\mathbf{A}x^k\) where r is known as the residual. This residual will be used to compute mutually conjugate vectors p (see the sequence diagram below) which will be used as a basis to find a new approximated solution \(x^{k+1}\).
Note: the CGLinearSolver in SOFA assumes that the right hand side (RHS) vector b is already computed. The computation of b is usually called in the integration scheme through the function computeForce().
Sequence diagram
Usage
The CGLinearSolver requires the use (above in the scene graph) of an integration scheme, and (below in the scene graph) of a MechanicalObject storing the state information that the CGLinearSolver will access.
When using a CGLinearSolver, make sure you carefully chose the value of the free data field iterations, tolerance and threshold. Both tolerance and threshold data must be chosen in accordance with the dimension of the degrees of freedom (DOFs). Usually, the value of these two data is close to the square of the expected error on the DOFs.
Remember that when using an iterative linear solver like the CGLinearSolver, no exact solution can be found. The accuracy of your solution will always depend on the conditioning of your system and your input data (iterations, tolerance and threshold).
Linear system solver using the conjugate gradient iterative algorithm
CompressedRowSparseMatrixMat2x2d
Templates:
- CompressedRowSparseMatrixMat2x2d
Target: Sofa.Component.LinearSolver.Iterative
namespace: sofa::component::linearsolver::iterative
parents:
- MatrixLinearSolver
Data
| Name | Description | Default value |
|---|---|---|
| name | object name | unnamed |
| printLog | if true, emits extra messages at runtime. | 0 |
| tags | list of the subsets the objet belongs to | |
| bbox | this object bounding box | |
| componentState | The state of the component among (Dirty, Valid, Undefined, Loading, Invalid). | Undefined |
| listening | if true, handle the events, otherwise ignore the events | 0 |
| parallelInverseProduct | Parallelize the computation of the product J*M^{-1}*J^T where M is the matrix of the linear system and J is any matrix with compatible dimensions | 0 |
| iterations | Maximum number of iterations after which the iterative descent of the Conjugate Gradient must stop | 25 |
| tolerance | Desired accuracy of the Conjugate Gradient solution evaluating: |r|²/|b|² (ratio of current residual norm over initial residual norm) | 1e-05 |
| threshold | Minimum value of the denominator (pT A p)^ in the conjugate Gradient solution | 1e-05 |
| warmStart | Use previous solution as initial solution, which may improve the initial guess if your system is evolving smoothly | 0 |
| graph | Graph of residuals at each iteration |
Links
| Name | Description | Destination type name |
|---|---|---|
| context | Graph Node containing this object (or BaseContext::getDefault() if no graph is used) | BaseContext |
| slaves | Sub-objects used internally by this object | BaseObject |
| master | nullptr for regular objects, or master object for which this object is one sub-objects | BaseObject |
| linearSystem | The linear system to solve | TypedMatrixLinearSystem<CompressedRowSparseMatrixMat2x2d> |
CompressedRowSparseMatrixMat3x3d
Templates:
- CompressedRowSparseMatrixMat3x3d
Target: Sofa.Component.LinearSolver.Iterative
namespace: sofa::component::linearsolver::iterative
parents:
- MatrixLinearSolver
Data
| Name | Description | Default value |
|---|---|---|
| name | object name | unnamed |
| printLog | if true, emits extra messages at runtime. | 0 |
| tags | list of the subsets the objet belongs to | |
| bbox | this object bounding box | |
| componentState | The state of the component among (Dirty, Valid, Undefined, Loading, Invalid). | Undefined |
| listening | if true, handle the events, otherwise ignore the events | 0 |
| parallelInverseProduct | Parallelize the computation of the product J*M^{-1}*J^T where M is the matrix of the linear system and J is any matrix with compatible dimensions | 0 |
| iterations | Maximum number of iterations after which the iterative descent of the Conjugate Gradient must stop | 25 |
| tolerance | Desired accuracy of the Conjugate Gradient solution evaluating: |r|²/|b|² (ratio of current residual norm over initial residual norm) | 1e-05 |
| threshold | Minimum value of the denominator (pT A p)^ in the conjugate Gradient solution | 1e-05 |
| warmStart | Use previous solution as initial solution, which may improve the initial guess if your system is evolving smoothly | 0 |
| graph | Graph of residuals at each iteration |
Links
| Name | Description | Destination type name |
|---|---|---|
| context | Graph Node containing this object (or BaseContext::getDefault() if no graph is used) | BaseContext |
| slaves | Sub-objects used internally by this object | BaseObject |
| master | nullptr for regular objects, or master object for which this object is one sub-objects | BaseObject |
| linearSystem | The linear system to solve | TypedMatrixLinearSystem<CompressedRowSparseMatrixMat3x3d> |
CompressedRowSparseMatrixMat4x4d
Templates:
- CompressedRowSparseMatrixMat4x4d
Target: Sofa.Component.LinearSolver.Iterative
namespace: sofa::component::linearsolver::iterative
parents:
- MatrixLinearSolver
Data
| Name | Description | Default value |
|---|---|---|
| name | object name | unnamed |
| printLog | if true, emits extra messages at runtime. | 0 |
| tags | list of the subsets the objet belongs to | |
| bbox | this object bounding box | |
| componentState | The state of the component among (Dirty, Valid, Undefined, Loading, Invalid). | Undefined |
| listening | if true, handle the events, otherwise ignore the events | 0 |
| parallelInverseProduct | Parallelize the computation of the product J*M^{-1}*J^T where M is the matrix of the linear system and J is any matrix with compatible dimensions | 0 |
| iterations | Maximum number of iterations after which the iterative descent of the Conjugate Gradient must stop | 25 |
| tolerance | Desired accuracy of the Conjugate Gradient solution evaluating: |r|²/|b|² (ratio of current residual norm over initial residual norm) | 1e-05 |
| threshold | Minimum value of the denominator (pT A p)^ in the conjugate Gradient solution | 1e-05 |
| warmStart | Use previous solution as initial solution, which may improve the initial guess if your system is evolving smoothly | 0 |
| graph | Graph of residuals at each iteration |
Links
| Name | Description | Destination type name |
|---|---|---|
| context | Graph Node containing this object (or BaseContext::getDefault() if no graph is used) | BaseContext |
| slaves | Sub-objects used internally by this object | BaseObject |
| master | nullptr for regular objects, or master object for which this object is one sub-objects | BaseObject |
| linearSystem | The linear system to solve | TypedMatrixLinearSystem<CompressedRowSparseMatrixMat4x4d> |
CompressedRowSparseMatrixMat6x6d
Templates:
- CompressedRowSparseMatrixMat6x6d
Target: Sofa.Component.LinearSolver.Iterative
namespace: sofa::component::linearsolver::iterative
parents:
- MatrixLinearSolver
Data
| Name | Description | Default value |
|---|---|---|
| name | object name | unnamed |
| printLog | if true, emits extra messages at runtime. | 0 |
| tags | list of the subsets the objet belongs to | |
| bbox | this object bounding box | |
| componentState | The state of the component among (Dirty, Valid, Undefined, Loading, Invalid). | Undefined |
| listening | if true, handle the events, otherwise ignore the events | 0 |
| parallelInverseProduct | Parallelize the computation of the product J*M^{-1}*J^T where M is the matrix of the linear system and J is any matrix with compatible dimensions | 0 |
| iterations | Maximum number of iterations after which the iterative descent of the Conjugate Gradient must stop | 25 |
| tolerance | Desired accuracy of the Conjugate Gradient solution evaluating: |r|²/|b|² (ratio of current residual norm over initial residual norm) | 1e-05 |
| threshold | Minimum value of the denominator (pT A p)^ in the conjugate Gradient solution | 1e-05 |
| warmStart | Use previous solution as initial solution, which may improve the initial guess if your system is evolving smoothly | 0 |
| graph | Graph of residuals at each iteration |
Links
| Name | Description | Destination type name |
|---|---|---|
| context | Graph Node containing this object (or BaseContext::getDefault() if no graph is used) | BaseContext |
| slaves | Sub-objects used internally by this object | BaseObject |
| master | nullptr for regular objects, or master object for which this object is one sub-objects | BaseObject |
| linearSystem | The linear system to solve | TypedMatrixLinearSystem<CompressedRowSparseMatrixMat6x6d> |
CompressedRowSparseMatrixMat8x8d
Templates:
- CompressedRowSparseMatrixMat8x8d
Target: Sofa.Component.LinearSolver.Iterative
namespace: sofa::component::linearsolver::iterative
parents:
- MatrixLinearSolver
Data
| Name | Description | Default value |
|---|---|---|
| name | object name | unnamed |
| printLog | if true, emits extra messages at runtime. | 0 |
| tags | list of the subsets the objet belongs to | |
| bbox | this object bounding box | |
| componentState | The state of the component among (Dirty, Valid, Undefined, Loading, Invalid). | Undefined |
| listening | if true, handle the events, otherwise ignore the events | 0 |
| parallelInverseProduct | Parallelize the computation of the product J*M^{-1}*J^T where M is the matrix of the linear system and J is any matrix with compatible dimensions | 0 |
| iterations | Maximum number of iterations after which the iterative descent of the Conjugate Gradient must stop | 25 |
| tolerance | Desired accuracy of the Conjugate Gradient solution evaluating: |r|²/|b|² (ratio of current residual norm over initial residual norm) | 1e-05 |
| threshold | Minimum value of the denominator (pT A p)^ in the conjugate Gradient solution | 1e-05 |
| warmStart | Use previous solution as initial solution, which may improve the initial guess if your system is evolving smoothly | 0 |
| graph | Graph of residuals at each iteration |
Links
| Name | Description | Destination type name |
|---|---|---|
| context | Graph Node containing this object (or BaseContext::getDefault() if no graph is used) | BaseContext |
| slaves | Sub-objects used internally by this object | BaseObject |
| master | nullptr for regular objects, or master object for which this object is one sub-objects | BaseObject |
| linearSystem | The linear system to solve | TypedMatrixLinearSystem<CompressedRowSparseMatrixMat8x8d> |
CompressedRowSparseMatrixd
Templates:
- CompressedRowSparseMatrixd
Target: Sofa.Component.LinearSolver.Iterative
namespace: sofa::component::linearsolver::iterative
parents:
- MatrixLinearSolver
Data
| Name | Description | Default value |
|---|---|---|
| name | object name | unnamed |
| printLog | if true, emits extra messages at runtime. | 0 |
| tags | list of the subsets the objet belongs to | |
| bbox | this object bounding box | |
| componentState | The state of the component among (Dirty, Valid, Undefined, Loading, Invalid). | Undefined |
| listening | if true, handle the events, otherwise ignore the events | 0 |
| parallelInverseProduct | Parallelize the computation of the product J*M^{-1}*J^T where M is the matrix of the linear system and J is any matrix with compatible dimensions | 0 |
| iterations | Maximum number of iterations after which the iterative descent of the Conjugate Gradient must stop | 25 |
| tolerance | Desired accuracy of the Conjugate Gradient solution evaluating: |r|²/|b|² (ratio of current residual norm over initial residual norm) | 1e-05 |
| threshold | Minimum value of the denominator (pT A p)^ in the conjugate Gradient solution | 1e-05 |
| warmStart | Use previous solution as initial solution, which may improve the initial guess if your system is evolving smoothly | 0 |
| graph | Graph of residuals at each iteration |
Links
| Name | Description | Destination type name |
|---|---|---|
| context | Graph Node containing this object (or BaseContext::getDefault() if no graph is used) | BaseContext |
| slaves | Sub-objects used internally by this object | BaseObject |
| master | nullptr for regular objects, or master object for which this object is one sub-objects | BaseObject |
| linearSystem | The linear system to solve | TypedMatrixLinearSystem<CompressedRowSparseMatrixd> |
FullMatrix
Templates:
- FullMatrix
Target: Sofa.Component.LinearSolver.Iterative
namespace: sofa::component::linearsolver::iterative
parents:
- MatrixLinearSolver
Data
| Name | Description | Default value |
|---|---|---|
| name | object name | unnamed |
| printLog | if true, emits extra messages at runtime. | 0 |
| tags | list of the subsets the objet belongs to | |
| bbox | this object bounding box | |
| componentState | The state of the component among (Dirty, Valid, Undefined, Loading, Invalid). | Undefined |
| listening | if true, handle the events, otherwise ignore the events | 0 |
| parallelInverseProduct | Parallelize the computation of the product J*M^{-1}*J^T where M is the matrix of the linear system and J is any matrix with compatible dimensions | 0 |
| iterations | Maximum number of iterations after which the iterative descent of the Conjugate Gradient must stop | 25 |
| tolerance | Desired accuracy of the Conjugate Gradient solution evaluating: |r|²/|b|² (ratio of current residual norm over initial residual norm) | 1e-05 |
| threshold | Minimum value of the denominator (pT A p)^ in the conjugate Gradient solution | 1e-05 |
| warmStart | Use previous solution as initial solution, which may improve the initial guess if your system is evolving smoothly | 0 |
| graph | Graph of residuals at each iteration |
Links
| Name | Description | Destination type name |
|---|---|---|
| context | Graph Node containing this object (or BaseContext::getDefault() if no graph is used) | BaseContext |
| slaves | Sub-objects used internally by this object | BaseObject |
| master | nullptr for regular objects, or master object for which this object is one sub-objects | BaseObject |
| linearSystem | The linear system to solve | TypedMatrixLinearSystem<FullMatrix> |
GraphScattered
Templates:
- GraphScattered
Target: Sofa.Component.LinearSolver.Iterative
namespace: sofa::component::linearsolver::iterative
parents:
- MatrixLinearSolver
Data
| Name | Description | Default value |
|---|---|---|
| name | object name | unnamed |
| printLog | if true, emits extra messages at runtime. | 0 |
| tags | list of the subsets the objet belongs to | |
| bbox | this object bounding box | |
| componentState | The state of the component among (Dirty, Valid, Undefined, Loading, Invalid). | Undefined |
| listening | if true, handle the events, otherwise ignore the events | 0 |
| parallelInverseProduct | Parallelize the computation of the product J*M^{-1}*J^T where M is the matrix of the linear system and J is any matrix with compatible dimensions | 0 |
| iterations | Maximum number of iterations after which the iterative descent of the Conjugate Gradient must stop | 25 |
| tolerance | Desired accuracy of the Conjugate Gradient solution evaluating: |r|²/|b|² (ratio of current residual norm over initial residual norm) | 1e-05 |
| threshold | Minimum value of the denominator (pT A p)^ in the conjugate Gradient solution | 1e-05 |
| warmStart | Use previous solution as initial solution, which may improve the initial guess if your system is evolving smoothly | 0 |
| graph | Graph of residuals at each iteration |
Links
| Name | Description | Destination type name |
|---|---|---|
| context | Graph Node containing this object (or BaseContext::getDefault() if no graph is used) | BaseContext |
| slaves | Sub-objects used internally by this object | BaseObject |
| master | nullptr for regular objects, or master object for which this object is one sub-objects | BaseObject |
| linearSystem | The linear system to solve | TypedMatrixLinearSystem<GraphScattered> |
SparseMatrix
Templates:
- SparseMatrix
Target: Sofa.Component.LinearSolver.Iterative
namespace: sofa::component::linearsolver::iterative
parents:
- MatrixLinearSolver
Data
| Name | Description | Default value |
|---|---|---|
| name | object name | unnamed |
| printLog | if true, emits extra messages at runtime. | 0 |
| tags | list of the subsets the objet belongs to | |
| bbox | this object bounding box | |
| componentState | The state of the component among (Dirty, Valid, Undefined, Loading, Invalid). | Undefined |
| listening | if true, handle the events, otherwise ignore the events | 0 |
| parallelInverseProduct | Parallelize the computation of the product J*M^{-1}*J^T where M is the matrix of the linear system and J is any matrix with compatible dimensions | 0 |
| iterations | Maximum number of iterations after which the iterative descent of the Conjugate Gradient must stop | 25 |
| tolerance | Desired accuracy of the Conjugate Gradient solution evaluating: |r|²/|b|² (ratio of current residual norm over initial residual norm) | 1e-05 |
| threshold | Minimum value of the denominator (pT A p)^ in the conjugate Gradient solution | 1e-05 |
| warmStart | Use previous solution as initial solution, which may improve the initial guess if your system is evolving smoothly | 0 |
| graph | Graph of residuals at each iteration |
Links
| Name | Description | Destination type name |
|---|---|---|
| context | Graph Node containing this object (or BaseContext::getDefault() if no graph is used) | BaseContext |
| slaves | Sub-objects used internally by this object | BaseObject |
| master | nullptr for regular objects, or master object for which this object is one sub-objects | BaseObject |
| linearSystem | The linear system to solve | TypedMatrixLinearSystem<SparseMatrix> |
Examples
CGLinearSolver.scn
<?xml version="1.0"?>
<Node name="root" dt="0.02" gravity="0 -10 0">
<RequiredPlugin name="Sofa.Component.Constraint.Projective"/> <!-- Needed to use components [FixedProjectiveConstraint] -->
<RequiredPlugin name="Sofa.Component.LinearSolver.Iterative"/> <!-- Needed to use components [CGLinearSolver] -->
<RequiredPlugin name="Sofa.Component.Mass"/> <!-- Needed to use components [UniformMass] -->
<RequiredPlugin name="Sofa.Component.ODESolver.Backward"/> <!-- Needed to use components [EulerImplicitSolver] -->
<RequiredPlugin name="Sofa.Component.SolidMechanics.FEM.Elastic"/> <!-- Needed to use components [HexahedronFEMForceField] -->
<RequiredPlugin name="Sofa.Component.StateContainer"/> <!-- Needed to use components [MechanicalObject] -->
<RequiredPlugin name="Sofa.Component.Topology.Container.Grid"/> <!-- Needed to use components [RegularGridTopology] -->
<RequiredPlugin name="Sofa.Component.Visual"/> <!-- Needed to use components [VisualStyle] -->
<VisualStyle displayFlags="showBehaviorModels showForceFields" />
<DefaultAnimationLoop/>
<Node>
<EulerImplicitSolver name="eulerimplicit_odesolver" printLog="false" rayleighStiffness="0.1" rayleighMass="0.1" />
<CGLinearSolver iterations="100" tolerance="1e-20" threshold="1e-20" warmStart="1" />
<MechanicalObject />
<UniformMass vertexMass="1" />
<RegularGridTopology nx="4" ny="4" nz="4" xmin="-9" xmax="-6" ymin="0" ymax="3" zmin="0" zmax="3" />
<FixedProjectiveConstraint indices="0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15" />
<HexahedronFEMForceField name="FEM" youngModulus="4000" poissonRatio="0.3" method="large" />
</Node>
</Node>
def createScene(root_node):
root = root_node.addChild('root', dt="0.02", gravity="0 -10 0")
root.addObject('RequiredPlugin', name="Sofa.Component.Constraint.Projective")
root.addObject('RequiredPlugin', name="Sofa.Component.LinearSolver.Iterative")
root.addObject('RequiredPlugin', name="Sofa.Component.Mass")
root.addObject('RequiredPlugin', name="Sofa.Component.ODESolver.Backward")
root.addObject('RequiredPlugin', name="Sofa.Component.SolidMechanics.FEM.Elastic")
root.addObject('RequiredPlugin', name="Sofa.Component.StateContainer")
root.addObject('RequiredPlugin', name="Sofa.Component.Topology.Container.Grid")
root.addObject('RequiredPlugin', name="Sofa.Component.Visual")
root.addObject('VisualStyle', displayFlags="showBehaviorModels showForceFields")
root.addObject('DefaultAnimationLoop', )
node = root.addChild('node')
node.addObject('EulerImplicitSolver', name="eulerimplicit_odesolver", printLog="false", rayleighStiffness="0.1", rayleighMass="0.1")
node.addObject('CGLinearSolver', iterations="100", tolerance="1e-20", threshold="1e-20", warmStart="1")
node.addObject('MechanicalObject', )
node.addObject('UniformMass', vertexMass="1")
node.addObject('RegularGridTopology', nx="4", ny="4", nz="4", xmin="-9", xmax="-6", ymin="0", ymax="3", zmin="0", zmax="3")
node.addObject('FixedProjectiveConstraint', indices="0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15")
node.addObject('HexahedronFEMForceField', name="FEM", youngModulus="4000", poissonRatio="0.3", method="large")