Class LagrangianDS

Defined in Program listing for file kernel/src/modelingTools/LagrangianDS.hpp

class LagrangianDS : public DynamicalSystem

Lagrangian non linear dynamical systems - \(M(q,z) \dot v = F(v, q, t, z) + p \).

This class defines and computes a generic ndof-dimensional Lagrangian Non Linear Dynamical System of the form :

\[\begin{split}\begin{cases} M(q,z) \dot v + F_{gyr}(v, q, z) + F_{int}(v , q , t, z) = F_{ext}(t, z) + p \\ \dot q = v \end{cases}\end{split}\]

where

  • \(q \in R^{ndof} \) is the set of the generalized coordinates,
  • \( \dot q =v \in R^{ndof} \) the velocity, i. e. the time derivative of the generalized coordinates (Lagrangian systems).
  • \( \ddot q =\dot v \in R^{ndof} \) the acceleration, i. e. the second time derivative of the generalized coordinates.
  • \( p \in R^{ndof} \) the reaction forces due to the Non Smooth Interaction.
  • \( M(q) \in R^{ndof \times ndof} \) is the inertia term saved in the SiconosMatrix mass().
  • \( F_{gyr}(\dot q, q) \in R^{ndof}\) is the non linear inertia term saved in the SiconosVector fGyr().
  • \( F_{int}(\dot q , q , t) \in R^{ndof} \) are the internal forces saved in the SiconosVector fInt().
  • \( F_{ext}(t) \in R^{ndof} \) are the external forces saved in the SiconosVector fExt().
  • \( z \in R^{zSize}\) is a vector of arbitrary algebraic variables, some sort of discrete state.

The equation of motion is also shortly denoted as \( M(q,z) \dot v = F(v, q, t, z) + p\)

where \(F(v, q, t, z) \in R^{ndof} \) collects the total forces acting on the system, that is \( F(v, q, t, z) = F_{ext}(t, z) - F_{gyr}(v, q, z) + F_{int}(v, q , t, z) \) This vector is stored in the SiconosVector forces().

q[i] is the derivative number i of q. Thus: q[0]= \( q \), global coordinates, q[1]= \( \dot q\), velocity, q[2]= \( \ddot q \), acceleration.

The following operators (and their jacobians) can be plugged, in the usual way (see User Guide, ‘User-defined plugins’)

If required (e.g. for Event-Driven like simulation), reformulation as a first-order system (DynamicalSystem) is possible, with:

  • \( n= 2 ndof \)
  • \( x = \left[\begin{array}{c}q \\ \dot q\end{array}\right]\)
  • rhs given by:

\[\begin{split}\dot x = \left[\begin{array}{c} \dot q\\ \ddot q = M^{-1}(q)\left[F(v, q , t, z) + p \right]\\ \end{array}\right]\end{split}\]

  • jacobian of the rhs, with respect to x

\[\begin{split}\nabla_{x}rhs(x,t) = \left[\begin{array}{cc} 0 & I \\ \nabla_{q}(M^{-1}(q)F(v, q , t, z)) & \nabla_{\dot q}(M^{-1}(q)F(v, q , t, z)) \\ \end{array}\right]\end{split}\]

  • input due to the non smooth law:

\[\begin{split}\left[\begin{array}{c} 0 \\ p \end{array}\right]\end{split}\]

Subclassed by LagrangianLinearDiagonalDS, LagrangianLinearTIDS

Right-hand side computation

void resetToInitialState()

reset the state to the initial state

void initRhs(double time)

allocate (if needed) and compute rhs and its jacobian.

Parameters
  • time: of initialization

void initializeNonSmoothInput(unsigned int level)

set nonsmooth input to zero

Parameters
  • level: input-level to be initialized.

virtual void computeRhs(double time)

update right-hand side for the current state

Parameters
  • time: of interest

virtual void computeJacobianRhsx(double time)

update \(\nabla_x rhs\) for the current state

Parameters
  • time: of interest

