Rod Interactions¶
Base class¶
This is the visible class that is output of the factory function.
-
class
BaseRodInteraction
: public mirheo::Interaction¶ Base class to manage rod interactions.
Rod interactions must be used with a RodVector. They are internal forces, meaning that halo() does not compute anything.
Subclassed by mirheo::RodInteraction< Nstates, StateParameters >
Public Functions
-
BaseRodInteraction
(const MirState *state, const std::string &name)¶ Construct a
BaseRodInteraction
.- Parameters
state
: The global state of the systemname
: Name of the interaction
-
void
halo
(ParticleVector *pv1, ParticleVector *pv2, CellList *cl1, CellList *cl2, cudaStream_t stream)¶ Compute interactions between bulk particles and halo particles.
The result of the interaction is
added to the corresponding channel of the ParticleVector. In general, the following interactions will be computed: pv1->halo() <> pv2->local() and pv2->halo() <> pv1->local().- Parameters
pv1
: first interacting ParticleVectorpv2
: second interacting ParticleVector. If it is the same as the pv1, self interactions will be computed.cl1
: cell-list built for the appropriate cut-off radius for pv1cl2
: cell-list built for the appropriate cut-off radius for pv2stream
: Execution stream
-
bool
isSelfObjectInteraction
() const¶ This is useful to know if we need exchange / cell-lists for that interaction. Example: membrane interactions are internal, all particles of a membrane are always on the same rank thus it does not need halo particles.
- Return
- boolean describing if the interaction is an internal interaction.
-
Implementation¶
The factory instantiates one of this templated class.
-
template <int Nstates, class StateParameters>
classRodInteraction
: public mirheo::BaseRodInteraction¶ Generic implementation of rod forces.
- Template Parameters
Nstates
: Number of polymorphic statesStateParameters
: parameters associated to the polymorphic state model
Public Functions
-
RodInteraction
(const MirState *state, const std::string &name, RodParameters parameters, StateParameters stateParameters, bool saveEnergies)¶ Construct a RodInteraction object.
- Parameters
state
: The global state of the systemname
: The name of the interactionparameters
: The common parameters from all kernel forcesstateParameters
: Parameters related to polymorphic states transitionsaveEnergies
:true
if the user wants to also compute the energies. In this case, energies will be saved in thechannel_names::energies
bisegment channel.
-
void
setPrerequisites
(ParticleVector *pv1, ParticleVector *pv2, CellList *cl1, CellList *cl2)¶ Add needed properties to the given ParticleVectors for future interactions.
Must be called before any other method of this class.
- Parameters
pv1
: One ParticleVector of the interactionpv2
: The other ParticleVector of that will interactcl1
: CellList of pv1cl2
: CellList of pv2
-
void
local
(ParticleVector *pv1, ParticleVector *pv2, CellList *cl1, CellList *cl2, cudaStream_t stream)¶ Compute interactions between bulk particles.
The result of the interaction is
added to the corresponding channel of the ParticleVector. The order of pv1 and pv2 may change the performance of the interactions.- Parameters
pv1
: first interacting ParticleVectorpv2
: second interacting ParticleVector. If it is the same as the pv1, self interactions will be computed.cl1
: cell-list built for the appropriate cut-off radius for pv1cl2
: cell-list built for the appropriate cut-off radius for pv2stream
: Execution stream
Kernels¶
The following support structure is used to compute the elastic energy:
-
template <int Nstates>
structBiSegment
¶ Helper class to compute elastic forces and energy on a bisegment.
- Template Parameters
Nstates
: Number of polymorphic states
Public Functions
-
__device__
BiSegment
(const RVview &view, int start)¶ Fetch bisegment data and prepare helper quantities.
-
__device__ rReal3
applyGrad0Bicur
(const rReal3 &v) const¶ compute gradient of the bicurvature w.r.t. r0 times v
-
__device__ rReal3
applyGrad2Bicur
(const rReal3 &v) const¶ compute gradient of the bicurvature w.r.t. r0 times v
-
__device__ void
computeBendingForces
(int state, const GPU_RodBiSegmentParameters<Nstates> ¶ms, rReal3 &fr0, rReal3 &fr2, rReal3 &fpm0, rReal3 &fpm1) const¶ compute the bending forces acting on the bisegment particles
This metho will add bending foces to the given variables. The other forces (on e.g. r1) can be computed by the symmetric nature of the model.
-
__device__ void
computeTwistForces
(int state, const GPU_RodBiSegmentParameters<Nstates> ¶ms, rReal3 &fr0, rReal3 &fr2, rReal3 &fpm0, rReal3 &fpm1) const¶ compute the torsion forces acting on the bisegment particles
This metho will add twist foces to the given variables. The other forces (on e.g. r1) can be computed by the symmetric nature of the model.
-
__device__ void
computeCurvatures
(rReal2 &kappa0, rReal2 &kappa1) const¶ Compute the curvatures along the material frames on each segment.
- Parameters
kappa0
: Curvature on the first segmentkappa1
: Curvature on the second segment
-
__device__ void
computeTorsion
(rReal &tau) const¶ Compute the torsion along the bisegment.
- Parameters
tau
: Torsion
-
__device__ void
computeCurvaturesGradients
(rReal3 &gradr0x, rReal3 &gradr0y, rReal3 &gradr2x, rReal3 &gradr2y, rReal3 &gradpm0x, rReal3 &gradpm0y, rReal3 &gradpm1x, rReal3 &gradpm1y) const¶ compute gradients of curvature term w.r.t. particle positions (see drivers)
-
__device__ void
computeTorsionGradients
(rReal3 &gradr0, rReal3 &gradr2, rReal3 &gradpm0, rReal3 &gradpm1) const¶ compute gradients of torsion term w.r.t. particle positions (see drivers)
-
__device__ rReal
computeEnergy
(int state, const GPU_RodBiSegmentParameters<Nstates> ¶ms) const¶ Compute the energy of the bisegment.
Public Members
-
rReal3
e0
¶ first segment
-
rReal3
e1
¶ second segment
-
rReal3
t0
¶ first segment direction
-
rReal3
t1
¶ second segment direction
-
rReal3
dp0
¶ first material frame direction
-
rReal3
dp1
¶ second material frame direction
-
rReal3
bicur
¶ bicurvature
-
rReal
bicurFactor
¶ helper scalar to compute bicurvature
-
rReal
e0inv
¶ 1 / length of first segment
-
rReal
e1inv
¶ 1 / length of second segment
-
rReal
linv
¶ 1 / l
-
rReal
l
¶ average of the lengths of the two segments