siconos.kernel.Simulation (Python class)

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

Bases: object

Description of the simulation process (integrators, time discretisation and so on).

!!! This is an abstract class !!!

The available simulations are TimeStepping, EventDriven and TimeSteppingD1Minus.

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

Constructors

Simulation()

default constructor, for serialization

Simulation(NonSmoothDynamicalSystem nsds, TimeDiscretisation td)

default constructor

Parameters:
  • nsds – current nonsmooth dynamical system
  • td – the timeDiscretisation for this Simulation
Simulation(TimeDiscretisation td)

constructor with only a TimeDiscretisation

Parameters:td – the timeDiscretisation for this Simulation
advanceToEvent() → void[source]

advanceToEvent()=0 -> None

step from current event to next event of EventsManager

associate(OneStepIntegrator osi, DynamicalSystem ds) → None[source]

associate an OSI with a DS

clear() → None[source]

clear all maps.

This function should not exist, but there is a cycle with the shared_ptr: the OneStepIntegrator and OneStepNSProblem have both a link to the Simulation, and here we have all the OneStepIntegrator and OneStepNSProblem in maps. Then the memory is never freed. The clumsy way to deal with it is to call this function from the Model destructor to free the maps and then the cycle is broken

Warning: do not call this yourself, it is meant to be called from the
desctructor of the Model
clearNSDSChangeLog() → None[source]

clear the NSDS changelog up to current position.

If you have a particularly dynamic simulation (DS and Interactions created and destroyed frequently), then it is important to call this periodically.

computeOneStepNSProblem(int nb) → int[source]

computes a one step NS problem

Parameters:nb – the id of the OneStepNSProblem to be computed
Returns:information about the solver convergence.
computeResiduR() → bool[source]
computeResiduY() → bool[source]
currentTimeStep() → double[source]

returns current timestep

eventsManager() → EventsManager[source]

returns a pointer to the EventsManager

getPrintStat() → bool[source]

get printStat value

Returns:true if stats are activated
getTk() → double[source]

returns time instant k of the time discretisation

getTkp1() → double[source]

get time instant k+1 of the time discretisation

Warning: : this instant may be different from nextTime(), if for example
some non-smooth events or some sensor events are present
Returns:a double. If the simulation is near the end (t_{k+1} > T), it returns NaN.
getTkp2() → double[source]

get time instant k+2 of the time discretisation

Warning: : this instant may be different from nextTime(), if for example
some non-smooth events or some sensor events are present
Returns:a double. If the simulation is near the end (t_{k+2} > T), it returns NaN.
hasNextEvent() → bool[source]

true if a future event is to be treated or not (ie if some events remain in the eventsManager).

indexSet(int i) → InteractionsGraph[source]

get a pointer to indexSets[i]

Parameters:i – number of the required index set
Returns:a graph of interactions
initOSNS() → void[source]

initOSNS()=0 -> None

initialisation for OneStepNSProblem.

initialize() → None[source]

Complete initialisation of the Simulation (OneStepIntegrators, OneStepNSProblem, TImediscretisation).

initializeInteraction(double time, Interaction inter) → None[source]

Initialize a single Interaction for this Simulation, used for dynamic topology updates.

insertIntegrator(OneStepIntegrator osi) → None[source]

insert an Integrator into the simulation list of integrators

Parameters:osi – the OneStepIntegrator to add
insertInteractionManager(InteractionManager manager) → None[source]

Set an object to automatically manage interactions during the simulation.

Parameters:manager
insertNonSmoothProblem(OneStepNSProblem osns, int Id=SICONOS_OSNSP_DEFAULT) → None[source]

add a OneStepNSProblem in the Simulation

Parameters:
  • osns – the OneStepNSProblem to insert
  • Id – its id: default is SICONOS_OSNSP_DEFAULT, at impact level SICONOS_OSNSP_ED_IMPACT, at acceleration level SICONOS_OSNSP_ED_ACCELERATION
lambda_(level: unsigned int = 0, coor: unsigned int = 0) → SP::SiconosVector[source]

lambda( int level=0, int coor=0) -> array_like (np.float64, 1D)

return input lambda[level](coor) for all the interactions

Parameters:
  • level – lambda min order to be computed
  • coor – the coordinate of interest
Returns:

a SP::SiconosVector that contains the concatenated value

Add a new Interaction between one or a pair of DSs.