void resetAllNonSmoothParts()

reset non-smooth part of the rhs (i.e.

p), for all ‘levels’

void resetNonSmoothPart(unsigned int level)

set nonsmooth part of the rhs (i.e.

p) to zero for a given level

Parameters
  • level:

void setRhs(const SiconosVector &newValue)

set the value of the right-hand side, \( \dot x \)

Parameters

void setRhsPtr(SP::SiconosVector newPtr)

set right-hand side, \( \dot x \) (pointer link)

Parameters
  • newPtr: SP::SiconosVector

virtual void computeForces(double time, SP::SiconosVector q, SP::SiconosVector velocity)

function to compute \(F(v,q,t,z)\) for the current state

Parameters
  • time: the current timeCompute \(F(v,q,t,z)\)
  • time: the current time
  • q: SP::SiconosVector: pointers on q
  • velocity: SP::SiconosVector: pointers on velocity

virtual void computeJacobianqForces(double time)

Compute \(\nabla_qF(v,q,t,z)\) for current \(q,v\) Default function to compute forces.

Parameters
  • time: the current time

virtual void computeJacobianqDotForces(double time)

Compute \(\nabla_{\dot q}F(v,q,t,z)\) for current \(q,v\).

Parameters
  • time: the current time

Attributes access

virtual unsigned int dimension() const

return the number of degrees of freedom of the system

Return
an unsigned int.

SP::SiconosVector q() const

generalized coordinates of the system (vector of size dimension())

Return
pointer on a SiconosVector

void setQ(const SiconosVector &newValue)

set value of generalized coordinates vector (copy)

Parameters
  • newValue:

void setQPtr(SP::SiconosVector newPtr)

set value of generalized coordinates vector (pointer link)

Parameters
  • newPtr:

SP::SiconosVector q0() const

return initial state of the system

Return
pointer on a SiconosVector

void setQ0(const SiconosVector &newValue)

set initial state (copy)

Parameters
  • newValue:

void setQ0Ptr(SP::SiconosVector newPtr)

set initial state (pointer link)

Parameters
  • newPtr:

SP::SiconosVector velocity() const

get velocity vector (pointer link)

Return
pointer on a SiconosVector

void setVelocity(const SiconosVector &newValue)

set velocity vector (copy)

Parameters
  • newValue:

void setVelocityPtr(SP::SiconosVector newPtr)

set velocity vector (pointer link)

Parameters
  • newPtr:

SP::SiconosVector velocity0() const

get initial velocity (pointer)

Return
pointer on a SiconosVector

void setVelocity0(const SiconosVector &newValue)

set initial velocity (copy)

Parameters
  • newValue:

void setVelocity0Ptr(SP::SiconosVector newPtr)

set initial velocity (pointer link)

Parameters
  • newPtr:

SP::SiconosVector acceleration() const

get acceleration (pointer link)

Return
pointer on a SiconosVector

SP::SiconosVector p(unsigned int level) const

get input due to nonsmooth interaction, \(p\)

Return
pointer on a SiconosVector
Parameters
  • level: unsigned int, required level for p

void setP(const SiconosVector &newValue, unsigned int level)

set nonsmooth input (copy)

Parameters
  • newValue:
  • level: required level for p

void setPPtr(SP::SiconosVector newPtr, unsigned int level)

set nonsmooth input (pointer link)

Parameters
  • newPtr:
  • level: required level for p

SP::SiconosMatrix mass() const

get mass matrix (pointer link)

Return
SP::SiconosMatrix

SP::SimpleMatrix inverseMass() const

get (pointer) LU-factorization of the mass, used for LU-forward-backward computation

Return
pointer SP::SimpleMatrix

void setMassPtr(SP::SiconosMatrix newPtr)

set mass to pointer newPtr

Parameters
  • newPtr: a plugged matrix SP

SP::SiconosVector fInt() const

get $F_{int}$ (pointer link)

Return
pointer on a plugged vector

