Integrators

See also the user interface.

Base class

class Integrator : public mirheo::MirSimulationObject

Advance ParticleVector objects in time.

Integrator objects are responsible to advance the state of ParticleVector objects on the device. After executed, the CellList of the ParticleVector object might be outdated; in that case, the Integrator is invalidates the current cell-lists, halo and redistribution status on the ParticleVector.

Subclassed by mirheo::IntegratorConstOmega, mirheo::IntegratorMinimize, mirheo::IntegratorOscillate, mirheo::IntegratorShear, mirheo::IntegratorShearPolChain, mirheo::IntegratorSubStep, mirheo::IntegratorSubStepShardlowSweep, mirheo::IntegratorTranslate, mirheo::IntegratorVV< ForcingTerm >, mirheo::IntegratorVVPolChain, mirheo::IntegratorVVRigid

Public Functions

Integrator(const MirState *state, const std::string &name)

Construct a Integrator object.

Parameters
  • state: The global state of the system. The time step and domain used during the execution are passed through this object.
  • name: The name of the integrator.

virtual void setPrerequisites(ParticleVector *pv)

Setup conditions on the ParticledVector.

Set specific properties to pv that will be modified during

execute(). Default: ask nothing. Must be called before execute() with the same pv.
Parameters

virtual void execute(ParticleVector *pv, cudaStream_t stream) = 0

Advance the ParticledVector for one time step.

Parameters
  • pv: The ParticleVector that will be advanced in time.
  • stream: The stream used for execution.

Derived classes

class IntegratorConstOmega : public mirheo::Integrator

Rotate ParticleVector objects with a constant angular velocity.

The center of rotation is defined in the global coordinate system.

Public Functions

IntegratorConstOmega(const MirState *state, const std::string &name, real3 center, real3 omega)

Construct a IntegratorConstOmega object.

Parameters
  • state: The global state of the system. The time step and domain used during the execution are passed through this object.
  • name: The name of the integrator.
  • center: The center of rotation, in global coordinates system.
  • omega: The angular velocity of rotation.

void execute(ParticleVector *pv, cudaStream_t stream)

Advance the ParticledVector for one time step.

Parameters
  • pv: The ParticleVector that will be advanced in time.
  • stream: The stream used for execution.

class IntegratorMinimize : public mirheo::Integrator

Energy minimization integrator.

Updates positions using a force-based gradient descent, without affecting or reading velocities.

Public Functions

IntegratorMinimize(const MirState *state, const std::string &name, real maxDisplacement)

Parameters
  • state: The global state of the system. The time step and domain used during the execution are passed through this object.
  • name: The name of the integrator.
  • maxDisplacement: Maximum particle displacement per time step.

void execute(ParticleVector *pv, cudaStream_t stream)

Advance the ParticledVector for one time step.

Parameters
  • pv: The ParticleVector that will be advanced in time.
  • stream: The stream used for execution.

class IntegratorOscillate : public mirheo::Integrator

Restrict ParticleVector velocities to a sine wave.

Set velocities to follow a sine wave:

\[v(t) = v \cos\left(\frac{ 2 \pi t} {T}\right)\]

The positions are integrated with forwards euler from the above velocities.

Public Functions

IntegratorOscillate(const MirState *state, const std::string &name, real3 vel, real period)

Parameters
  • state: The global state of the system. The time step and domain used during the execution are passed through this object.
  • name: The name of the integrator.
  • vel: Velocity magnitude.
  • period: The time taken for one oscillation.

void execute(ParticleVector *pv, cudaStream_t stream)

Advance the ParticledVector for one time step.

Parameters
  • pv: The ParticleVector that will be advanced in time.
  • stream: The stream used for execution.

class IntegratorVVRigid : public mirheo::Integrator

Integrate RigidObjectVector objects given torque and force.

Advance the RigidMotion and the frozen particles of the RigidObjectVector objects. The particles of each object are given the velocities corresponding to the rigid object motion.

Public Functions

IntegratorVVRigid(const MirState *state, const std::string &name)

Parameters
  • state: The global state of the system. The time step and domain used during the execution are passed through this object.
  • name: The name of the integrator.

void setPrerequisites(ParticleVector *pv)

Setup conditions on the ParticledVector.

Set specific properties to pv that will be modified during

execute(). Default: ask nothing. Must be called before execute() with the same pv.
Parameters

void execute(ParticleVector *pv, cudaStream_t stream)

