siconos.kernel.NewtonEulerDS (Python class)

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

Bases: siconos.kernel.DynamicalSystem

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

NewtonEulerDS()

Default constructor.

NewtonEulerDS(array_like (np.float64, 1D) position, array_like (np.float64, 1D) twist, double mass, array_like (np.float64, 2D) inertia)

constructor from a minimum set of data

Parameters:
  • position – initial coordinates of this DynamicalSystem
  • twist – initial twist of this DynamicalSystem
  • mass – the mass
  • inertia – the inertia matrix
T() -> array_like (np.float64, 2D)[source]
Tdot() -> array_like (np.float64, 2D)[source]
addExtForceAtPos(array_like (np.float64, 1D) force, bool forceAbsRef, array_like (np.float64, 1D) pos=array_like (np.float64, 1D)(), bool posAbsRef=False) → None[source]

Adds a force/torque impulse to a body’s FExt and MExt vectors in either absolute (inertial) or relative (body) frame.

Modifies contents of _fExt and _mExt! Therefore these must have been set as constant vectors using setFExtPtr and setMExtPtr prior to calling this function. Adjustments to _mExt will take into account the value of _isMextExpressedInInertialFrame.

Parameters:
  • force – A force vector to be added.
  • forceAbsRef – If true, force is in inertial frame, otherwise it is in body frame.
  • pos – A position at which force should be applied. If NULL, the center of mass is assumed.
  • posAbsRef – If true, pos is in inertial frame, otherwise it is in body frame.
angularVelocity(*args)[source]

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

angularVelocity(bool absoluteRef) -> array_like (np.float64, 1D)[source]

Get the angular velocity in the absolute (inertial) or relative (body) frame of reference.

Parameters:absoluteRef – If true, velocity is returned in the inertial frame, otherwise velocity is returned in the body frame.
Returns:
  • A SiconosVector of size 3 containing the angular velocity of this dynamical
  • system.
angularVelocity(bool absoluteRef, array_like (np.float64, 1D) w) → None[source]

Fill a SiconosVector with the angular velocity in the absolute (inertial) or relative (body) frame of reference.

Parameters:
  • absoluteRef – If true, velocity is returned in the inertial frame, otherwise velocity is returned in the body frame.
  • w – A SiconosVector of size 3 to receive the angular velocity.
boundaryConditions() → BoundaryCondition[source]

get Boundary Conditions

Returns:SP::BoundaryCondition pointer on a BoundaryConditions
computeFExt(*args)[source]

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

computeFExt(double time) → None[source]

function to compute the external forces

Parameters:time – the current time
computeFExt(double time, array_like (np.float64, 1D) fExt) → None[source]

default function to compute the external forces

Parameters:
  • time – the current time
  • fExt – the computed external force (in-out param)
computeFInt(*args)[source]

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

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

default function to compute the internal forces

Parameters:
  • time – the current time function to compute the internal forces
  • time – the current time
  • q
  • v
computeFInt(double time, array_like (np.float64, 1D) q, array_like (np.float64, 1D) v, array_like (np.float64, 1D) fInt) → None[source]

default function to compute the internal forces

Parameters:
  • time – the current time
  • q
  • v
  • fInt – the computed internal force vector
computeForces(*args)[source]

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

computeForces(double time) → None[source]

Default function to compute forces.

Parameters:time – double, the current time
computeForces(double time, array_like (np.float64, 1D) q, array_like (np.float64, 1D) twist) → None[source]

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

Parameters:
  • time – double : the current time
  • q – SP::SiconosVector: pointers on q
  • twist – SP::SiconosVector: pointers on twist
computeJacobianFIntq(*args)[source]

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

computeJacobianFIntq(double time) → None[source]

Default function to compute the jacobian following v of mGyr.

Parameters:
  • time – the current timeTo compute the jacobian w.r.t q of the internal forces
  • time – double : the current time
computeJacobianFIntq(double time, array_like (np.float64, 1D) position, array_like (np.float64, 1D) twist) → None[source]

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

Parameters:
  • time – double
  • position – SP::SiconosVector
  • twist – SP::SiconosVector
