siconos.kernel.LagrangianDS (Python class)

class siconos.kernel.LagrangianDS(*args)[source]

Bases: siconos.kernel.SecondOrderDS

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 :

\[\]

M(q,z) dot v + F_{gyr}(v, q, z) + F_{int}(v , q , t, z) = F_{ext}(t, z) + p \ dot q = v

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 (access : mass() method).

  • \(F_{gyr}(\dot q, q) \in R^{ndof}\) is the non linear inertia term (access fGyr() method).

  • \(F_{int}(\dot q , q , t) \in R^{ndof}\) are the internal forces (access fInt() method).

  • \(F_{ext}(t) \in R^{ndof}\) are the external forces (access fExt() method).

  • \(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 saved and may be accessed using forces() method.

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’)

  • \(M(q)\) (computeMass())

  • \(F_{gyr}(v, q, z)\) (computeFGyr())

  • \(F_{int}(v , q , t, z)\) (computeFInt())

  • \(F_{ext}(t, z)\) (computeFExt())

If required (e.g. for Event-Driven like simulation), formulation as a first- order system is also available, and writes:

  • \(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}\]

with the input due to the non smooth law:

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

In that case, use the following methods:

  • initRhs() to allocate/initialize memory for these new operators,

  • rhs() to get the rhs vector

  • computeRhs(), computeJacobianRhsx() …, to update the content of rhs, its jacobians …

Generated class (swig), based on C++ header Program listing for file kernel/src/modelingTools/LagrangianDS.hpp.

Constructors

LagrangianDS(array_like (np.float64, 1D) position, array_like (np.float64, 1D) velocity)

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

Parameters
  • position – SiconosVector : initial coordinates of this DynamicalSystem

  • velocity – SiconosVector : initial velocity of this DynamicalSystem

LagrangianDS(array_like (np.float64, 1D) position, array_like (np.float64, 1D) velocity, array_like (np.float64, 2D) mass)

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

Parameters
  • position – SiconosVector : initial coordinates of this DynamicalSystem

  • velocity – SiconosVector : initial velocity of this DynamicalSystem

  • mass – SiconosMatrix : mass matrix

LagrangianDS(array_like (np.float64, 1D) position, array_like (np.float64, 1D) velocity, str plugin)

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

Parameters
  • position – SiconosVector : initial coordinates of this DynamicalSystem

  • velocity – SiconosVector : initial velocity of this DynamicalSystem

  • plugin – std::string: plugin path to compute mass matrix

acceleration() -> array_like (np.float64, 1D)[source]

get acceleration (pointer link)

Returns

pointer on a SiconosVector

allocateFExt() → None[source]
allocateFInt() → None[source]
allocateJacobianFIntq() → None[source]
allocateJacobianFIntqDot() → None[source]
allocateMass() → None[source]
computeFExt(double time) → None[source]

default function to compute the external strengths

Parameters

time – the current time

