File kernel/src/modelingTools/NewtonEulerDS.hpp

Go to the source code of this file

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:

SP::SiconosVector forces() const

get forces

Return
pointer on a SiconosVector

SP::SimpleMatrix jacobianqForces() const

get JacobianqForces

Return
pointer on a SiconosMatrix

SP::SimpleMatrix jacobianvForces() const

get JacobianvForces

Return
pointer on a SiconosMatrix

Attributes access

virtual unsigned int dimension() const

return the number of degrees of freedom of the system

Return
an unsigned int.

virtual unsigned int getqDim() const

Returns dimension of vector q.

SP::SiconosVector q() const

get q (pointer link)

Return
pointer on a SiconosVector

SP::SiconosVector q0() const

get initial state (pointer link)

Return
pointer on a SiconosVector

SP::SiconosVector twist() const

get twist

Return
pointer on a SiconosVector

SP::SiconosVector velocity() const

get twist

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

SP::SiconosVector twist0() const
SP::SiconosVector linearVelocity(bool absoluteRef) const

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

Return
A SiconosVector of size 3 containing the linear velocity of this dynamical system.
Parameters
  • absoluteRef: If true, velocity is returned in the inertial frame, otherwise velocity is returned in the body frame.

void linearVelocity(bool absoluteRef, SiconosVector &v) const

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.

SP::SiconosVector angularVelocity(bool absoluteRef) const

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

Return
A SiconosVector of size 3 containing the angular velocity of this dynamical system.
Parameters
  • absoluteRef: If true, velocity is returned in the inertial frame, otherwise velocity is returned in the body frame.

void angularVelocity(bool absoluteRef, SiconosVector &w) const

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.

SP::SiconosVector p(unsigned int level = 2) const

get p

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

double scalarMass() const

get mass value

Return
a double

void setScalarMass(double mass)

Modify the scalar mass.

SP::SiconosMatrix inertia()
void setInertia(SP::SiconosMatrix newInertia)
void setInertia(double ix, double iy, double iz)
void updateMassMatrix()

to be called after scalar mass or inertia matrix have changed

SP::SiconosVector fExt() const

get fExt

Return
pointer on a plugged vector

void setFExtPtr(SP::SiconosVector newPtr)

set fExt to pointer newPtr

Parameters
  • newPtr: a SP to a Simple vector

void setMExtPtr(SP::SiconosVector newPtr)

set mExt to pointer newPtr

Parameters
  • newPtr: a SP to a Simple vector

SP::SiconosVector mGyr() const

get mGyr

Return
pointer on a plugged vector

SP::SimpleMatrix mass()
SP::SimpleMatrix inverseMass() const

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

Return
pointer SP::SimpleMatrix

SP::SimpleMatrix T()
SP::SimpleMatrix Tdot()
SP::SiconosVector dotq()
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

Memory vectors management

const SiconosMemory &qMemory()

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

Return
a memory

const SiconosMemory &twistMemory()

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

Return
a memory

void initMemory(unsigned int steps)

initialize the SiconosMemory objects with a positive size.

Parameters
  • steps: the size of the SiconosMemory (i)

void swapInMemory()

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

const SiconosMemory &forcesMemory()
const SiconosMemory &dotqMemory()

Miscellaneous public methods

double computeKineticEnergy()

To compute the kinetic energy.

void display() const

print the data to the screen

void setIsMextExpressedInInertialFrame(bool value)
void setNullifyMGyr(bool value)
virtual void normalizeq()
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.

Plugins management

void setComputeJacobianFIntqByFD(bool value)
void setComputeJacobianFIntvByFD(bool value)
void setComputeJacobianMIntqByFD(bool value)
void setComputeJacobianMIntvByFD(bool value)
void setComputeFExtFunction(const std::string &pluginPath, const std::string &functionName)

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

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

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

void setComputeFExtFunction(FExt_NE fct)

set a specified function to compute _fExt

Parameters
  • fct: a pointer on the plugin function

void setComputeMExtFunction(FExt_NE fct)

set a specified function to compute _mExt

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: the complete path to the plugin
  • functionName: the name of the function to use in this plugin

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

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

void setComputeFIntFunction(FInt_NE fct)

set a specified function to compute _fInt

Parameters
  • fct: a pointer on the plugin function