Advance the ParticledVector for one time step.

Parameters
  • pv: The ParticleVector that will be advanced in time.
  • stream: The stream used for execution.

class IntegratorSubStepShardlowSweep : public mirheo::Integrator

Advance one MembraneVector associated with internal forces with smaller time step, similar to IntgratorSubStep.

We distinguish slow forces, which are computed outside of this class, from fast forces, computed only inside this class. Each time step given by the simulation is split into n sub time steps. Each of these sub time step advances the object using the non updated slow forces and the updated fast forces n times using the Shardlow method.

This was motivated by the separation of time scale of membrane viscosity (fast forces) and solvent viscosity (slow forces) in blood.

Warning

The fast forces should NOT be registered in the c Simulation. Otherwise, it will be executed twice (by the simulation and by this class).

Public Functions

IntegratorSubStepShardlowSweep(const MirState *state, const std::string &name, int substeps, BaseMembraneInteraction *fastForces, real gammaC, real kBT, int nsweeps)

construct a IntegratorSubStepShardlowSweep object.

Parameters
  • state: The global state of the system. The time step and domain used during the execution are passed through this object.
  • name: The name of the integrator.
  • substeps: Number of sub steps. Must be at least 1.
  • fastForces: Internal interactions executed at each sub step.
  • gammaC: The dissipation coefficient.
  • kBT: The temperature in energy units.
  • nsweeps: The number of iterations spent for each viscous update. Must be at least 1.

void execute(ParticleVector *pv, cudaStream_t stream)

Advance the ParticledVector for one time step.

Parameters
  • pv: The ParticleVector that will be advanced in time.
  • stream: The stream used for execution.

void setPrerequisites(ParticleVector *pv)

Setup conditions on the ParticledVector.

Set specific properties to pv that will be modified during

execute(). Default: ask nothing. Must be called before execute() with the same pv.
Parameters

class IntegratorSubStep : public mirheo::Integrator

Advance one ObjectVector associated with internal forces with smaller time step.

We distinguish slow forces, which are computed outside of this class, from fast forces, computed only inside this class. Each time step given by the simulation is split into n sub time steps. Each of these sub time step advances the object using the non updated slow forces and the updated fast forces n times.

This was motivated by the separation of time scale of membrane viscosity (fast forces) and solvent viscosity (slow forces) in blood.

Warning

The fast forces should NOT be registered in the c Simulation. Otherwise, it will be executed twice (by the simulation and by this class).

Public Functions

IntegratorSubStep(const MirState *state, const std::string &name, int substeps, const std::vector<Interaction *> &fastForces)

construct a IntegratorSubStep object.

This constructor will die if the fast forces need to exchange ghost particles with other ranks.

Parameters
  • state: The global state of the system. The time step and domain used during the execution are passed through this object.
  • name: The name of the integrator.
  • substeps: Number of sub steps
  • fastForces: Internal interactions executed at each sub step.

void execute(ParticleVector *pv, cudaStream_t stream)

Advance the ParticledVector for one time step.

Parameters
  • pv: The ParticleVector that will be advanced in time.
  • stream: The stream used for execution.

void setPrerequisites(ParticleVector *pv)

Setup conditions on the ParticledVector.

Set specific properties to pv that will be modified during

execute(). Default: ask nothing. Must be called before execute() with the same pv.
Parameters

class IntegratorShear : public mirheo::Integrator

Set ParticleVector velocities a linear shear.

The positions are integrated with forwards euler.

Public Functions

IntegratorShear(const MirState *state, const std::string &name, std::array<real, 9> shear, real3 origin)

Parameters
  • state: The global state of the system. The time step and domain used during the execution are passed through this object.
  • name: The name of the integrator.
  • origin: a point with zero flow (in global coordinates)
  • shear: shear rate tensor.

void execute(ParticleVector *pv, cudaStream_t stream)

Advance the ParticledVector for one time step.

Parameters
  • pv: The ParticleVector that will be advanced in time.
  • stream: The stream used for execution.

class IntegratorShearPolChain : public mirheo::Integrator

Set ParticleVector velocities a linear shear.

The positions are integrated with forwards Euler. Evolve the polymeric chain vectors with forwards Euler.

Public Functions

IntegratorShearPolChain(const MirState *state, const std::string &name, std::array<real, 9> shear, real3 origin)

Parameters
  • state: The global state of the system. The time step and domain used during the execution are passed through this object.
  • name: The name of the integrator.
  • origin: a point with zero flow (in global coordinates)
  • shear: shear rate tensor.

