OpenMPCD
Namespaces | Functions
OpenMPCD::ImplementationDetails Namespace Reference

Namespace for implementation details that do not have to concern the users of the main functionality. More...

Namespaces

 Configuration
 Namespace for implementation details of OpenMPCD::Configuration.
 

Functions

unsigned int getCollisionCellIndex (const Vector3D< MPCParticlePositionType > &position, const unsigned int boxSizeX, const unsigned int boxSizeY, const unsigned int boxSizeZ)
 Returns the collision cell index for the given position. More...
 
void sortIntoCollisionCellsLeesEdwards (const unsigned int particleID, const MPCParticlePositionType *const gridShift_, const FP shearRate, const FP mpcTime, const MPCParticlePositionType *const positions, MPCParticleVelocityType *const velocities, MPCParticleVelocityType *const velocityCorrections, unsigned int *const collisionCellIndices, const unsigned int boxSizeX, const unsigned int boxSizeY, const unsigned int boxSizeZ)
 Sorts the given particle into the collision cells, temporarily applying Lees-Edwards boundary conditions. More...
 
void collisionCellContributions (const unsigned int particleCount, const MPCParticleVelocityType *const velocities, const unsigned int *const collisionCellIndices, MPCParticleVelocityType *const collisionCellMomenta, FP *const collisionCellMasses, const FP particleMass)
 Computes the collision cell mass and momentum contributions by the given particles. More...
 
Vector3D< MPCParticleVelocityTypegetCollisionCellCenterOfMassMomentum (const unsigned int collisionCellIndex, const MPCParticleVelocityType *const collisionCellMomenta)
 Returns the center-of-mass momentum of a collision cell. More...
 
Vector3D< MPCParticleVelocityTypegetCollisionCellCenterOfMassVelocity (const unsigned int collisionCellIndex, const MPCParticleVelocityType *const collisionCellMomenta, const FP *const collisionCellMasses)
 Returns the center-of-mass velocity of a collision cell. More...
 
void collisionCellStochasticRotationStep1 (const unsigned int particleCount, MPCParticleVelocityType *const velocities, const FP particleMass, const unsigned int *const collisionCellIndices, const MPCParticleVelocityType *const collisionCellMomenta, const FP *const collisionCellMasses, const MPCParticlePositionType *const collisionCellRotationAxes, const FP collisionAngle, FP *const collisionCellFrameInternalKineticEnergies, unsigned int *const collisionCellParticleCounts)
 Applies the first step of the SRD rotation to the given particles. More...
 
void collisionCellStochasticRotationStep2 (const unsigned int particleCount, MPCParticleVelocityType *const velocities, const unsigned int *const collisionCellIndices, const MPCParticleVelocityType *const collisionCellMomenta, const FP *const collisionCellMasses, const FP *const collisionCellVelocityScalings)
 Applies the first step of the SRD rotation to the given particles. More...
 
void undoLeesEdwardsVelocityCorrections (const unsigned int particleCount, MPCParticleVelocityType *const velocities, const MPCParticleVelocityType *const velocityCorrections)
 Undoes the velocity corrections applied by sortIntoCollisionCellsLeesEdwards. More...
 

Detailed Description

Namespace for implementation details that do not have to concern the users of the main functionality.

Function Documentation

◆ collisionCellContributions()

void OpenMPCD::ImplementationDetails::collisionCellContributions ( const unsigned int  particleCount,
const MPCParticleVelocityType *const  velocities,
const unsigned int *const  collisionCellIndices,
MPCParticleVelocityType *const  collisionCellMomenta,
FP *const  collisionCellMasses,
const FP  particleMass 
)

Computes the collision cell mass and momentum contributions by the given particles.

Parameters
[in]particleCountThe number of particles there are.
[in]velocitiesThe particle velocities.
[in]collisionCellIndicesThe collision cell indices for the particles.
[in,out]collisionCellMomentaThe collision cell momenta. This array is expected to be long enough to store 3 coordinates for each collision cell, and is assumed to have each element set to 0 prior to calling this function.
[in,out]collisionCellMassesThe masses of the collision cells. This array is expected to be long enough to store 1 element for each collision cell, and is assumed to have each element set to 0 prior to calling this function.
[in]particleMassThe mass of any one particle.

