VariationalSymplecticSolver
Implicit time integrator which conserves linear momentum and mechanical energy. Implicit and Explicit time integrator using the Variational Symplectic Integrator as defined in : Kharevych, L et al. “Geometric, Variational Integrators for Computer Animation.” ACM SIGGRAPH Symposium on Computer Animation 4 (2006): 43–51. The current implementation for implicit integration assume alpha =0.5 (quadratic accuracy) and uses several Newton steps to estimate the velocity
The `VariationalSymplecticSolver` conserves linear momentum and mechanical energy using variational symplectic integrators for both implicit and explicit time integration schemes.
- module
- Sofa.Component.ODESolver.Backward
- namespace
- sofa::component::odesolver::backward
- include
- sofa/component/odesolver/backward/VariationalSymplecticSolver.h
- inherits
-
- LinearSolverAccessor
- OdeSolver
- description
The VariationalSymplecticSolver is an ODE solver within the SOFA framework designed for both implicit and explicit time integration using variational symplectic integrators as defined by Kharevych et al. (2006). This solver conserves linear momentum and mechanical energy, making it particularly suitable for simulations where such conservation laws are crucial.
Governing Equations
Implicit Integration Scheme
For implicit time integration with $\alpha = 0.5$ (quadratic accuracy), the VariationalSymplecticSolver solves the following equation iteratively using Newton steps:
where:
- $q_{k,i}$ is the estimate of the mid-point position at iteration $k$ and Newton step $i$.
- $
es_i$ is the position increment that minimizes the variational energy.
- $f(q_{k,i})$ represents internal forces, which include potential energy gradients.
- $K(q_{k,i})$ is the tangent stiffness matrix.
- $M$ is the mass matrix.
- $m{p}_k$ is momentum at time step $k$.
- $h$ is the time step size.
The solver computes $
es_i$, updates the position estimate, and checks for convergence based on the specified error tolerance (newtonError).
Explicit Integration Scheme
For explicit integration, the acceleration at each time step is computed as:
$$ \bm{a}_{k+1} = \frac{m{f}(q_k) + r_M M m{v}_k}{M} $$where $r_M$ is the Rayleigh damping coefficient related to mass.
The updated velocity and position are then:
and
$$ q_{k+1} = q_k + h m{v}_{k+1} $$This scheme does not require iterative solving and is less computationally expensive.
Constitutive Laws and Kinematic Models
The component uses Rayleigh damping to account for energy dissipation, characterized by the following damping forces:
$$ f_d = r_M M m{v} + r_K K(q) q $$where $r_K$ is the Rayleigh stiffness coefficient.
The solver also handles constraints and projections using mop.solveConstraint to ensure that the updated positions, velocities, and momenta satisfy any imposed constraints.
Role in FEM Pipeline
The VariationalSymplecticSolver fits into the broader variational mechanics framework by conserving key physical quantities such as linear momentum and mechanical energy. It operates within the time integration phase of the simulation pipeline:
1. Initialization: The solver initializes the system, including allocating memory for storing momentum.
2. Time Integration: The solver computes forces and accelerations, updating positions and velocities using either explicit or implicit schemes depending on user settings (d_explicit).
3. Nonlinear Solve: For implicit integration, it solves a nonlinear system iteratively using Newton steps to find the updated position that minimizes energy.
4. Energy Computation: The solver computes kinetic and potential energies at each time step for monitoring purposes or saving in a CSV file if d_saveEnergyInFile is enabled.
Numerical Methods and Discretization Choices
The component uses variational symplectic integrators, which are designed to preserve the Hamiltonian structure of mechanical systems. This ensures that energy and momentum are conserved over long-time simulations. The choice of $\alpha = 0.5$ for implicit integration provides quadratic accuracy in time step size.
The solver supports Rayleigh damping, which is a form of proportional damping where the damping matrix is a linear combination of the mass and stiffness matrices.
The VariationalSymplecticSolver also leverages linear solvers (via LinearSolverAccessor) to solve the resulting systems of equations efficiently.
Summary
The VariationalSymplecticSolver is designed for accurate time integration with strong conservation properties. It integrates seamlessly into SOFA's modular framework, providing a robust and efficient solution for simulating mechanical systems.
Data Fields
| Name | Type | Default | Help |
|---|---|---|---|
d_newtonError |
SReal | |
Error tolerance for Newton iterations |
d_newtonSteps |
unsigned int | |
Maximum number of Newton steps |
d_rayleighStiffness |
SReal | |
Rayleigh damping coefficient related to stiffness, > 0 |
d_rayleighMass |
SReal | |
Rayleigh damping coefficient related to mass, > 0 |
d_saveEnergyInFile |
bool | |
If kinetic and potential energies should be dumped in a CSV file at each iteration |
d_explicit |
bool | |
Use explicit integration scheme |
d_computeHamiltonian |
bool | |
Compute hamiltonian |
d_hamiltonianEnergy |
SReal | |
hamiltonian energy |
d_useIncrementalPotentialEnergy |
bool | |
use real potential energy, if false use approximate potential energy |
d_threadSafeVisitor |
bool | |
If true, do not use realloc and free visitors in fwdInteractionForceField. |
Methods
void
init
()
virtual
void
solve
(const core::ExecParams * params, SReal dt, sofa::core::MultiVecCoordId xResult, sofa::core::MultiVecDerivId vResult)
virtual
SReal
getVelocityIntegrationFactor
()
virtual
SReal
getPositionIntegrationFactor
()
virtual
SReal
getIntegrationFactor
(int , int )
virtual
SReal
getSolutionIntegrationFactor
(int )
virtual
void
parse
(core::objectmodel::BaseObjectDescription * arg)
virtual
{
"name": "VariationalSymplecticSolver",
"namespace": "sofa::component::odesolver::backward",
"module": "Sofa.Component.ODESolver.Backward",
"include": "sofa/component/odesolver/backward/VariationalSymplecticSolver.h",
"doc": "Implicit time integrator which conserves linear momentum and mechanical energy.\n\nImplicit and Explicit time integrator using the Variational Symplectic Integrator as defined in :\nKharevych, L et al. “Geometric, Variational Integrators for Computer Animation.” ACM SIGGRAPH Symposium on Computer Animation 4 (2006): 43–51.\nThe current implementation for implicit integration assume alpha =0.5 (quadratic accuracy) and uses\nseveral Newton steps to estimate the velocity",
"inherits": [
"LinearSolverAccessor",
"OdeSolver"
],
"templates": [],
"data_fields": [
{
"name": "d_newtonError",
"type": "SReal",
"xmlname": "newtonError",
"help": "Error tolerance for Newton iterations"
},
{
"name": "d_newtonSteps",
"type": "unsigned int",
"xmlname": "steps",
"help": "Maximum number of Newton steps"
},
{
"name": "d_rayleighStiffness",
"type": "SReal",
"xmlname": "rayleighStiffness",
"help": "Rayleigh damping coefficient related to stiffness, > 0"
},
{
"name": "d_rayleighMass",
"type": "SReal",
"xmlname": "rayleighMass",
"help": "Rayleigh damping coefficient related to mass, > 0"
},
{
"name": "d_saveEnergyInFile",
"type": "bool",
"xmlname": "saveEnergyInFile",
"help": "If kinetic and potential energies should be dumped in a CSV file at each iteration"
},
{
"name": "d_explicit",
"type": "bool",
"xmlname": "explicitIntegration",
"help": "Use explicit integration scheme"
},
{
"name": "d_computeHamiltonian",
"type": "bool",
"xmlname": "computeHamiltonian",
"help": "Compute hamiltonian"
},
{
"name": "d_hamiltonianEnergy",
"type": "SReal",
"xmlname": "hamiltonianEnergy",
"help": "hamiltonian energy"
},
{
"name": "d_useIncrementalPotentialEnergy",
"type": "bool",
"xmlname": "useIncrementalPotentialEnergy",
"help": "use real potential energy, if false use approximate potential energy"
},
{
"name": "d_threadSafeVisitor",
"type": "bool",
"xmlname": "threadSafeVisitor",
"help": "If true, do not use realloc and free visitors in fwdInteractionForceField."
}
],
"links": [],
"methods": [
{
"name": "init",
"return_type": "void",
"params": [],
"is_virtual": true,
"is_pure_virtual": false,
"is_static": false,
"access": "public"
},
{
"name": "solve",
"return_type": "void",
"params": [
{
"name": "params",
"type": "const core::ExecParams *"
},
{
"name": "dt",
"type": "SReal"
},
{
"name": "xResult",
"type": "sofa::core::MultiVecCoordId"
},
{
"name": "vResult",
"type": "sofa::core::MultiVecDerivId"
}
],
"is_virtual": true,
"is_pure_virtual": false,
"is_static": false,
"access": "public"
},
{
"name": "getVelocityIntegrationFactor",
"return_type": "SReal",
"params": [],
"is_virtual": true,
"is_pure_virtual": false,
"is_static": false,
"access": "public"
},
{
"name": "getPositionIntegrationFactor",
"return_type": "SReal",
"params": [],
"is_virtual": true,
"is_pure_virtual": false,
"is_static": false,
"access": "public"
},
{
"name": "getIntegrationFactor",
"return_type": "SReal",
"params": [
{
"name": "",
"type": "int"
},
{
"name": "",
"type": "int"
}
],
"is_virtual": true,
"is_pure_virtual": false,
"is_static": false,
"access": "public"
},
{
"name": "getSolutionIntegrationFactor",
"return_type": "SReal",
"params": [
{
"name": "",
"type": "int"
}
],
"is_virtual": true,
"is_pure_virtual": false,
"is_static": false,
"access": "public"
},
{
"name": "parse",
"return_type": "void",
"params": [
{
"name": "arg",
"type": "core::objectmodel::BaseObjectDescription *"
}
],
"is_virtual": true,
"is_pure_virtual": false,
"is_static": false,
"access": "public"
}
],
"description": "The `VariationalSymplecticSolver` is an ODE solver within the SOFA framework, specifically designed for implicit and explicit time integration using variational symplectic integrators as defined by Kharevych et al. (2006). This solver conserves linear momentum and mechanical energy. It is part of the `sofa::component::odesolver::backward` namespace and belongs to the Sofa.Component.ODESolver.Backward module. The component inherits from both `LinearSolverAccessor` and `OdeSolver`, integrating mechanical systems with a focus on accuracy through multiple Newton steps for implicit integration.\n\nIt interacts closely with other components such as `core::ExecParams` for execution parameters, `sofa::core::MultiVecCoordId` and `sofa::core::MultiVecDerivId` to handle positions and velocities, and leverages linear solvers via the `LinearSolverAccessor` interface. The solver supports error tolerance settings (`newtonError`) and a maximum number of Newton steps (`steps`). It also allows for Rayleigh damping through `rayleighStiffness` and `rayleighMass`, and can compute Hamiltonian energy with options to save it in a CSV file.\n\nThe component provides methods like `init()` for initialization, `solve()` for solving the ODE system over time steps, and other integration factor getters (`getVelocityIntegrationFactor()`, `getPositionIntegrationFactor()`) which are essential for its operation within the simulation loop. Additional fields allow toggling between explicit and implicit schemes, controlling energy computation and file output, and specifying thread-safe behavior.",
"maths": "The `VariationalSymplecticSolver` is an ODE solver within the SOFA framework designed for both implicit and explicit time integration using variational symplectic integrators as defined by Kharevych et al. (2006). This solver conserves linear momentum and mechanical energy, making it particularly suitable for simulations where such conservation laws are crucial.\n\n### Governing Equations\n\n#### Implicit Integration Scheme\nFor implicit time integration with \\(\\alpha = 0.5\\) (quadratic accuracy), the `VariationalSymplecticSolver` solves the following equation iteratively using Newton steps:\n\\[\nK(q_{k,i-1}) \res_i + f(q_{k,i-1}) + \\frac{4}{h^2}M\res_i - \\frac{2}{h}\bm{p}_k = 0\n\\]\nwhere:\n- \\(q_{k,i}\\) is the estimate of the mid-point position at iteration \\(k\\) and Newton step \\(i\\).\n- \\(\res_i\\) is the position increment that minimizes the variational energy.\n- \\(f(q_{k,i})\\) represents internal forces, which include potential energy gradients.\n- \\(K(q_{k,i})\\) is the tangent stiffness matrix.\n- \\(M\\) is the mass matrix.\n- \\(\bm{p}_k\\) is momentum at time step \\(k\\).\n- \\(h\\) is the time step size.\n\nThe solver computes \\(\res_i\\), updates the position estimate, and checks for convergence based on the specified error tolerance (`newtonError`).\n\n#### Explicit Integration Scheme\nFor explicit integration, the acceleration at each time step is computed as:\n\\[\n\\bm{a}_{k+1} = \\frac{\bm{f}(q_k) + r_M M \bm{v}_k}{M}\n\\]\nwhere \\(r_M\\) is the Rayleigh damping coefficient related to mass.\nThe updated velocity and position are then:\n\\[\n\bm{v}_{k+1} = \bm{v}_k + h \bm{a}_{k+1}\n\\]\nand\n\\[\nq_{k+1} = q_k + h \bm{v}_{k+1}\n\\]\nThis scheme does not require iterative solving and is less computationally expensive.\n\n### Constitutive Laws and Kinematic Models\nThe component uses Rayleigh damping to account for energy dissipation, characterized by the following damping forces:\n\\[\nf_d = r_M M \bm{v} + r_K K(q) q\n\\]\nwhere \\(r_K\\) is the Rayleigh stiffness coefficient.\nThe solver also handles constraints and projections using `mop.solveConstraint` to ensure that the updated positions, velocities, and momenta satisfy any imposed constraints.\n\n### Role in FEM Pipeline\nThe `VariationalSymplecticSolver` fits into the broader variational mechanics framework by conserving key physical quantities such as linear momentum and mechanical energy. It operates within the time integration phase of the simulation pipeline:\n1. **Initialization:** The solver initializes the system, including allocating memory for storing momentum.\n2. **Time Integration:** The solver computes forces and accelerations, updating positions and velocities using either explicit or implicit schemes depending on user settings (`d_explicit`).\n3. **Nonlinear Solve:** For implicit integration, it solves a nonlinear system iteratively using Newton steps to find the updated position that minimizes energy.\n4. **Energy Computation:** The solver computes kinetic and potential energies at each time step for monitoring purposes or saving in a CSV file if `d_saveEnergyInFile` is enabled.\n\n### Numerical Methods and Discretization Choices\nThe component uses variational symplectic integrators, which are designed to preserve the Hamiltonian structure of mechanical systems. This ensures that energy and momentum are conserved over long-time simulations. The choice of \\(\\alpha = 0.5\\) for implicit integration provides quadratic accuracy in time step size.\n\nThe solver supports Rayleigh damping, which is a form of proportional damping where the damping matrix is a linear combination of the mass and stiffness matrices.\nThe `VariationalSymplecticSolver` also leverages linear solvers (via `LinearSolverAccessor`) to solve the resulting systems of equations efficiently.\n\n### Summary\nThe `VariationalSymplecticSolver` is designed for accurate time integration with strong conservation properties. It integrates seamlessly into SOFA's modular framework, providing a robust and efficient solution for simulating mechanical systems.",
"abstract": "The `VariationalSymplecticSolver` conserves linear momentum and mechanical energy using variational symplectic integrators for both implicit and explicit time integration schemes.",
"sheet": "# VariationalSymplecticSolver\n\n## Overview\nThe `VariationalSymplecticSolver` is an ODE solver within the SOFA framework designed to conserve linear momentum and mechanical energy. It supports both implicit and explicit time integration using variational symplectic integrators, as defined by Kharevych et al. (2006). The component inherits from `LinearSolverAccessor` and `OdeSolver`, providing methods for initialization (`init()`), solving the ODE system (`solve()`), and computing integration factors.\n\n## Mathematical Model\n### Implicit Integration Scheme\nFor implicit time integration with \\(\\alpha = 0.5\\) (quadratic accuracy), the solver solves the following equation iteratively using Newton steps:\n\begin{equation}\nK(q_{k,i-1}) \\\\Delta s_i + f(q_{k,i-1}) + \\frac{4}{h^2}M\\\\Delta s_i - \\frac{2}{h}\\mathbf{p}_k = 0\n\\end{equation}\nwhere:\n- \\(q_{k,i}\\) is the estimate of the mid-point position at iteration \\(k\\) and Newton step \\(i\\).\n- \\\\Delta s_i is the position increment that minimizes the variational energy.\n- \\(f(q_{k,i})\\) represents internal forces, which include potential energy gradients.\n- \\(K(q_{k,i})\\) is the tangent stiffness matrix.\n- \\(M\\) is the mass matrix.\n- \\(\\mathbf{p}_k\\) is momentum at time step \\(k\\).\n- \\(h\\) is the time step size.\nThe solver computes \\\\Delta s_i, updates the position estimate, and checks for convergence based on the specified error tolerance (`newtonError`).\n\n### Explicit Integration Scheme\nFor explicit integration, the acceleration at each time step is computed as:\n\begin{equation}\n\\mathbf{a}_{k+1} = \\frac{\\mathbf{f}(q_k) + r_M M \\mathbf{v}_k}{M}\n\\end{equation}\nwhere \\(r_M\\) is the Rayleigh damping coefficient related to mass.\nThe updated velocity and position are then:\n\begin{equation}\n\\mathbf{v}_{k+1} = \\mathbf{v}_k + h \\mathbf{a}_{k+1}\n\\end{equation}\nand\n\begin{equation}\nq_{k+1} = q_k + h \\mathbf{v}_{k+1}\n\\end{equation}\nThis scheme does not require iterative solving and is less computationally expensive.\n\n## Parameters and Data\n- **newtonError (SReal)**: Error tolerance for Newton iterations.\n- **steps (unsigned int)**: Maximum number of Newton steps.\n- **rayleighStiffness (SReal)**: Rayleigh damping coefficient related to stiffness, > 0.\n- **rayleighMass (SReal)**: Rayleigh damping coefficient related to mass, > 0.\n- **saveEnergyInFile (bool)**: If kinetic and potential energies should be dumped in a CSV file at each iteration.\n- **explicitIntegration (bool)**: Use explicit integration scheme.\n- **computeHamiltonian (bool)**: Compute hamiltonian.\n- **hamiltonianEnergy (SReal)**: Hamiltonian energy.\n- **useIncrementalPotentialEnergy (bool)**: Use real potential energy, if false use approximate potential energy.\n- **threadSafeVisitor (bool)**: If true, do not use realloc and free visitors in fwdInteractionForceField."
}