void execute(ParticleVector *pv, cudaStream_t stream)

Advance the ParticledVector for one time step.

Parameters
  • pv: The ParticleVector that will be advanced in time.
  • stream: The stream used for execution.

class IntegratorTranslate : public mirheo::Integrator

Restrict ParticleVector velocities to a constant.

The positions are integrated with forwards euler with a constant velocity.

Public Functions

IntegratorTranslate(const MirState *state, const std::string &name, real3 vel)

Parameters
  • state: The global state of the system. The time step and domain used during the execution are passed through this object.
  • name: The name of the integrator.
  • vel: Velocity magnitude.

void execute(ParticleVector *pv, cudaStream_t stream)

Advance the ParticledVector for one time step.

Parameters
  • pv: The ParticleVector that will be advanced in time.
  • stream: The stream used for execution.

template <class ForcingTerm>
class IntegratorVV : public mirheo::Integrator

Advance individual particles with Velocity-Verlet scheme.

Template Parameters
  • ForcingTerm: a functor that can add additional force to the particles depending on their position

Public Functions

IntegratorVV(const MirState *state, const std::string &name, ForcingTerm forcingTerm)

Parameters
  • state: The global state of the system. The time step and domain used during the execution are passed through this object.
  • name: The name of the integrator.
  • forcingTerm: Additional force added to the particles.

void execute(ParticleVector *pv, cudaStream_t stream)

Advance the ParticledVector for one time step.

Parameters
  • pv: The ParticleVector that will be advanced in time.
  • stream: The stream used for execution.

class IntegratorVVPolChain : public mirheo::Integrator

Advance individual particles with Velocity-Verlet scheme, and evolve their polymeric chain vector.

Public Functions

IntegratorVVPolChain(const MirState *state, const std::string &name)

Parameters
  • state: The global state of the system.
  • name: The name of the integrator.

void execute(ParticleVector *pv, cudaStream_t stream)

Advance the ParticledVector for one time step.

Parameters
  • pv: The ParticleVector that will be advanced in time.
  • stream: The stream used for execution.

Forcing terms

The forcing terms must follow the same interface. Currently implemented forcing terms:

class ForcingTermNone

No forcing term.

Public Functions

void setup(ParticleVector *pv, real t)

Initialize internal state.

This method must be called at every time step

Parameters
  • pv: the ParticleVector that will be updated
  • t: Current simulation time

real3 operator()(real3 original, Particle p) const

Add the additional force to the current one on a particle.

Return
The total force that must be applied to the particle
Parameters
  • original: Original force acting on the particle
  • p: Particle on which to apply the additional force

class ForcingTermConstDP

Apply a constant force independently of the position.

Public Functions

ForcingTermConstDP(real3 extraForce)

Construct a ForcingTermConstDP object.

Parameters
  • extraForce: The force to add to every particle

void setup(ParticleVector *pv, real t)

Initialize internal state.

This method must be called at every time step.

Parameters
  • pv: the ParticleVector that will be updated
  • t: Current simulation time

real3 operator()(real3 original, Particle p) const

Add the additional force to the current one on a particle.

Return
The total force that must be applied to the particle
Parameters
  • original: Original force acting on the particle
  • p: Particle on which to apply the additional force

class ForcingTermPeriodicPoiseuille

Apply equal but opposite forces in two halves of the global domain.

\[\begin{split}f_x = \begin{cases} F, & y_p > L_y / 2 \\ -F, & y_p \leqslant L_y / 2 \end{cases}\end{split}\]

Similarly, if the force is parallel to the y axis, its sign will depend on z; parallel to z it will depend on x.

Public Types

enum Direction

Encode directions.

Values:

x
y
z

Public Functions

ForcingTermPeriodicPoiseuille(real magnitude, Direction dir)

Construct a ForcingTermPeriodicPoiseuille object.

Parameters
  • magnitude: force magnitude to be applied.
  • dir: The force will be applied parallel to the specified axis.

void setup(ParticleVector *pv, real t)

Initialize internal state.

This method must be called at every time step.

Parameters
  • pv: the ParticleVector that will be updated
  • t: Current simulation time

real3 operator()(real3 original, Particle p) const

Add the additional force to the current one on a particle.

Return
The total force that must be applied to the particle
Parameters
  • original: Original force acting on the particle
  • p: Particle on which to apply the additional force