siconos.kernel.MoreauJeanGOSI (Python class)

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

Bases: siconos.kernel.OneStepIntegrator

One Step time Integrator for First Order Dynamical Systems for mechanical Systems (LagrangianDS and NewtonEulerDS)

This integrator is the work horse of the event–capturing time stepping schemes for mechanical systems. It is mainly based on the pioneering works of M. Jean and J.J. Moreau for the time integration of mechanical systems with unilateral contact, impact and Coulomb’s friction with \(\theta\) scheme

For the linear Lagrangina system, the scheme reads as

\begin{cases} M (v_{k+1}-v_k) + h K q_{k+\theta} + h C v_{k+\theta} - h F_{k+\theta} = p_{k+1} = G P_{k+1},\label{eq:MoreauTS-motion}\\[1mm] q_{k+1} = q_{k} + h v_{k+\theta}, \quad \\[1mm] U_{k+1} = G^\top\, v_{k+1}, \\[1mm] \begin{array}{lcl} 0 \leq U^\alpha_{k+1} + e U^\alpha_{k} \perp P^\alpha_{k+1} \geq 0,& \quad&\alpha \in \mathcal I_1, \\[1mm] P^\alpha_{k+1} =0,&\quad& \alpha \in \mathcal I \setminus \mathcal I_1, \end{array} \end{cases}

with \(\theta \in [0,1]\). The index set \(\mathcal I_1\) is the discrete equivalent to the rule that allows us to apply the Signorini condition at the velocity level. In the numerical practice, we choose to define this set by begin{equation} label{eq:index-set1} mathcal I_1 = {alpha in mathcal I mid G^top (q_{k} + h v_{k}) + w leq 0text{ and } U_k leq 0 }. end{equation}.

For more details, we refer to

M. Jean and J.J. Moreau. Dynamics in the presence of unilateral contacts and dry friction: a numerical approach. In G. Del Pietro and F. Maceri, editors, Unilateral problems in structural analysis. II, pages 151–196. CISM 304, Spinger Verlag, 1987.

J.J. Moreau. Unilateral contact and dry friction in finite freedom dynamics. In J.J. Moreau and Panagiotopoulos P.D., editors, Nonsmooth Mechanics and Applications, number 302 in CISM, Courses and lectures, pages 1–82. CISM 302, Spinger Verlag, Wien- New York, 1988a.

J.J. Moreau. Numerical aspects of the sweeping process. Computer Methods in Applied Mechanics and Engineering, 177:329–349, 1999.

M. Jean. The non smooth contact dynamics method. Computer Methods in Applied Mechanics and Engineering, 177:235–257, 1999.

and for a review :

V. Acary and B. Brogliato. Numerical Methods for Nonsmooth Dynamical Systems: Applications in Mechanics and Electronics, volume 35 of Lecture Notes in Applied and Computational Mechanics. Springer Verlag, 2008.

MoreauJeanGOSI class is used to define some time-integrators methods for a list of dynamical systems. A MoreauJeanGOSI instance is defined by the value of theta and the list of concerned dynamical systems.

Each DynamicalSystem is associated to a SiconosMatrix, named “W”, the “iteration” matrix” W matrices are initialized and computed in initializeIterationMatrixW and computeW. Depending on the DS type, they may depend on time t and DS state x.

For mechanical systems, the implementation uses _p for storing the the input due to the nonsmooth law. This MoreauJeanGOSI scheme assumes that the relative degree is two.

For Lagrangian systems, the implementation uses _p[1] for storing the discrete impulse.

Main functions:

  • computeFreeState(): computes xfree (or vfree), dynamical systems state without taking non-smooth part into account
  • updateState(): computes x (q,v), the complete dynamical systems states. See User’s guide for details.

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

Constructors

MoreauJeanGOSI(double theta=0.5, double gamma=numeric_limits< double >::quiet_NaN())

constructor from theta value only

Parameters:
  • theta – value for all linked DS (default = 0.5).
  • gamma – value for all linked DS (default = NaN and gamma is not used).
NSLcontrib(Interaction inter, OneStepNSProblem osnsp) → None[source]

Compute the nonsmooth law contribution.

Parameters:
  • inter – the interaction (for y_k)
  • osnsp – the non-smooth integrator
addInteractionInIndexSet(Interaction inter, int i) → bool[source]

Apply the rule to one Interaction to known if is it should be included in the IndexSet of level i.

Parameters:
  • inter – the Interaction to test
  • i – level of the IndexSet
Returns:

Boolean

computeFreeState() → None[source]

Perform the integration of the dynamical systems linked to this integrator without taking into account the nonsmooth input (_r or _p)

computeInitialNewtonState() → None[source]

compute the initial state of the Newton loop.

computeResidu() → double[source]

return the maximum of all norms for the “MoreauJeanGOSI-discretized” residus of DS

Returns:a double
computeW(double time, DynamicalSystem ds, array_like (np.float64, 2D) W) → None[source]

compute W MoreauJeanGOSI matrix at time t

Parameters:
  • time – (double)
  • ds – a pointer to DynamicalSystem
  • W – the matrix to compute
computeWBoundaryConditions(DynamicalSystem ds) → None[source]

compute WBoundaryConditionsMap[ds] MoreauJeanGOSI matrix at time t

Parameters:ds – a pointer to DynamicalSystem
display() → None[source]

Displays the data of the MoreauJeanGOSI’s integrator.

initializeIterationMatrixW(double time, DynamicalSystem ds) → None[source]

initialize iteration matrix W MoreauJeanGOSI matrix at time t

Parameters:
  • time
  • ds – a pointer to DynamicalSystem
initializeIterationMatrixWBoundaryConditions(DynamicalSystem ds) → None[source]

initialize iteration matrix WBoundaryConditionsMap[ds] MoreauJeanGOSI

Parameters:ds – a pointer to DynamicalSystem
initializeWorkVectorsForDS(double t, DynamicalSystem ds) → None[source]

initialization of the work vectors and matrices (properties) related to one dynamical system on the graph and needed by the osi

Parameters:
  • t – time of initialization
  • ds – the dynamical system
initializeWorkVectorsForInteraction(Interaction inter, InteractionProperties interProp, DynamicalSystemsGraph DSG) → None[source]

initialization of the work vectors and matrices (properties) related to one interaction on the graph and needed by the osi

Parameters:
  • inter – the interaction
  • interProp – the properties on the graph
  • DSG – the dynamical systems graph
initialize_nonsmooth_problems() → None[source]

Initialization process of the nonsmooth problems linked to this OSI.

integrate(double tinit, double tend, double tout, int notUsed) → None[source]

integrate the system, between tinit and tend (->iout=true), with possible stop at tout (->iout=false)

Parameters:
  • tinit – the initial time
  • tend – the end time
  • tout – the real end time
  • notUsed – useless flag (for MoreauJeanGOSI, used in LsodarOSI)
numberOfIndexSets() → int[source]

get the number of index sets required for the simulation

Returns:unsigned int
prepareNewtonIteration(double time) → None[source]

method to prepare the fist Newton iteration

Parameters:time
removeInteractionFromIndexSet(Interaction inter, int i) → bool[source]

Apply the rule to one Interaction to known if is it should be removed in the IndexSet of level i.

Parameters:
  • inter – the Interaction to test
  • i – level of the IndexSet
Returns:

Boolean

updatePosition(DynamicalSystem ds) → None[source]

update the state of the dynamical systems

Parameters:ds – the dynamical to update
updateState(int level) → None[source]

update the state of the dynamical systems

Parameters:level – the level of interest for the dynamics: not used at the time