computeJacobianFIntqByFD(double time, array_like (np.float64, 1D) position, array_like (np.float64, 1D) twist) → None[source]

To compute the jacobian w.r.t q of the internal forces by forward finite difference.

Parameters:
  • time – double
  • position – SP::SiconosVector
  • twist – SP::SiconosVector
computeJacobianFIntv(*args)[source]

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

computeJacobianFIntv(double time) → None[source]

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

Parameters:time – double : the current time
computeJacobianFIntv(double time, array_like (np.float64, 1D) position, array_like (np.float64, 1D) twist) → None[source]

To compute the jacobian w.r.t.

v of the internal forces

Parameters:
  • time – double: the current time
  • position – SP::SiconosVector
  • twist – SP::SiconosVector
computeJacobianFIntvByFD(double time, array_like (np.float64, 1D) position, array_like (np.float64, 1D) twist) → None[source]

To compute the jacobian w.r.t v of the internal forces by forward finite difference.

Parameters:
  • time – double
  • position – SP::SiconosVector
  • twist – SP::SiconosVector
computeJacobianMExtqExpressedInInertialFrame(double time, array_like (np.float64, 1D) q) → None[source]
computeJacobianMExtqExpressedInInertialFrameByFD(double time, array_like (np.float64, 1D) q) → None[source]
computeJacobianMGyrtwist(double time) → None[source]

Default function to compute the jacobian following q of mGyr.

Parameters:time – the current time
computeJacobianMGyrtwistByFD(double time, array_like (np.float64, 1D) q, array_like (np.float64, 1D) twist) → None[source]

Default function to compute the jacobian following q of mGyr by forward finite difference.

Parameters:
  • time – the current time
  • q – current state
  • twist – pointer to twist vector
computeJacobianMIntq(*args)[source]

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

computeJacobianMIntq(double time) → None[source]

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

Parameters:time – double : the current time
computeJacobianMIntq(double time, array_like (np.float64, 1D) position, array_like (np.float64, 1D) twist) → None[source]

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

Parameters:
  • time – double : the current time,
  • position – SP::SiconosVector
  • twist – SP::SiconosVector
computeJacobianMIntqByFD(double time, array_like (np.float64, 1D) position, array_like (np.float64, 1D) twist) → None[source]

To compute the jacobian w.r.t q of the internal moments by forward finite difference.

Parameters:
  • time – double
  • position – SP::SiconosVector
  • twist – SP::SiconosVector
computeJacobianMIntv(*args)[source]

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

computeJacobianMIntv(double time) → None[source]

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

Parameters:time – double : the current time
computeJacobianMIntv(double time, array_like (np.float64, 1D) position, array_like (np.float64, 1D) twist) → None[source]

To compute the jacobian w.r.t.

v of the internal forces

Parameters:
  • time – double: the current time
  • position – SP::SiconosVector
  • twist – SP::SiconosVector
computeJacobianMIntvByFD(double time, array_like (np.float64, 1D) position, array_like (np.float64, 1D) twist) → None[source]

To compute the jacobian w.r.t v of the internal moments by forward finite difference.

Parameters:
  • time – double
  • position – SP::SiconosVector
  • twist – SP::SiconosVector
computeJacobianRhsx(time: double) → void[source]

computeJacobianRhsx(double time)=0 -> None

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

Parameters:time – of interest
computeJacobianqForces(double time) → None[source]

Default function to compute the jacobian w.r.t.

q of forces

Parameters:time – double, the current time
computeJacobianvForces(double time) → None[source]

Default function to compute the jacobian w.r.t.

v of forces

Parameters:time – double, the current time
computeKineticEnergy() → double[source]

To compute the kinetic energy.

computeMExt(*args)[source]

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

computeMExt(double time, array_like (np.float64, 1D) mExt) → None[source]

function to compute the external moments The external moments are expressed by default in the body frame, since the Euler equation for Omega is written in the body–fixed frame.

Nevertheless, if _isMextExpressedInInertialFrame) is set to true, we assume that the external moment is given in the inertial frame and we perform the rotation afterwards

Parameters:
  • time – the current time
  • mExt – the computed external moment (in-out param)