computeFGyr(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

computeFGyr() → None[source]

default function to compute the inertia

computeFGyr(array_like (np.float64, 1D) position, array_like (np.float64, 1D) velocity) → None[source]

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

computeFInt(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

computeFInt(double time) → None[source]

default function to compute the internal strengths

Parameters

time – the current time

computeFInt(double time, array_like (np.float64, 1D) position, array_like (np.float64, 1D) velocity) → None[source]

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

computeForces(double time, array_like (np.float64, 1D) q, array_like (np.float64, 1D) velocity) → None[source]

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

computeJacobianFGyrq(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

computeJacobianFGyrq() → None[source]

function to compute the jacobian w.r.t.

q of the inertia forces

computeJacobianFGyrq(array_like (np.float64, 1D) position, array_like (np.float64, 1D) velocity) → None[source]

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

computeJacobianFGyrqDot(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

computeJacobianFGyrqDot() → None[source]

function to compute the jacobian w.r.t.

qDot of the inertia forces

computeJacobianFGyrqDot(array_like (np.float64, 1D) position, array_like (np.float64, 1D) velocity) → None[source]

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

computeJacobianFIntq(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

computeJacobianFIntq(double time) → None[source]

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

Parameters

time – the current time

computeJacobianFIntq(double time, array_like (np.float64, 1D) position, array_like (np.float64, 1D) velocity) → None[source]

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

computeJacobianFIntqDot(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

computeJacobianFIntqDot(double time) → None[source]

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

Parameters

time – the current time

computeJacobianFIntqDot(double time, array_like (np.float64, 1D) position, array_like (np.float64, 1D) velocity) → None[source]

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

computeJacobianRhsx(double time) → None[source]

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

Parameters

time – of interest

computeJacobianqDotForces(double time) → None[source]

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

Parameters

time – the current time

computeJacobianqForces(double time) → None[source]

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

Parameters

time – the current time

computeJacobianvForces(double time) → None[source]

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

Parameters

time – the current time

computeKineticEnergy() → double[source]

To compute the kinetic energy.

computeMass(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

computeMass() → None[source]

default function to compute the mass

computeMass(array_like (np.float64, 1D) position) → None[source]

function to compute the mass

Parameters

position – value used to evaluate the mass matrix

computePostImpactVelocity() → None[source]

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

Used in EventDriven (LsodarOSI->updateState)

computeRhs(double time) → None[source]

update right-hand side for the current state

Parameters

time – of interest

display(bool brief=True) → None[source]

print the data of the dynamical system on the standard output

fExt() -> array_like (np.float64, 1D)[source]

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

Returns

pointer on a plugged vector

fGyr() -> array_like (np.float64, 1D)[source]

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

Returns

pointer on a plugged vector

fInt() -> array_like (np.float64, 1D)[source]

get $F_{int}$ (pointer link)

Returns

pointer on a plugged vector

forces() -> array_like (np.float64, 1D)[source]

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

Returns

pointer on a SiconosVector

forcesMemory() → SiconosMemory[source]

get forces in memory buff

Returns

pointer on a SiconosMemory

initMemory(int size) → None[source]

initialize the SiconosMemory objects with a positive size.

Parameters

size – the size of the SiconosMemory. must be >= 0

initRhs(double time) → None[source]

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

Parameters

time – of initialization

init_forces() → None[source]

Allocate memory for forces and its jacobian.

init_generalized_coordinates(int level) → None[source]

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

Parameters

level – the required level

init_inverse_mass() → None[source]

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

Useful for some integrators with system inversion involving the mass

initializeNonSmoothInput(int level) → None[source]

set nonsmooth input to zero

Parameters

level – input-level to be initialized.

jacobianFGyrq() -> array_like (np.float64, 2D)[source]

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

Returns

pointer on a SiconosMatrix

jacobianFGyrqDot() -> array_like (np.float64, 2D)[source]

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

Returns

pointer on a SiconosMatrix

jacobianFIntq() -> array_like (np.float64, 2D)[source]

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

Returns

pointer on a SiconosMatrix

jacobianFIntqDot() -> array_like (np.float64, 2D)[source]

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

Returns

pointer on a SiconosMatrix

jacobianqForces() -> array_like (np.float64, 2D)[source]

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

Returns

pointer on a SiconosMatrix

jacobianvForces() -> array_like (np.float64, 2D)[source]

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

Returns

pointer on a SiconosMatrix

pMemory(int level) → SiconosMemory[source]

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

Parameters

level

Returns

a memory

q() -> array_like (np.float64, 1D)[source]

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

Returns

pointer on a SiconosVector

qMemory() → SiconosMemory[source]

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

note: not const due to SchatzmanPaoliOSI::initializeWorkVectorsForDS

Returns

a memory

resetAllNonSmoothParts() → None[source]

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

p), for all ‘levels’

resetNonSmoothPart(int level) → None[source]

set nonsmooth part of the rhs (i.e.

  1. to zero for a given level

Parameters

level

resetToInitialState() → None[source]

reset the state to the initial state

setComputeFExtFunction(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

setComputeFExtFunction(str pluginPath, str functionName) → None[source]

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

setComputeFExtFunction(VectorFunctionOfTime fct) → None[source]

set a specified function to compute fExt

Parameters

fct – a pointer on the plugin function

setComputeFGyrFunction(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

setComputeFGyrFunction(str pluginPath, str functionName) → None[source]

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

setComputeFGyrFunction(FPtr5 fct) → None[source]

set a specified function to compute FGyr

Parameters

fct – a pointer on the plugin function

setComputeFIntFunction(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

setComputeFIntFunction(str pluginPath, str functionName) → None[source]

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

setComputeFIntFunction(FPtr6 fct) → None[source]

set a specified function to compute fInt

Parameters

fct – a pointer on the plugin function

setComputeJacobianFGyrqDotFunction(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

setComputeJacobianFGyrqDotFunction(str pluginPath, str functionName) → None[source]

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

setComputeJacobianFGyrqDotFunction(FPtr5 fct) → None[source]

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

Parameters

fct – a pointer on the plugin function

setComputeJacobianFGyrqFunction(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

setComputeJacobianFGyrqFunction(str pluginPath, str functionName) → None[source]

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

setComputeJacobianFGyrqFunction(FPtr5 fct) → None[source]

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

Parameters

fct – a pointer on the plugin function

setComputeJacobianFIntqDotFunction(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

setComputeJacobianFIntqDotFunction(str pluginPath, str functionName) → None[source]

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

setComputeJacobianFIntqDotFunction(FPtr6 fct) → None[source]

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

Parameters

fct – a pointer on the plugin function

setComputeJacobianFIntqFunction(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

setComputeJacobianFIntqFunction(str pluginPath, str functionName) → None[source]

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

setComputeJacobianFIntqFunction(FPtr6 fct) → None[source]

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

Parameters

fct – a pointer on the plugin function

setComputeMassFunction(*args)[source]

Warning - Overloaded function : multiple signatures available, check prototypes below.

setComputeMassFunction(str pluginPath, str functionName) → None[source]

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

setComputeMassFunction(FPtr7 fct) → None[source]

set a specified function to compute Mass

Parameters

fct – a pointer on the plugin function

setFExtPtr(array_like (np.float64, 1D) newPtr) → None[source]

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

Parameters

newPtr – a SP to a Simple vector

setFGyrPtr(array_like (np.float64, 1D) newPtr) → None[source]

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

Parameters

newPtr – a SP to plugged vector

setFIntPtr(array_like (np.float64, 1D) newPtr) → None[source]

set $F_{int}$ (pointer link)

Parameters

newPtr – a SP to plugged vector

setJacobianFGyrqDotPtr(array_like (np.float64, 2D) newPtr) → None[source]

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

Parameters

newPtr – a SP SiconosMatrix

setJacobianFGyrqPtr(array_like (np.float64, 2D) newPtr) → None[source]

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

Parameters

newPtr – a SP SiconosMatrix

setJacobianFIntqDotPtr(array_like (np.float64, 2D) newPtr) → None[source]

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

Parameters

newPtr – a SP SiconosMatrix

setJacobianFIntqPtr(array_like (np.float64, 2D) newPtr) → None[source]

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

Parameters

newPtr – a SP SiconosMatrix

setQ(array_like (np.float64, 1D) newValue) → None[source]

set value of generalized coordinates vector (copy)

Parameters

newValue

setQ0(array_like (np.float64, 1D) newValue) → None[source]

set initial state (copy)

Parameters

newValue

setQ0Ptr(array_like (np.float64, 1D) newPtr) → None[source]

set initial state (pointer link)

Parameters

newPtr

setQPtr(array_like (np.float64, 1D) newPtr) → None[source]

set value of generalized coordinates vector (pointer link)

Parameters

newPtr

setRhs(array_like (np.float64, 1D) newValue) → None[source]

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

Parameters

newValue – SiconosVector

setRhsPtr(array_like (np.float64, 1D) newPtr) → None[source]

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

Parameters

newPtr – SP::SiconosVector

setVelocity(array_like (np.float64, 1D) newValue) → None[source]

set velocity vector (copy)

Parameters

newValue

setVelocity0(array_like (np.float64, 1D) newValue) → None[source]

set initial velocity (copy)

Parameters

newValue

setVelocity0Ptr(array_like (np.float64, 1D) newPtr) → None[source]

set initial velocity (pointer link)

Parameters

newPtr

setVelocityPtr(array_like (np.float64, 1D) newPtr) → None[source]

set velocity vector (pointer link)

Parameters

newPtr

swapInMemory() → None[source]

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

updatePlugins(double time) → None[source]

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

Parameters

time – the current time

update_inverse_mass() → None[source]

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

velocity() -> array_like (np.float64, 1D)[source]

get velocity vector (pointer link)

Returns

pointer on a SiconosVector

velocity0() -> array_like (np.float64, 1D)[source]

get initial velocity (pointer)

Returns

pointer on a SiconosVector

velocityMemory() → SiconosMemory[source]

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

note: not const due to SchatzmanPaoliOSI::initializeWorkVectorsForDS

Returns

a memory