void setFIntPtr(SP::SiconosVector newPtr)

set $F_{int}$ (pointer link)

Parameters
  • newPtr: a SP to plugged vector

SP::SiconosVector fExt() const

get \(F_{ext}\), (pointer link)

Return
pointer on a plugged vector

void setFExtPtr(SP::SiconosVector newPtr)

set \(F_{ext}\), (pointer link)

Parameters
  • newPtr: a SP to a Simple vector

SP::SiconosVector fGyr() const

get \(F_{gyr}\), (pointer link)

Return
pointer on a plugged vector

void setFGyrPtr(SP::SiconosVector newPtr)

set \(F_{gyr}\), (pointer link)

Parameters
  • newPtr: a SP to plugged vector

SP::SiconosMatrix jacobianFIntq() const

get \(\nabla_qF_{int}\), (pointer link)

Return
pointer on a SiconosMatrix

SP::SiconosMatrix jacobianFIntqDot() const

get \(\nabla_{\dot q}F_{int}\), (pointer link)

Return
pointer on a SiconosMatrix

void setJacobianFIntqPtr(SP::SiconosMatrix newPtr)

set \(\nabla_{q}F_{int}\), (pointer link)

Parameters
  • newPtr: a SP SiconosMatrix

void setJacobianFIntqDotPtr(SP::SiconosMatrix newPtr)

set \(\nabla_{\dot q}F_{int}\), (pointer link)

Parameters
  • newPtr: a SP SiconosMatrix

SP::SiconosMatrix jacobianFGyrq() const

get \(\nabla_{q}F_{gyr}\), (pointer link)

Return
pointer on a SiconosMatrix

SP::SiconosMatrix jacobianFGyrqDot() const

get \(\nabla_{\dot q}F_{gyr}\), (pointer link)

Return
pointer on a SiconosMatrix

void setJacobianFGyrqPtr(SP::SiconosMatrix newPtr)

get \(\nabla_{q}F_{gyr}\), (pointer link)

Parameters
  • newPtr: a SP SiconosMatrix

void setJacobianFGyrqDotPtr(SP::SiconosMatrix newPtr)

get \(\nabla_{\dot q}F_{gyr}\), (pointer link)

Parameters
  • newPtr: a SP SiconosMatrix

SP::SiconosVector forces() const

get \( F(v,q,t,z)\) (pointer link)

Return
pointer on a SiconosVector

virtual SP::SiconosMatrix jacobianqForces() const

get \( \nabla_qF(v,q,t,z)\) (pointer link)

Return
pointer on a SiconosMatrix

virtual SP::SiconosMatrix jacobianqDotForces() const

get \( \nabla_{\dot q}F(v,q,t,z)\) (pointer link)

Return
pointer on a SiconosMatrix

Memory vectors management

SiconosMemory &qMemory()

get all the values of the state vector q stored in memory.

note: not const due to SchatzmanPaoliOSI::initializeWorkVectorsForDS

Return
a memory

SiconosMemory &velocityMemory()

get all the values of the state vector velocity stored in memory.

note: not const due to SchatzmanPaoliOSI::initializeWorkVectorsForDS

Return
a memory

const SiconosMemory &pMemory(unsigned int level)

get all the values of the state vector p stored in memory

Return
a memory
Parameters
  • level:

const SiconosMemory &forcesMemory()

get forces in memory buff

Return
pointer on a SiconosMemory

void initMemory(unsigned int size)

initialize the SiconosMemory objects with a positive size.

Parameters
  • size: the size of the SiconosMemory. must be >= 0

void swapInMemory()

push the current values of x, q and r in the stored previous values xMemory, qMemory, rMemory,

Plugins management

void setComputeMassFunction(const std::string &pluginPath, const std::string &functionName)

allow to set a specified function to compute the mass

Parameters
  • pluginPath: std::string : the complete path to the plugin
  • functionName: std::string : the name of the function to use in this plugin

void setComputeMassFunction(FPtr7 fct)