void setComputeMIntFunction(FInt_NE fct)

set a specified function to compute _mInt

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 setComputeJacobianFIntvFunction(const std::string &pluginPath, const std::string &functionName)

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

void setComputeJacobianFIntqFunction(FInt_NE fct)

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

Parameters
  • fct: a pointer on the plugin function

void setComputeJacobianFIntvFunction(FInt_NE fct)

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

Parameters
  • fct: a pointer on the plugin function

void setComputeJacobianMIntqFunction(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 setComputeJacobianMIntvFunction(const std::string &pluginPath, const std::string &functionName)

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

void setComputeJacobianMIntqFunction(FInt_NE fct)

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

Parameters
  • fct: a pointer on the plugin function

void setComputeJacobianMIntvFunction(FInt_NE fct)

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

Parameters
  • fct: a pointer on the plugin function

virtual void computeFExt(double time)

function to compute the external forces

Parameters
  • time: the current time

virtual void computeFExt(double time, SP::SiconosVector fExt)

default function to compute the external forces

Parameters
  • time: the current time
  • fExt: the computed external force (in-out param)

virtual void computeMExt(double time, SP::SiconosVector mExt)

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

virtual void computeMExt(double time)
void addExtForceAtPos(SP::SiconosVector force, bool forceAbsRef, SP::SiconosVector pos = SP::SiconosVector(), bool posAbsRef = false)

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.

void computeJacobianMExtqExpressedInInertialFrameByFD(double time, SP::SiconosVector q)
void computeJacobianMExtqExpressedInInertialFrame(double time, SP::SiconosVector q)
void computeFInt(double time, SP::SiconosVector q, SP::SiconosVector v)

default function to compute the internal forces

Parameters
  • time: the current time function to compute the internal forces
  • time: the current time
  • q:
  • v:

virtual void computeFInt(double time, SP::SiconosVector q, SP::SiconosVector v, SP::SiconosVector fInt)

default function to compute the internal forces

Parameters
  • time: the current time
  • q:
  • v:
  • fInt: the computed internal force vector

void computeMInt(double time, SP::SiconosVector q, SP::SiconosVector v)

default function to compute the internal moments

Parameters
  • time: the current time
  • q:
  • v:

virtual void computeMInt(double time, SP::SiconosVector q, SP::SiconosVector v, SP::SiconosVector mInt)

default function to compute the internal moments

Parameters
  • time: the current time
  • q:
  • v:
  • mInt: the computed internal moment vector

virtual void updatePlugins(double time)

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

Parameters
  • time: the current time

virtual void computeForces(double time)

Default function to compute forces.

Parameters
  • time: double, the current time

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

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

virtual void computeJacobianqForces(double time)

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

q of forces

Parameters
  • time: double, the current time

virtual void computeJacobianvForces(double time)

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

v of forces

Parameters
  • time: double, the current time

virtual void computeMGyr(SP::SiconosVector twist)

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

virtual void computeMGyr(SP::SiconosVector twist, SP::SiconosVector mGyr)

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

virtual void computeJacobianMGyrtwist(double time)

Default function to compute the jacobian following q of mGyr.

Parameters
  • time: the current time

virtual void computeJacobianMGyrtwistByFD(double time, SP::SiconosVector q, SP::SiconosVector twist)

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

void computeJacobianFIntq(double time)

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

void computeJacobianFIntv(double time)

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

Parameters
  • time: double : the current time

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

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

Parameters
  • time: double
  • position: SP::SiconosVector
  • twist: SP::SiconosVector

void computeJacobianFIntqByFD(double time, SP::SiconosVector position, SP::SiconosVector twist)

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

virtual void computeJacobianFIntv(double time, SP::SiconosVector position, SP::SiconosVector twist)

To compute the jacobian w.r.t.

v of the internal forces

Parameters
  • time: double: the current time
  • position: SP::SiconosVector
  • twist: SP::SiconosVector

void computeJacobianFIntvByFD(double time, SP::SiconosVector position, SP::SiconosVector twist)

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

virtual void computeJacobianMIntq(double time)

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

Parameters
  • time: double : the current time

virtual void computeJacobianMIntv(double time)

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

Parameters
  • time: double : the current time

virtual void computeJacobianMIntq(double time, SP::SiconosVector position, SP::SiconosVector twist)

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

Parameters
  • time: double : the current time,
  • position: SP::SiconosVector
  • twist: SP::SiconosVector

void computeJacobianMIntqByFD(double time, SP::SiconosVector position, SP::SiconosVector twist)

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

virtual void computeJacobianMIntv(double time, SP::SiconosVector position, SP::SiconosVector twist)

To compute the jacobian w.r.t.

v of the internal forces

Parameters
  • time: double: the current time
  • position: SP::SiconosVector
  • twist: SP::SiconosVector

void computeJacobianMIntvByFD(double time, SP::SiconosVector position, SP::SiconosVector twist)

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

virtual void computeT()
virtual void computeTdot()

Typedefs

typedef void (*FExt_NE)(double t, double *f, unsigned int size_z, double *z)
typedef void (*FInt_NE)(double t, double *q, double *v, double *f, unsigned int size_z, double *z)

Pointer to function for plug-in.

Enums

enum NewtonEulerDSRhsMatrices

Values:

jacobianXBloc00
jacobianXBloc01
jacobianXBloc10
jacobianXBloc11
zeroMatrix
zeroMatrixqDim
numberOfRhsMatrices

Functions

virtual ~NewtonEulerDS()

destructor

void _zeroPlugin()

build all _plugin…

PluggedObject

ACCEPT_STD_VISITORS()
double axisAngleFromQuaternion(double q0, double q1, double q2, double q3, SP::SiconosVector axis)
double axisAngleFromQuaternion(SP::SiconosVector q, SP::SiconosVector axis)
void changeFrameAbsToBody(const SiconosVector &q, SiconosVector &v)
void changeFrameAbsToBody(SP::SiconosVector q, SP::SiconosVector v)
void changeFrameAbsToBody(SP::SiconosVector q, SP::SimpleMatrix m)
void changeFrameBodyToAbs(const SiconosVector &q, SiconosVector &v)
void changeFrameBodyToAbs(SP::SiconosVector q, SP::SiconosVector v)
void changeFrameBodyToAbs(SP::SiconosVector q, SP::SimpleMatrix m)
void computeExtForceAtPos(SP::SiconosVector q, bool isMextExpressedInInertialFrame, SP::SiconosVector force, bool forceAbsRef, SP::SiconosVector pos, bool posAbsRef, SP::SiconosVector fExt, SP::SiconosVector mExt, bool accumulate)

Compute the force and moment vectors applied to a body with state q from a force vector at a given position.

void computeRotationMatrix(double q0, double q1, double q2, double q3, SP::SimpleMatrix rotationMatrix)
void computeRotationMatrix(SP::SiconosVector q, SP::SimpleMatrix rotationMatrix)
void computeRotationMatrixTransposed(SP::SiconosVector q, SP::SimpleMatrix rotationMatrix)
void computeT(SP::SiconosVector q, SP::SimpleMatrix T)
NewtonEulerDS()

Default constructor.

NewtonEulerDS(SP::SiconosVector position, SP::SiconosVector twist, double mass, SP::SiconosMatrix inertia)

constructor from a minimum set of data

Parameters

void normalizeq(SP::SiconosVector q)
void quaternionFromAxisAngle(SP::SiconosVector axis, double angle, SP::SiconosVector q)
void quaternionFromRotationVector(SP::SiconosVector rotationVector, SP::SiconosVector q)
void rotateAbsToBody(SP::SiconosVector q, SP::SiconosVector v)
void rotateAbsToBody(SP::SiconosVector q, SP::SimpleMatrix m)
void rotateAbsToBody(double q0, double q1, double q2, double q3, SiconosVector &v)
void rotateAbsToBody(double q0, double q1, double q2, double q3, SP::SiconosVector v)
void rotateAbsToBody(double q0, double q1, double q2, double q3, SP::SimpleMatrix m)
void rotationVectorFromQuaternion(double q0, double q1, double q2, double q3, SP::SiconosVector rotationVector)
void rotationVectorFromQuaternion(SP::SiconosVector q, SP::SiconosVector rotationVector)

Variables

SP::SiconosVector _acceleration

_acceleration contains the time derivative of the twist

SP::BoundaryCondition _boundaryConditions

Boundary condition applied to a dynamical system.

bool _computeJacobianFIntqByFD

If true, we compute the missing Jacobian by forward finite difference.

bool _computeJacobianFInttwistByFD

If true, we compute the missing Jacobian by forward finite difference.

bool _computeJacobianMIntqByFD

If true, we compute the missing Jacobian by forward finite difference.

bool _computeJacobianMInttwistByFD

If true, we compute the missing Jacobian by forward finite difference.

SP::SiconosVector _dotq

NewtonEuler non linear dynamical systems.

The equations of motion in the Newton-Euler formalism can be stated as

\left\{\begin{array}{rcl} M \dot v + F_{int}(q,v, \Omega, t)&=& F_{ext}(t), \ \ I \dot \Omega + \Omega \wedge I\Omega + M_{int}(q,v, \Omega, t) &=& M_{ext}(t), \ \ \dot q &=& T(q) [ v, \Omega] \ \ \dot R &=& R \tilde \Omega,\quad R^{-1}=R^T,\quad \det(R)=1 . \end{array}\right.

endrst

with

<ul> <li> f$x_G,v_Gf$ position and velocity of the center of mass expressed in a inertial frame of reference (world frame) </li> <li> f$Omegaf$ angular velocity vector expressed in the body-fixed frame (frame attached to the object) </li> <li> f$Rf$ rotation matrix form the inertial frame to the body-fixed frame f$R^{-1}=R^T, det(R)=1f$, i.e f$ Rin SO^+(3)f$ </li> <li> f$M=m,I_{3times 3}f$ diagonal mass matrix with f$m in mathbb{R}f$ the scalar mass </li> <li> f$If$ constant inertia matrix </li> <li> f$F_{ext}f$ and f$ M_{ext}f$ are the external applied forces and moment </li> </ul>

In the current implementation, f$Rf$ is parametrized by a unit quaternion.

*/ class NewtonEulerDS : public DynamicalSystem { protected:

/* serialization hooks */ ACCEPT_SERIALIZATION(NewtonEulerDS);

/** Common code for constructors
should be replaced in C++11 by delegating constructors

*/

void _init();

// &ndash; MEMBERS &ndash;

/** Dimension of _twist, _ndof=6 */ unsigned int _ndof;

/** _twist contains the twist of the Newton Euler dynamical system.
_twist[0:2] : f$v_G in RR^3 f$ velocity of the center of mass in
the inertial frame of reference (world frame).
_twist[3:5] : f$OmegainRR^3f$ angular velocity expressed in the body-fixed frame

*/

SP::SiconosVector _twist;

/** Initial twist */ SP::SiconosVector _twist0;

/** _q contains the representation of the system
In the current implementation, we have
_q[0:2] : the coordinates of the center of mass expressed
in the inertial frame of reference (world frame)
_q[3:6] : an unit quaternion representing the orientation of the solid.
This unit quaternion encodes the rotation mapping from the inertial frame of reference to the body-fixed frame

*/

SP::SiconosVector _q;

//SP::SiconosVector _deltaq;

/** Initial position */ SP::SiconosVector _q0;

/** Dimension of _q, is not necessary equal to _ndof. In our case, _qDim = 7 and _ndof =6*/ unsigned int _qDim;

/** The time derivative of f$qf$, f$dot qf$

SiconosMemory _dotqMemory
double _epsilonFD

value of the step in finite difference

SP::SiconosVector _fExt

external forces of the system

SP::SiconosVector _fInt

internal forces of the system

SiconosMemory _forcesMemory
bool _hasConstantFExt

boolean if _fext is constant (set thanks to setFExtPtr for instance) false by default

bool _hasConstantMExt

boolean if _mext is constant (set thanks to setMExtPtr for instance) false by default

SP::SiconosMatrix _I

Inertial matrix.

SP::SimpleMatrix _inverseMass

inverse or factorization of the mass of the system

bool _isMextExpressedInInertialFrame

if true, we assume that mExt is given in inertial frame (default false)

SP::SimpleMatrix _jacobianFIntq

jacobian_q FInt w.r.t q

SP::SimpleMatrix _jacobianFInttwist

jacobian_twist FInt w.r.t the twist

SP::SimpleMatrix _jacobianMExtq

jacobian_q MExt w.r.t q

SP::SimpleMatrix _jacobianMGyrtwist

jacobian_twist of mGyr w.r.t the twist

SP::SimpleMatrix _jacobianMIntq

jacobian_q MInt w.r.t q

SP::SimpleMatrix _jacobianMInttwist

jacobian_twist MInt w.r.t the twist

SP::SimpleMatrix _jacobianWrenchq

jacobian_q forces

SP::SimpleMatrix _jacobianWrenchTwist

jacobian_{twist} forces

SP::SimpleMatrix _massMatrix

used for concatenate _I and _scalarMass.I_3

SP::SiconosVector _mExt

external moment expressed in the inertial frame

SP::SiconosVector _mGyr

gyroscpical moment

SP::SiconosVector _mInt

external moment expressed in the body-fixed frame

internal moment of the forces

bool _nullifyMGyr

if true, we set the gyroscopic forces equal to 0 (default false)

std::vector<SP::SiconosVector> _p

“Reaction” due to the non smooth law - The index corresponds to the dynamic levels.

SP::PluggedObject _pluginFExt

Plugin to compute strength of external forces.

SP::PluggedObject _pluginFInt

Plugin to compute strength of internal forces.

SP::PluggedObject _pluginJacqFInt

The following code is commented because the jacobian of _mInt and _fInt are not yet used by the numerical scheme.

Will be needed by a fully implicit scheme for instance.jacobian_q jacobian_{qDot} NewtonEulerDS plug-in to compute \(\nabla_qF_{Int}(\dot q, q, t)\), id = “jacobianFIntq”

Parameters
  • time: : current time
  • sizeOfq: : size of vector q
  • q: : pointer to the first element of q
  • twist: : pointer to the first element of twist
  • jacob: : pointer to the first element of the jacobian
  • size: of vector z
  • z: : a vector of user-defined parameters

SP::PluggedObject _pluginJacqMInt

NewtonEulerDS plug-in to compute \(\nabla_qM_{Int}(\dot q, q, t)\), id = “jacobianMInttwist”.

Parameters
  • time: : current time
  • sizeOfq: : size of vector q
  • q: : pointer to the first element of q
  • twist: : pointer to the first element of twist
  • jacob: : pointer to the first element of the jacobian
  • size: of vector z
  • z: : a vector of user-defined parameters

SP::PluggedObject _pluginJactwistFInt

NewtonEulerDS plug-in to compute \(\nabla_{\dot q}F_{Int}(\dot q, q, t)\), id = “jacobianFIntTwist”.

Parameters
  • time: : current time
  • sizeOfq: : size of vector q
  • q: : pointer to the first element of q
  • twist: : pointer to the first element of twist
  • jacob: : pointer to the first element of the jacobian
  • size: of vector z
  • z: : a vector of user-defined parameters

SP::PluggedObject _pluginJactwistMInt

NewtonEulerDS plug-in to compute \(\nabla_{\dot q}M_{Int}(\dot q, q, t)\), id = “jacobianMInttwist”.

Parameters
  • time: : current time
  • sizeOfq: : size of vector q
  • q: : pointer to the first element of q
  • twist: : pointer to the first element of twist
  • jacob: : pointer to the first element of the jacobian
  • size: of vector z
  • z: : a vector of user-defined parameters

SP::PluggedObject _pluginMExt

Plugin to compute the external moment expressed in the inertial frame.

SP::PluggedObject _pluginMInt

Plugin to compute moments of internal forces.

SiconosMemory _qMemory
SP::SiconosVector _reactionToBoundaryConditions

Reaction to an applied boundary condition.

VectorOfSMatrices _rhsMatrices

A container of matrices to save matrices that are involed in first order from of NewtonEulerDS system values (jacobianXBloc10, jacobianXBloc11, zeroMatrix, idMatrix) No get-set functions at the time.

Only used as a protected member.

double _scalarMass

Scalar mass of the system.

SP::SimpleMatrix _T

Matrix depending on the parametrization of the orientation \(v = T(q) \dot q\).

SP::SimpleMatrix _Tdot

Time derivative of T.

\(\dot v = \dot T(q) \dot q + T(q) \ddot q\)

SiconosMemory _twistMemory

Memory vectors that stores the values within the timestep.

SP::SiconosVector _wrench

wrench (q,twist,t)= [ fExt - fInt ; mExtBodyFrame - mGyr - mInt ]^T