computeMExt(double time) → None[source]
computeMGyr(*args)[source]

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

computeMGyr(array_like (np.float64, 1D) twist) → None[source]

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

Parameters:twist – SP::SiconosVector: pointers on twist vector
computeMGyr(array_like (np.float64, 1D) twist, array_like (np.float64, 1D) mGyr) → None[source]

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

Parameters:
  • twist – pointer to twist vector
  • mGyr – pointer to gyroscopic forces
computeMInt(*args)[source]

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

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

default function to compute the internal moments

Parameters:
  • time – the current time
  • q
  • v
computeMInt(double time, array_like (np.float64, 1D) q, array_like (np.float64, 1D) v, array_like (np.float64, 1D) mInt) → None[source]

default function to compute the internal moments

Parameters:
  • time – the current time
  • q
  • v
  • mInt – the computed internal moment vector
computeRhs(time: double) → void[source]

computeRhs(double time)=0 -> None

update right-hand side for the current state

Parameters:time – of interest
computeT(*args)[source]

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

computeT() → None[source]
computeT(array_like (np.float64, 1D) q, array_like (np.float64, 2D) T) → None[source]
computeTdot() → None[source]
dimension() → int[source]

returns the dimension of the system (n for first order, ndof for Lagrangian).

Useful to avoid if(typeOfDS) when size is required.

display() → void[source]

display() =0 -> None

print the data of the dynamical system on the standard output

dotq() -> array_like (np.float64, 1D)[source]
dotqMemory() → SiconosMemory[source]
fExt() -> array_like (np.float64, 1D)[source]

get fExt

Returns:pointer on a plugged vector
forces() -> array_like (np.float64, 1D)[source]

get forces

Returns:pointer on a SiconosVector
forcesMemory() → SiconosMemory[source]
getqDim() → int[source]

Returns dimension of vector q.

inertia() -> array_like (np.float64, 2D)[source]
initMemory(int steps) → None[source]

initialize the SiconosMemory objects: reserve memory for i vectors in memory and reset all to zero.

Parameters:steps – the size of the SiconosMemory (i)
initRhs(time: double) → void[source]

initRhs(double time)=0 -> None

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

Parameters:time – of initialization
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(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) LU-factorization of the mass, used for LU-forward-backward computation

Returns:pointer SP::SimpleMatrix
jacobianqForces() -> array_like (np.float64, 2D)[source]

get JacobianqForces

Returns:pointer on a SiconosMatrix
jacobianvForces() -> array_like (np.float64, 2D)[source]

get JacobianvForces

Returns:pointer on a SiconosMatrix
linearVelocity(*args)[source]

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

linearVelocity(bool absoluteRef) -> array_like (np.float64, 1D)[source]

Get the linear velocity in the absolute (inertial) or relative (body) frame of reference.

Parameters:absoluteRef – If true, velocity is returned in the inertial frame, otherwise velocity is returned in the body frame.
Returns:
  • A SiconosVector of size 3 containing the linear velocity of this dynamical
  • system.
linearVelocity(bool absoluteRef, array_like (np.float64, 1D) v) → None[source]

Fill a SiconosVector with the linear velocity in the absolute (inertial) or relative (body) frame of reference.

Parameters:
  • absoluteRef – If true, velocity is returned in the inertial frame, otherwise velocity is returned in the body frame.
  • v – A SiconosVector of size 3 to receive the linear velocity.
mGyr() -> array_like (np.float64, 1D)[source]

get mGyr

Returns:pointer on a plugged vector
mass() -> array_like (np.float64, 2D)[source]
normalizeq(*args)[source]

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

normalizeq() → None[source]
normalizeq(array_like (np.float64, 1D) q) → None[source]
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() -> array_like (np.float64, 1D)[source]

get q (pointer link)

Returns:pointer on a SiconosVector
q0() -> array_like (np.float64, 1D)[source]

get initial state (pointer link)

Returns:pointer on a SiconosVector
qMemory() → SiconosMemory[source]

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

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 nonsmooth part of the rhs, for all ‘levels’

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

resetNonSmoothPart( int level)=0 -> None