set a specified function to compute Mass

Parameters
  • fct: a pointer on the plugin function

void setComputeFIntFunction(const std::string &pluginPath, const std::string &functionName)

allow to set a specified function to compute FInt

Parameters
  • pluginPath: std::string : the complete path to the plugin
  • functionName: std::string : the name of the function to use in this plugin

void setComputeFIntFunction(FPtr6 fct)

set a specified function to compute fInt

Parameters
  • fct: a pointer on the plugin function

void setComputeFExtFunction(const std::string &pluginPath, const std::string &functionName)

allow to set a specified function to compute Fext

Parameters
  • pluginPath: std::string : the complete path to the plugin
  • functionName: std::string : the name of the function to use in this plugin

void setComputeFExtFunction(VectorFunctionOfTime fct)

set a specified function to compute fExt

Parameters
  • fct: a pointer on the plugin function

void setComputeFGyrFunction(const std::string &pluginPath, const std::string &functionName)

allow to set a specified function to compute the inertia

Parameters
  • pluginPath: std::string : the complete path to the plugin
  • functionName: std::string : the name of the function to use in this plugin

void setComputeFGyrFunction(FPtr5 fct)

set a specified function to compute FGyr

Parameters
  • fct: a pointer on the plugin function

void setComputeJacobianFIntqFunction(const std::string &pluginPath, const std::string &functionName)

allow to set a specified function to compute the jacobian w.r.t q of the internal forces

Parameters
  • pluginPath: std::string : the complete path to the plugin
  • functionName: std::string : the name of the function to use in this plugin

void setComputeJacobianFIntqDotFunction(const std::string &pluginPath, const std::string &functionName)

allow to set a specified function to compute the jacobian following qDot of the internal forces w.r.t.

Parameters
  • pluginPath: std::string : the complete path to the plugin
  • functionName: std::string : the name of the function to use in this plugin

void setComputeJacobianFIntqFunction(FPtr6 fct)

set a specified function to compute jacobian following q of the FInt

Parameters
  • fct: a pointer on the plugin function

void setComputeJacobianFIntqDotFunction(FPtr6 fct)

set a specified function to compute jacobian following qDot of the FInt

Parameters
  • fct: a pointer on the plugin function

void setComputeJacobianFGyrqFunction(const std::string &pluginPath, const std::string &functionName)

allow to set a specified function to compute the jacobian w.r.t q of the the external forces

Parameters
  • pluginPath: std::string : the complete path to the plugin
  • functionName: std::string : the name of the function to use in this plugin

void setComputeJacobianFGyrqDotFunction(const std::string &pluginPath, const std::string &functionName)

allow to set a specified function to compute the jacobian w.r.t qDot of the the external strength

Parameters
  • pluginPath: std::string : the complete path to the plugin
  • functionName: std::string : the name of the function to use in this plugin

void setComputeJacobianFGyrqFunction(FPtr5 fct)

set a specified function to compute the jacobian following q of FGyr

Parameters
  • fct: a pointer on the plugin function

void setComputeJacobianFGyrqDotFunction(FPtr5 fct)

set a specified function to compute the jacobian following qDot of FGyr

Parameters
  • fct: a pointer on the plugin function

virtual void computeMass()

default function to compute the mass

virtual void computeMass(SP::SiconosVector position)

function to compute the mass

Parameters
  • position: value used to evaluate the mass matrix

virtual void computeFInt(double time)

default function to compute the internal strengths

Parameters
  • time: the current time

virtual void computeFInt(double time, SP::SiconosVector position, SP::SiconosVector velocity)

function to compute the internal strengths with some specific values for position and velocity (ie not those of the current state).

Parameters
  • time: the current time,
  • position: value used to evaluate the internal forces
  • velocity: value used to evaluate the internal forces

virtual void computeFExt(double time)

default function to compute the external strengths

Parameters
  • time: the current time

virtual void computeFGyr()

default function to compute the inertia