Definition at line 74 of file Simulation.cpp.

◆ collisionCellStochasticRotationStep1()

void OpenMPCD::ImplementationDetails::collisionCellStochasticRotationStep1 ( const unsigned int  particleCount,
MPCParticleVelocityType *const  velocities,
const FP  particleMass,
const unsigned int *const  collisionCellIndices,
const MPCParticleVelocityType *const  collisionCellMomenta,
const FP *const  collisionCellMasses,
const MPCParticlePositionType *const  collisionCellRotationAxes,
const FP  collisionAngle,
FP *const  collisionCellFrameInternalKineticEnergies,
unsigned int *const  collisionCellParticleCounts 
)

Applies the first step of the SRD rotation to the given particles.

Parameters
[in]particleCountThe number of particles there are.
[in,out]velocitiesThe particle velocities as input; after the function call, the values in this array will contain the velocities in the collision cell's center-of-mass frame.
[in]particleMassThe mass of any one particle.
[in]collisionCellIndicesAn array which holds, for each particle, the index of the collision cell it is currently in.
[in]collisionCellMomentaAn array holding, for each collision cell, its three cartesian momentum coordinates.
[in]collisionCellMassesAn array holding, for each collision cell, its total mass content.
[in]collisionCellRotationAxesAn array holding, for each collision cell, the random rotation axis it is to use during the collision. Each axis consists of the x, y, and z coordinate of a normalized vector.
[in]collisionAngleThe rotation angle to use for the SRD collision.
[out]collisionCellFrameInternalKineticEnergiesAn array, long enough to hold an element for each collision cell. The respective elements are set to the sum of the kinetic energies of the particles in that collision cell, as measured in that collision cell's center-of-mass frame. The array is not set to 0 by this function prior to execution of the algorithm.
[out]collisionCellParticleCountsAn array, long enough to hold an element for each collision cell. The respective elements are set to the number of particles in that collision cell. The array is assumed to contain only values 0 prior to calling this function.

Definition at line 131 of file Simulation.cpp.

◆ collisionCellStochasticRotationStep2()

void OpenMPCD::ImplementationDetails::collisionCellStochasticRotationStep2 ( const unsigned int  particleCount,
MPCParticleVelocityType *const  velocities,
const unsigned int *const  collisionCellIndices,
const MPCParticleVelocityType *const  collisionCellMomenta,
const FP *const  collisionCellMasses,
const FP *const  collisionCellVelocityScalings 
)

Applies the first step of the SRD rotation to the given particles.

Parameters
[in]particleCountThe number of particles there are.
[in,out]velocitiesThe particle new velocities as output; as input, the values in this array are to contain the rotated velocities in the collision cell's center-of-mass frame.
[in]collisionCellIndicesAn array which holds, for each particle, the index of the collision cell it is currently in.
[in]collisionCellMomentaAn array holding, for each collision cell, its three cartesian momentum coordinates.
[in]collisionCellMassesAn array holding, for each collision cell, its total mass content.
[in]collisionCellVelocityScalingsFor each collision cell, stores the factor by which to scale the relative velocities.

Definition at line 176 of file Simulation.cpp.

◆ getCollisionCellCenterOfMassMomentum()

Vector3D< MPCParticleVelocityType > OpenMPCD::ImplementationDetails::getCollisionCellCenterOfMassMomentum ( const unsigned int  collisionCellIndex,
const MPCParticleVelocityType *const  collisionCellMomenta 
)

Returns the center-of-mass momentum of a collision cell.

Parameters
[in]collisionCellIndexThe index of the collision cell.
[in]collisionCellMomentaThe momenta of the collision cells.

Definition at line 102 of file Simulation.cpp.

◆ getCollisionCellCenterOfMassVelocity()

