siconos.kernel.SecondOrderDS (Python class)

class siconos.kernel.SecondOrderDS(*args, **kwargs)[source]

Bases: siconos.kernel.DynamicalSystem

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

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

\[\]

M(q,z) dot v = F(v, q, t, z) + p \ dot q = G(q,v)

where

  • \(q \in R^{ndof}\) is the set of the coordinates,

  • \(\dot q =v \in R^{ndof}\) the velocity,

  • \(\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(\dot q , q , t) \in R^{ndof}\) are the forces (access forces() method).

  • \(z \in R^{zSize}\) is a vector of arbitrary algebraic variables, some sort of discrete state.

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(v , q , t, z)\) (computeF()) 3 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/SecondOrderDS.hpp.

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

get acceleration (pointer link)

Returns

pointer on a SiconosVector

boundaryConditions() → BoundaryCondition[source]

get Boundary Conditions

Returns

SP::BoundaryCondition pointer on a BoundaryConditions

computeForces(time: double, q: SP::SiconosVector, velocity: SP::SiconosVector) → void[source]

computeForces(double time, array_like (np.float64, 1D) q, array_like (np.float64, 1D) velocity)=0 -> None

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

computeJacobianRhsx(time: double) → void[source]

computeJacobianRhsx(double time)=0 -> None

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

Parameters

time – of interest

computeJacobianqForces(time: double) → void[source]

computeJacobianqForces(double time)=0 -> None

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

Parameters

time – the current time

computeJacobianvForces(time: double) → void[source]

computeJacobianvForces(double time)=0 -> None

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

Parameters

time – the current time

computeMass(*args)[source]

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

computeMass()=0 -> None

default function to compute the mass

computeMass(array_like (np.float64, 1D) position)=0 -> None

function to compute the mass

Parameters

position – value used to evaluate the mass matrix

computeRhs(time: double) → void[source]

computeRhs(double time)=0 -> None

update right-hand side for the current state

Parameters

time – of interest

dimension() → int[source]

return the number of degrees of freedom of the system

Returns

an unsigned int.

display(brief: bool = True) → void[source]

display(bool brief=True) =0 -> None

print the data of the dynamical system on the standard output

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

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

Returns

pointer on a SiconosVector

forcesMemory() → SiconosMemory const &[source]

forcesMemory()=0 -> SiconosMemory

get forces in memory buff

Returns

pointer on a SiconosMemory

initMemory(size: unsigned int) → void[source]

initMemory( int size)=0 -> None

initialize the SiconosMemory objects with a positive size.

Parameters

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

initRhs(time: double) → void[source]

initRhs(double time)=0 -> None

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

Parameters

time – of initialization

init_forces() → void[source]

init_forces()=0 -> None

Allocate memory for forces and its jacobian.

init_inverse_mass() → void[source]

init_inverse_mass()=0 -> None

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

Useful for some integrators with system inversion involving the mass

initializeNonSmoothInput(level: unsigned int) → void[source]

initializeNonSmoothInput( int level)=0 -> None

set nonsmooth input to zero

Parameters

level – input-level to be initialized.

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

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

Returns

pointer SP::SimpleMatrix

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

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

Returns

pointer on a SiconosMatrix

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

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

Returns

pointer on a SiconosMatrix

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

get mass matrix (pointer link)

Returns

SP::SiconosMatrix

p(int level=2) -> array_like (np.float64, 1D)[source]

get p

Parameters

level – unsigned int, required level for p, default = 2

Returns

pointer on a SiconosVector

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

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

Returns

pointer on a SiconosVector

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

get initial state (pointer link)

Returns

pointer on a SiconosVector

qMemory() → SiconosMemory const &[source]

qMemory()=0 -> SiconosMemory

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

note: not const due to SchatzmanPaoliOSI::initializeWorkVectorsForDS

Returns

a memory

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

get Reaction to Boundary Conditions

Returns

pointer on a BoundaryConditions

resetAllNonSmoothParts() → void[source]

resetAllNonSmoothParts()=0 -> None

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

p), for all ‘levels’

resetNonSmoothPart(level: unsigned int) → void[source]

resetNonSmoothPart( int level)=0 -> None

set nonsmooth part of the rhs (i.e.

  1. to zero for a given level

Parameters

level

resetToInitialState() → void[source]

resetToInitialState()=0 -> None

reset the state to the initial state

setBoundaryConditions(BoundaryCondition newbd) → None[source]

set Boundary Conditions

Parameters

newbd – BoundaryConditions

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

set mass to pointer newPtr

Parameters

newPtr – a plugged matrix SP

setQ(newValue: SiconosVector) → void[source]

setQ( array_like (np.float64, 1D) newValue)=0 -> None

set value of generalized coordinates vector (copy)

Parameters

newValue

setQ0(newValue: SiconosVector) → void[source]

setQ0( array_like (np.float64, 1D) newValue)=0 -> None

set initial state (copy)

Parameters

newValue

setQ0Ptr(newPtr: SP::SiconosVector) → void[source]

setQ0Ptr(array_like (np.float64, 1D) newPtr)=0 -> None

set initial state (pointer link)

Parameters

newPtr

setQPtr(newPtr: SP::SiconosVector) → void[source]

setQPtr(array_like (np.float64, 1D) newPtr)=0 -> None

set value of generalized coordinates vector (pointer link)

Parameters

newPtr

setReactionToBoundaryConditions(array_like (np.float64, 1D) newrbd) → None[source]

set Reaction to Boundary Conditions

Parameters

newrbd – BoundaryConditions pointer

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(newValue: SiconosVector) → void[source]

setVelocity( array_like (np.float64, 1D) newValue)=0 -> None

set velocity vector (copy)

Parameters

newValue

setVelocity0(newValue: SiconosVector) → void[source]

setVelocity0( array_like (np.float64, 1D) newValue)=0 -> None

set initial velocity (copy)

Parameters

newValue

setVelocity0Ptr(newPtr: SP::SiconosVector) → void[source]

setVelocity0Ptr(array_like (np.float64, 1D) newPtr)=0 -> None

set initial velocity (pointer link)

Parameters

newPtr

setVelocityPtr(newPtr: SP::SiconosVector) → void[source]

setVelocityPtr(array_like (np.float64, 1D) newPtr)=0 -> None

set velocity vector (pointer link)

Parameters

newPtr

swapInMemory() → void[source]

swapInMemory()=0 -> None

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

updatePlugins(time: double) → void[source]

updatePlugins(double time)=0 -> None

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

Parameters

time – the current time

update_inverse_mass() → void[source]

update_inverse_mass()=0 -> None

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

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

get velocity vector (pointer link)

Returns

pointer on a SiconosVector

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

get initial velocity (pointer)

Returns

pointer on a SiconosVector

velocityMemory() → SiconosMemory const &[source]

velocityMemory()=0 -> SiconosMemory

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

note: not const due to SchatzmanPaoliOSI::initializeWorkVectorsForDS

Returns

a memory