set nonsmooth part of the rhs to zero for a given level

Parameters:level
resetToInitialState() → None[source]

reset the state x() to the initial state x0

scalarMass() → double[source]

get mass value

Returns:a double
setBoundaryConditions(BoundaryCondition newbd) → None[source]

set Boundary Conditions

Parameters:newbd – BoundaryConditions
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 – the complete path to the plugin
  • functionName – the name of the function to use in this plugin
setComputeFExtFunction(FExt_NE fct) → None[source]

set a specified function to compute _fExt

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 – the complete path to the plugin
  • functionName – the name of the function to use in this plugin
setComputeFIntFunction(FInt_NE fct) → None[source]

set a specified function to compute _fInt

Parameters:fct – a pointer on the plugin function
setComputeJacobianFIntqByFD(bool value) → None[source]
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(FInt_NE fct) → None[source]

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

Parameters:fct – a pointer on the plugin function
setComputeJacobianFIntvByFD(bool value) → None[source]
setComputeJacobianFIntvFunction(*args)[source]

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

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

allow to set a specified function to compute the jacobian following v 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
setComputeJacobianFIntvFunction(FInt_NE fct) → None[source]

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

Parameters:fct – a pointer on the plugin function
setComputeJacobianMIntqByFD(bool value) → None[source]
setComputeJacobianMIntqFunction(*args)[source]

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

setComputeJacobianMIntqFunction(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
setComputeJacobianMIntqFunction(FInt_NE fct) → None[source]

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

Parameters:fct – a pointer on the plugin function
setComputeJacobianMIntvByFD(bool value) → None[source]
setComputeJacobianMIntvFunction(*args)[source]

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

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

allow to set a specified function to compute the jacobian following v 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
setComputeJacobianMIntvFunction(FInt_NE fct) → None[source]

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

Parameters:fct – a pointer on the plugin function
setComputeMExtFunction(*args)[source]

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

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

allow to set a specified function to compute _mExt

Parameters:
  • pluginPath – the complete path to the plugin
  • functionName – the name of the function to use in this plugin
setComputeMExtFunction(FExt_NE fct) → None[source]

set a specified function to compute _mExt

Parameters:fct – a pointer on the plugin function
setComputeMIntFunction(*args)[source]

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

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

allow to set a specified function to compute _mInt

Parameters:
  • pluginPath – the complete path to the plugin
  • functionName – the name of the function to use in this plugin
setComputeMIntFunction(FInt_NE fct) → None[source]

set a specified function to compute _mInt

Parameters:fct – a pointer on the plugin function
setFExtPtr(array_like (np.float64, 1D) newPtr) → None[source]

set fExt to pointer newPtr

Parameters:newPtr – a SP to a Simple vector
setInertia(*args)[source]

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

setInertia(array_like (np.float64, 2D) newInertia) → None[source]
setInertia(double ix, double iy, double iz) → None[source]
setIsMextExpressedInInertialFrame(bool value) → None[source]
setMExtPtr(array_like (np.float64, 1D) newPtr) → None[source]

set mExt to pointer newPtr

Parameters:newPtr – a SP to a Simple vector
setNullifyMGyr(bool value) → None[source]
setReactionToBoundaryConditions(array_like (np.float64, 1D) newrbd) → None[source]

set Reaction to Boundary Conditions

Parameters:newrbd – BoundaryConditions pointer
setScalarMass(double mass) → None[source]

Modify the scalar mass.

swapInMemory() → void[source]

swapInMemory()=0 -> None

push the current values of x and r in memory (index 0 of memory is the last inserted vector) xMemory and rMemory,

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

get twist

Returns:pointer on a SiconosVector
twist0() -> array_like (np.float64, 1D)[source]
twistMemory() → SiconosMemory[source]

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

Returns:a memory
updateMassMatrix() → None[source]

to be called after scalar mass or inertia matrix have changed

updatePlugins(time: double) → void[source]

updatePlugins(double time)=0 -> None

call all plugged functions for the current state

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 twist

Returns:
  • pointer on a SiconosVector this accessor is left to get a uniform access to
  • velocity. This should be removed with MechanicalDS class