virtual void computeFGyr(SP::SiconosVector position, SP::SiconosVector velocity)

function to compute the inertia with some specific values for q and velocity (ie not those of the current state).

Parameters
  • position: value used to evaluate the inertia forces
  • velocity: value used to evaluate the inertia forces

virtual void computeJacobianFIntq(double time)

To compute the jacobian w.r.t q of the internal forces.

Parameters
  • time: the current time

virtual void computeJacobianFIntqDot(double time)

To compute the jacobian w.r.t qDot of the internal forces.

Parameters
  • time: the current time

virtual void computeJacobianFIntq(double time, SP::SiconosVector position, SP::SiconosVector velocity)

To compute the jacobian w.r.t q of the internal forces.

Parameters
  • time: the current time
  • position: value used to evaluate the jacobian
  • velocity: value used to evaluate the jacobian

virtual void computeJacobianFIntqDot(double time, SP::SiconosVector position, SP::SiconosVector velocity)

To compute the jacobian w.r.t.

qDot of the internal forces

Parameters
  • time: the current time
  • position: value used to evaluate the jacobian
  • velocity: value used to evaluate the jacobian

virtual void computeJacobianFGyrq()

function to compute the jacobian w.r.t.

q of the inertia forces

virtual void computeJacobianFGyrqDot()

function to compute the jacobian w.r.t.

qDot of the inertia forces

virtual void computeJacobianFGyrq(SP::SiconosVector position, SP::SiconosVector velocity)

function to compute the jacobian w.r.t.

q of the inertia forces

Parameters
  • position: value used to evaluate the jacobian
  • velocity: value used to evaluate the jacobian

virtual void computeJacobianFGyrqDot(SP::SiconosVector position, SP::SiconosVector velocity)

function to compute the jacobian w.r.t.

qDot of the inertia forces

Parameters
  • position: value used to evaluate the jacobian
  • velocity: value used to evaluate the jacobian

virtual void updatePlugins(double time)

default function to update the plugins functions using a new time:

Parameters
  • time: the current time

Miscellaneous public methods

double computeKineticEnergy()

To compute the kinetic energy.

void display() const

print the data of the dynamical system on the standard output

void computePostImpactVelocity()

Computes post-impact velocity, using pre-impact velocity and impulse (p) value.

Used in EventDriven (LsodarOSI->updateState)

void setBoundaryConditions(SP::BoundaryCondition newbd)

set Boundary Conditions

Parameters
  • newbd: BoundaryConditions

SP::BoundaryCondition boundaryConditions()

get Boundary Conditions

Return
SP::BoundaryCondition pointer on a BoundaryConditions

void setReactionToBoundaryConditions(SP::SiconosVector newrbd)

set Reaction to Boundary Conditions

Parameters
  • newrbd: BoundaryConditions pointer

SP::SiconosVector reactionToBoundaryConditions()

get Reaction to Boundary Conditions

Return
pointer on a BoundaryConditions

void init_generalized_coordinates(unsigned int level)

Allocate memory for q[level], level > 1 Useful for some integrators that need q[2] or other coordinates vectors.

Parameters
  • level: the required level

void init_inverse_mass()

Allocate memory for the lu factorization of the mass of the system.

Useful for some integrators with system inversion involving the mass

void update_inverse_mass()

Update the content of the lu factorization of the mass of the system, if required.

void init_forces()

Allocate memory for forces and its jacobian.

Public Functions

LagrangianDS(SP::SiconosVector position, SP::SiconosVector velocity)

constructor from initial state only, \(dv = p \)

Parameters

LagrangianDS(SP::SiconosVector position, SP::SiconosVector velocity, SP::SiconosMatrix mass)

constructor from initial state and mass, \(Mdv = p \)

Parameters

LagrangianDS(SP::SiconosVector position, SP::SiconosVector velocity, const std::string &plugin)

constructor from initial state and mass (plugin) \(Mdv = p \)

Parameters

virtual ~LagrangianDS()

destructor