Class EulerMoreauOSI¶

Defined in Program listing for file kernel/src/simulationTools/EulerMoreauOSI.hpp

class EulerMoreauOSI : public OneStepIntegrator ¶

One Step time Integrator for First Order Dynamical Systems.

This integrator is the work horse of the event—capturing time stepping schemes for first order systems. It is mainly based on some extensions of the Backward Euler and \( \theta-\gamma \) schemes proposed in the pionnering work of J.J. Moreau for the sweeping process

J.J. Moreau. Evolution problem associated with a moving convex set in a Hilbert space. Journal of Differential Equations, 26, pp 347—374, 1977.

Variants are now used to integrate LCS, Relay systems, Higher order sweeping process see for instance

Consistency of a time-stepping method for a class of piecewise linear networks

M.K. Camlibel, W.P.M.H. Heemels, and J.M. Schumacher IEEE Transactions on Circuits and Systems I, 2002, 49(3):349—357

Numerical methods for nonsmooth dynamical systems: applications in mechanics and electronics

V Acary, B Brogliato Springer Verlag 2008

Convergence of time-stepping schemes for passive and extended linear complementarity systems L. Han, A. Tiwari, M.K. Camlibel, and J.-S. Pang SIAM Journal on Numerical Analysis 2009, 47(5):3768-3796

On preserving dissipativity properties of linear complementarity dynamical systems with the &theta-method

Greenhalgh Scott, Acary Vincent, Brogliato Bernard Numer. Math., , 2013.

Main time—integration schemes are based on the following \( \theta-\gamma \) scheme

\[\begin{split} \begin{cases} \label{eq:toto1} M x_{k+1} = M x_{k} +h\theta f(x_{k+1},t_{k+1})+h(1-\theta) f(x_k,t_k) + h \gamma r(t_{k+1}) + h(1-\gamma)r(t_k) \\[2mm] y_{k+1} = h(t_{k+1},x_{k+1},\lambda _{k+1}) \\[2mm] r_{k+1} = g(x_{k+1},\lambda_{k+1},t_{k+1})\\[2mm] \mbox{nslaw} ( y_{k+1} , \lambda_{k+1}) \end{cases} \end{split}\]