Vector3D< MPCParticleVelocityType > OpenMPCD::ImplementationDetails::getCollisionCellCenterOfMassVelocity ( const unsigned int  collisionCellIndex,
const MPCParticleVelocityType *const  collisionCellMomenta,
const FP *const  collisionCellMasses 
)

Returns the center-of-mass velocity of a collision cell.

Parameters
[in]collisionCellIndexThe index of the collision cell.
[in]collisionCellMomentaThe momenta of the collision cells.
[in]collisionCellMassesThe masses of the collision cells.

Definition at line 113 of file Simulation.cpp.

◆ getCollisionCellIndex()

unsigned int OpenMPCD::ImplementationDetails::getCollisionCellIndex ( const Vector3D< MPCParticlePositionType > &  position,
const unsigned int  boxSizeX,
const unsigned int  boxSizeY,
const unsigned int  boxSizeZ 
)

Returns the collision cell index for the given position.

The collision cell index is determined as follows: The three cartesian coordinates of the given particle position are rounded down towards the nearest integer. This triple of numbers then determine the coordinates of the collision cell the particle lies in. If they are called cellX, cellY, and cellZ, respectively, what is returned is cellX + cellY * boxSizeX + cellZ * boxSizeX * boxSizeY`.

Parameters
[in]positionThe position to consider, which must lie in the primary simulation volume, i.e. the x, y, and z coordinates must be non-negative and smaller than the simulation box size along the respective direction.
[in]boxSizeXThe number of collision cells in the primary simulation volume along the x axis.
[in]boxSizeYThe number of collision cells in the primary simulation volume along the y axis.
[in]boxSizeZThe number of collision cells in the primary simulation volume along the z axis.

Definition at line 12 of file Simulation.cpp.

◆ sortIntoCollisionCellsLeesEdwards()

void OpenMPCD::ImplementationDetails::sortIntoCollisionCellsLeesEdwards ( const unsigned int  particleID,
const MPCParticlePositionType *const  gridShift_,
const FP  shearRate,
const FP  mpcTime,
const MPCParticlePositionType *const  positions,
MPCParticleVelocityType *const  velocities,
MPCParticleVelocityType *const  velocityCorrections,
unsigned int *const  collisionCellIndices,
const unsigned int  boxSizeX,
const unsigned int  boxSizeY,
const unsigned int  boxSizeZ 
)

Sorts the given particle into the collision cells, temporarily applying Lees-Edwards boundary conditions.

See also
LeesEdwardsBoundaryConditions
Parameters
[in]particleIDThe particle ID.
[in]gridShift_The grid shift vector to temporarily add to the position.
[in]shearRateThe Lees-Edwards shear rate \( \dot{\gamma} \).
[in]mpcTimeThe MPC time that has passed since the start of the simulation.
[in]positionsThe particle positions.
[in,out]velocitiesThe particle velocities, which may be changed due to Lees-Edwards boundary conditions.
[out]velocityCorrectionsAn array, long enough to store one element per particle; at position particleID, the velocity corrections along the x axis applied due to Lees-Edwards boundary conditions will be saved.
[out]collisionCellIndicesThe collision cell indices for the particles.
[in]boxSizeXThe number of collision cells in the primary simulation volume along the x axis.
[in]boxSizeYThe number of collision cells in the primary simulation volume along the y axis.
[in]boxSizeZThe number of collision cells in the primary simulation volume along the z axis.

Definition at line 36 of file Simulation.cpp.

◆ undoLeesEdwardsVelocityCorrections()

void OpenMPCD::ImplementationDetails::undoLeesEdwardsVelocityCorrections ( const unsigned int  particleCount,
MPCParticleVelocityType *const  velocities,
const MPCParticleVelocityType *const  velocityCorrections 
)

Undoes the velocity corrections applied by sortIntoCollisionCellsLeesEdwards.

Parameters
[in]particleCountThe number of particles there are.
[in,out]velocitiesThe particle velocities.
[in]velocityCorrectionsThe velocity corrections as returned by sortIntoCollisionCellsLeesEdwards.

Definition at line 206 of file Simulation.cpp.