Parameters:
  • inter – the SP::Interaction to add
  • ds1 – the first SP::DynamicalSystem in the Interaction
  • ds2 – the second SP::DynamicalSystem in the Interaction, if any
name() → str[source]

get the name of the Simulation

Returns:std::string (the name of the Simulation)
nextTime() → double[source]

get “next time” (ie ending point for current integration, time of nextEvent of eventsManager.)

Returns:a double.
nonSmoothDynamicalSystem() → NonSmoothDynamicalSystem[source]

get the NonSmoothDynamicalSystem

Returns:NonSmoothDynamicalSystem
numberOfOSI() → int[source]

get the number of OSIs in the Simulation (ie the size of allOSI)

Returns:an unsigned int
numberOfOSNSProblems() → int[source]

get the number of OSNSP in the Simulation (ie the size of allNSProblems)

Returns:an unsigned int
oneStepIntegrators() → OSISet[source]

get all the Integrators of the Simulation

Returns:an OSISset
oneStepNSProblem(int id) → OneStepNSProblem[source]

get a OneStep nonsmooth problem of the simulation, identify with its number.

Parameters:id – number of the required osnspb
Returns:a pointer to OneStepNSProblem
oneStepNSProblems() → OneStepNSProblems[source]

get allNSProblems

Returns:a pointer to OneStepNSProblems object (container of SP::OneStepNSProblem)
processEvents() → None[source]

call eventsManager processEvents.

relativeConvergenceCriterionHeld() → bool[source]
Returns:true if the relative convergence criterion held.
relativeConvergenceTol() → double[source]
Returns:the relative convergence tolerence.
run() → None[source]

run the simulation, from t0 to T with default parameters if any particular settings has been done

setName(str newName) → None[source]

set the name of the Simulation

Parameters:newName – the new name
setNonSmoothDynamicalSystemPtr(NonSmoothDynamicalSystem newPtr) → None[source]

set the NonSmoothDynamicalSystem of the Simulation

Parameters:newPtr – a pointer on NonSmoothDynamicalSystem
setPrintStat(bool newVal) → None[source]

set printStat value: if true, print solver stats.

Parameters:newVal – true to activate stats
setRelativeConvergenceCriterionHeld(bool newVal) → None[source]
Parameters:newVal – a new relative convergence criterion
setRelativeConvergenceTol(double v) → None[source]

Set the relative convergence tolerance.

Parameters:v – tolerance value
setStaticLevels(bool b) → None[source]

set staticLevels

Parameters:b – decides whether levels should be computed at each iteration
setTolerance(double inputVal) → None[source]

set the value of offset for q dof vector in dynamical systems (to avoid events accumulation)

Parameters:inputVal – new tolerance
setUseRelativeConvergenceCriteron(bool use) → None[source]
startingTime() → double[source]

get “current time” (ie starting point for current integration, time of currentEvent of eventsManager.)

Returns:a double.
timeStep() → double[source]

get the current time step size (“next time”-“current time”)

Returns:a double.
tolerance() → double[source]

get tolerance

Returns:a double

Remove an Interaction from the simulation.

Parameters:inter – the SP::Interaction to remove
update(int level=0) → None[source]

update output, state, and input

Parameters:level – lambda order used to compute input level is set to 0 by default since in all time-stepping schemes we update all the state
updateIndexSet(level: unsigned int) → void[source]

updateIndexSet( int level)=0 -> None

update indexSets[i] of the topology, using current y and lambda values of Interactions.

Parameters:level – the number of the set to be updated
updateIndexSets() → None[source]

update all index sets of the topology, using current y and lambda values of Interactions.

updateInput(int level=0) → None[source]

update input

Parameters:level – lambda order used to compute input level is set to 0 by default since in all time-stepping schemes we update all the state
updateOutput(int level=0) → None[source]

update output

Parameters:level – lambda order used to compute output level is set to 0 by default since in all time-stepping schemes we update all the state
updateState(int level=0) → None[source]

update state of each dynamical system

updateT(double T) → None[source]

This updates the end of the Simulation.

Warning: this should be called only from the Model, to synchronise the 2
values
Parameters:T – the new final time
useRelativeConvergenceCriteron() → bool[source]
Returns:true if the relative convergence criterion is activated.
y(int level=0, int coor=0) -> array_like (np.float64, 1D)[source]

return output y[level](coor) for all the interactions

Parameters:
  • level – y min order to be computed
  • coor – the coordinate of interest
Returns:

a SP::SiconosVector that contains the concatenated value