where \( \theta = [0,1] \) and \( \gamma \in [0,1] \). As in Acary & Brogliato 2008, we call the previous problem the `one–step nonsmooth problem’.

Another variant can also be used (FullThetaGamma scheme)

\[\begin{split} \begin{cases} M x_{k+1} = M x_{k} +h f(x_{k+\theta},t_{k+1}) + h r(t_{k+\gamma})\\[2mm] y_{k+\gamma} = h(t_{k+\gamma},x_{k+\gamma},\lambda _{k+\gamma})\\[2mm] r_{k+\gamma} = g(x_{k+\gamma},\lambda_{k+\gamma},t_{k+\gamma})\\[2mm] \mbox{nslaw} ( y_{k+\gamma} , \lambda_{k+\gamma}) \end{cases} \end{split}\]

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

Each DynamicalSystem is associated to a SiconosMatrix, named “W”, which is 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 first order systems, the implementation uses _r for storing the the input due to the nonsmooth law. This EulerMoreauOSI scheme assumes that the relative degree is zero or one and one level for _r is sufficient

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.

Public Functions

EulerMoreauOSI(double theta)¶

constructor from theta value only

Parameters: theta – value for all DS.

EulerMoreauOSI(double theta, double gamma)¶

constructor from theta value only

Parameters

theta – value for all linked DS.
gamma – value for all linked DS.

inline virtual ~EulerMoreauOSI()¶: destructor

const SimpleMatrix getW(SP::DynamicalSystem ds = SP::DynamicalSystem())¶

get the value of W corresponding to DynamicalSystem ds

Parameters: ds – a pointer to DynamicalSystem, optional, default = nullptr. get W[0] in that case
Returns: SimpleMatrix

SP::SimpleMatrix W(SP::DynamicalSystem ds)¶

get W corresponding to DynamicalSystem ds

Parameters: ds – a pointer to DynamicalSystem
Returns: pointer to a SiconosMatrix

const SimpleMatrix getWBoundaryConditions(SP::DynamicalSystem ds = SP::DynamicalSystem())¶

get the value of WBoundaryConditions corresponding to DynamicalSystem ds

Parameters: ds – a pointer to DynamicalSystem, optional, default = nullptr. get WBoundaryConditions[0] in that case
Returns: SimpleMatrix

SP::SiconosMatrix WBoundaryConditions(SP::DynamicalSystem ds)¶

get WBoundaryConditions corresponding to DynamicalSystem ds

Parameters: ds – a pointer to DynamicalSystem, optional, default = nullptr. get WBoundaryConditions[0] in that case
Returns: pointer to a SiconosMatrix

inline double theta()¶

get theta

Returns: a double

inline void setTheta(double newTheta)¶

set the value of theta

Parameters: newTheta – a double

inline double gamma()¶

get gamma

Returns: a double

inline void setGamma(double newGamma)¶

set the value of gamma

Parameters: newGamma – a double

inline bool useGamma()¶

get bool useGamma

Returns: a bool

inline void setUseGamma(bool b)¶

set the boolean to indicate that we use gamma

Parameters: b – true if gamma has to be used, false otherwise

inline bool useGammaForRelation()¶

get bool gammaForRelation for the relation

Returns: a

inline void setUseGammaForRelation(bool newUseGammaForRelation)¶

set the boolean to indicate that we use gamma for the relation

Parameters: newUseGammaForRelation – a bool

virtual void initializeWorkVectorsForDS(double t, SP::DynamicalSystem ds) override¶

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

virtual void initializeWorkVectorsForInteraction(Interaction &inter, InteractionProperties &interProp, DynamicalSystemsGraph &DSG) override¶

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

inline virtual unsigned int numberOfIndexSets() const override¶

get the number of index sets required for the simulation

Returns: unsigned int

void initializeIterationMatrixW(double time, SP::DynamicalSystem ds)¶

initialize iteration matrix W EulerMoreauOSI matrix at time t

Parameters

time – the time (double)
ds – a pointer to DynamicalSystem

void computeW(double time, DynamicalSystem &ds, DynamicalSystemsGraph::VDescriptor &dsv, SiconosMatrix &W)¶

compute W EulerMoreauOSI matrix at time t

Parameters

time – the current time
ds – the DynamicalSystem
dsv – a descriptor of the ds on the graph (redundant to avoid invocation)
W – the matrix to compute

void computeWBoundaryConditions(SP::DynamicalSystem ds)¶

compute WBoundaryConditionsMap[ds] EulerMoreauOSI matrix at time t

Parameters: ds – a pointer to DynamicalSystem

void initializeIterationMatrixWBoundaryConditions(SP::DynamicalSystem ds)¶

initialize iteration matrix WBoundaryConditionsMap[ds] EulerMoreauOSI

Parameters: ds – a pointer to DynamicalSystem

virtual double computeResidu() override¶

Computes the residuFree and residu of all the DynamicalSystems.

Returns: the maximum of the 2-norm over all the residu

virtual void computeFreeState() override¶: Perform the integration of the dynamical systems linked to this integrator without taking into account the nonsmooth input r.

virtual void updateOutput(double time) override¶: update the output of the Interaction attached to this Integrator

virtual void updateInput(double time) override¶: update the input of the Interaction attached to this Integrator

virtual void updateOutput(double time, unsigned int level) override¶

update the output of the Interaction attached to this Integrator

Parameters

time – current time
level – level of interest for the dynamics

virtual void updateInput(double time, unsigned int level) override¶

update the input of the Interaction attached to this Integrator

Parameters

time – current time
level – level of interest for the dynamics

virtual double computeResiduOutput(double time, SP::InteractionsGraph indexSet) override¶

compute the residu of the output of the relation (y) This computation depends on the type of OSI

Parameters

time – time of computation
indexSet – the index set of the interaction that are concerned

virtual double computeResiduInput(double time, SP::InteractionsGraph indexSet) override¶

compute the residu of the input of the relation (R or p) This computation depends on the type of OSI

Parameters

time – time of computation
indexSet – the index set of the interaction that are concerned

virtual void computeFreeOutput(InteractionsGraph::VDescriptor &vertex_inter, OneStepNSProblem *osnsp) override¶

integrates the Interaction linked to this integrator, without taking non-smooth effects into account

Parameters

vertex_inter – of the interaction graph
osnsp – pointer to OneStepNSProblem

virtual void prepareNewtonIteration(double time) override¶

computes all the W matrices

Parameters: time – current time

virtual void integrate(double &tinit, double &tend, double &tout, int &useless) override¶

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

Parameters

tinit – initial time
tend – end time
tout – real end time
useless – flag (for EulerMoreauOSI, used in LsodarOSI)

virtual void updateState(const unsigned int level) override¶

updates the state of the Dynamical Systems

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

virtual void display() override¶: Displays the data of the EulerMoreauOSI’s integrator.