Open Qmin  0.8.0
GPU-accelerated Q-tensor-based liquid crystal simulations
Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | List of all members
multirankSimulation Class Reference

#include <multirankSimulation.h>

Inheritance diagram for multirankSimulation:
Inheritance graph
[legend]
Collaboration diagram for multirankSimulation:
Collaboration graph
[legend]

Public Member Functions

 multirankSimulation (int _myRank, int xDiv, int yDiv, int zDiv, bool _edges, bool _corners)
 
virtual void moveParticles (GPUArray< dVec > &displacements, scalar scale=1.0)
 move particles, and also communicate halo sites More...
 
void finalizeObjects ()
 transfer buffers and make sure sites on the skin of each rank have correct type More...
 
void createMultirankBoundaryObject (vector< int3 > &latticeSites, vector< dVec > &qTensors, boundaryType _type, scalar Param1, scalar Param2)
 the flexible base function...given lattice sites composing an object, determine which are on this rank and add the object More...
 
void createWall (int xyz, int plane, boundaryObject &bObj)
 make a wall with x, y, or z normal More...
 
void createSphericalColloid (scalar3 center, scalar radius, boundaryObject &bObj)
 make a simple sphere, setting all points within radius of center to be the object More...
 
void createSphericalCavity (scalar3 center, scalar radius, boundaryObject &bObj)
 make a simple sphere, setting all points farther than radius of center to be the object More...
 
void createCylindricalObject (scalar3 cylinderStart, scalar3 cylinderEnd, scalar radius, bool colloidOrCapillary, boundaryObject &bObj)
 make a cylindrical object, with either the inside or outside defined as the object More...
 
void createSpherocylinder (scalar3 cylinderStart, scalar3 cylinderEnd, scalar radius, boundaryObject &bObj)
 make a spherocylinder, defined by the start and end of the cylindrical section and the radius More...
 
virtual void createBoundaryFromFile (string fname, bool verbose=false)
 import a boundary object from a (carefully prepared) text file More...
 
void setDipolarField (scalar3 center, scalar ThetaD, scalar radius, scalar range, scalar S0)
 a function of convenience... make a dipolar field a la Lubensky et al. More...
 
void setDipolarField (scalar3 center, scalar3 direction, scalar radius, scalar range, scalar S0)
 a function of convenience... make a dipolar field a la Lubensky et al. Simpler, uses the ravnik expression More...
 
virtual void communicateHaloSitesRoutine ()
 handles calls to all necessary halo site transfer More...
 
virtual void synchronizeAndTransferBuffers ()
 synchronize mpi and make transfer buffers More...
 
shared_ptr< multirankSimulationgetPointer ()
 return a shared pointer to this Simulation More...
 
void setConfiguration (MConfigPtr _config)
 Pass in a reference to the configuration. More...
 
virtual void computeForces ()
 Call the force computer to compute the forces. More...
 
void addUpdater (UpdaterPtr _upd)
 Add an updater. More...
 
void addUpdater (UpdaterPtr _upd, MConfigPtr _config)
 Add an updater with a reference to a configuration. More...
 
virtual void addForce (ForcePtr _force)
 Add a force computer configuration. More...
 
virtual void addForce (ForcePtr _force, MConfigPtr _config)
 Add a force computer configuration. More...
 
void clearForceComputers ()
 Clear out the vector of forceComputes. More...
 
void clearUpdaters ()
 Clear out the vector of updaters. More...
 
scalar getMaxForce ()
 A utility function that just checks the first updater for a max force. More...
 
void performTimestep ()
 Call every updater to advance one time step. More...
 
virtual void sumUpdaterData (vector< scalar > &data)
 manipulate data from updaters More...
 
virtual scalar computePotentialEnergy (bool verbose=false)
 compute the potential energy associated with all of the forces More...
 
virtual scalar computeKineticEnergy (bool verbose=false)
 compute the kinetic energy More...
 
virtual scalar computeEnergy (bool verbose=false)
 compute the total energy More...
 
void setIntegrationTimestep (scalar dt)
 Set the simulation timestep. More...
 
void setCPUOperation (bool setcpu)
 turn on CPU-only mode for all components More...
 
void setReproducible (bool reproducible)
 Enforce reproducible dynamics. More...
 
void saveState (string fname, int latticeSkip=1, int defectType=0)
 save a file for each rank recording the expanded lattice; lattice skip controls the sparsity of saved sites More...
 
void setBox (BoxPtr _box)
 This changes the contents of the Box pointed to by Box to match that of _box. More...
 
void setSortPeriod (int sp)
 Set the time between spatial sorting operations. More...
 
virtual void setCurrentTime (scalar _cTime)
 reset the simulation clock More...
 
virtual void setCurrentTimestep (int _cTime)
 reset the simulation clock counter More...
 

Public Attributes

profiler p1 = profiler("total communication time")
 
profiler p4 = profiler("MPI recv")
 
profiler p3 = profiler("MPI send")
 
profiler p2 = profiler("GPU data buffering kernels time")
 
WeakMConfigPtr mConfiguration
 The configuration of latticeSites. More...
 
vector< WeakUpdaterPtrupdaters
 A vector of updaters that the simulation will loop through. More...
 
vector< WeakForcePtrforceComputers
 A vector of force computes the simulation will loop through. More...
 
int3 latticeMinPosition
 in multi-rank simulations, this stores the lowest (x,y,z) coordinate controlled by the current rank More...
 
WeakConfigPtr configuration
 The configuration of particles. More...
 
BoxPtr Box
 The domain of the simulation. More...
 
int integerTimestep
 An integer that keeps track of how often performTimestep has been called. More...
 
scalar Time
 The current simulation time. More...
 
scalar integrationTimestep
 The dt of a time step. More...
 
bool useGPU
 A flag controlling whether to use the GPU. More...
 
int myRank
 integer for this rank More...
 
int nRanks
 total number of ranks More...
 
int NActive = 0
 some measure of the number of active degrees of freedom More...
 

Protected Member Functions

void setRankTopology (int x, int y, int z)
 
void determineCommunicationPattern (bool _edges, bool _corners)
 

Protected Attributes

vector< int2 > communicationDirections
 
vector< bool > communicationDirectionParity
 
vector< int > communicationTargets
 
int3 rankTopology
 the number of ranks per {x,y,z} axis More...
 
int3 rankParity
 the local {x,y,z} rank coordinate More...
 
Index3D parityTest
 
bool edges
 do edges need to be communicated? More...
 
bool corners
 do corner sites? More...
 
bool transfersUpToDate
 have the halo sites been communicated? More...
 
MPI_Status mpiStatus
 
vector< MPI_Status > mpiStatuses
 
vector< MPI_Request > mpiRequests
 
vector< scalardataBuffer
 
int sortPeriod
 Determines how frequently the spatial sorter be called...once per sortPeriod Timesteps. When sortPeriod < 0 no sorting occurs. More...
 
bool spatialSortThisStep
 A flag that determins if a spatial sorting is due to occur this Timestep. More...
 

Constructor & Destructor Documentation

◆ multirankSimulation()

multirankSimulation::multirankSimulation ( int  _myRank,
int  xDiv,
int  yDiv,
int  zDiv,
bool  _edges,
bool  _corners 
)
inline

Member Function Documentation

◆ moveParticles()

void multirankSimulation::moveParticles ( GPUArray< dVec > &  displacements,
scalar  scale = 1.0 
)
virtual

move particles, and also communicate halo sites

Calls the configuration to displace the degrees of freedom, and communicates halo sites according to the rankTopology and boolean settings

Implements basicSimulation.

References communicateHaloSitesRoutine(), profiler::end(), mConfiguration, p1, profiler::start(), and transfersUpToDate.

◆ finalizeObjects()

void multirankSimulation::finalizeObjects ( )

transfer buffers and make sure sites on the skin of each rank have correct type

References communicateHaloSitesRoutine(), basicSimulation::NActive, and updaters.

◆ createMultirankBoundaryObject()

void multirankSimulation::createMultirankBoundaryObject ( vector< int3 > &  latticeSites,
vector< dVec > &  qTensors,
boundaryType  _type,
scalar  Param1,
scalar  Param2 
)

the flexible base function...given lattice sites composing an object, determine which are on this rank and add the object

References latticeMinPosition, mConfiguration, rankParity, rankTopology, and wrap().

Referenced by createBoundaryFromFile(), createCylindricalObject(), createSphericalCavity(), createSphericalColloid(), createSpherocylinder(), and createWall().

◆ createWall()

void multirankSimulation::createWall ( int  xyz,
int  plane,
boundaryObject bObj 
)

make a wall with x, y, or z normal

/param xyz int that is 0, 1, or 2 for wall normal x, y, and z, respectively

References boundaryObject::boundary, createMultirankBoundaryObject(), degeneratePlanar, homeotropic, mConfiguration, boundaryObject::P1, boundaryObject::P2, rankTopology, scalar, and UNWRITTENCODE.

◆ createSphericalColloid()

void multirankSimulation::createSphericalColloid ( scalar3  center,
scalar  radius,
boundaryObject bObj 
)

make a simple sphere, setting all points within radius of center to be the object

References boundaryObject::boundary, createMultirankBoundaryObject(), degeneratePlanar, homeotropic, norm(), boundaryObject::P1, boundaryObject::P2, qTensorFromDirector(), scalar, scalar3, and UNWRITTENCODE.

◆ createSphericalCavity()

void multirankSimulation::createSphericalCavity ( scalar3  center,
scalar  radius,
boundaryObject bObj 
)

◆ createCylindricalObject()

void multirankSimulation::createCylindricalObject ( scalar3  cylinderStart,
scalar3  cylinderEnd,
scalar  radius,
bool  colloidOrCapillary,
boundaryObject bObj 
)

make a cylindrical object, with either the inside or outside defined as the object

define a cylinder by its two endpoints and its prependicular radius.

Parameters
colloidOrCapillaryif true, points inside the radius are the object (i.e., a colloid), else points outside are the object (i.e., a surrounding capillary)

References boundaryObject::boundary, createMultirankBoundaryObject(), degeneratePlanar, homeotropic, mConfiguration, norm(), boundaryObject::P1, boundaryObject::P2, qTensorFromDirector(), rankTopology, scalar, scalar3, truncatedPointSegmentDistance(), and UNWRITTENCODE.

◆ createSpherocylinder()

void multirankSimulation::createSpherocylinder ( scalar3  cylinderStart,
scalar3  cylinderEnd,
scalar  radius,
boundaryObject bObj 
)

◆ createBoundaryFromFile()

void multirankSimulation::createBoundaryFromFile ( string  fname,
bool  verbose = false 
)
virtual

import a boundary object from a (carefully prepared) text file

Reads a carefully prepared text file to create a new boundary object... For convenience this information is copied into the main README.md documentation file, so both should be updated together if needed.

The first line must be a single integer specifying the number of objects to be read in. This is meant in the colloquial English sense, not in the zero-indexed counting sense. So, if you want your file to specify one object, make sure the first line is the number 1.

Each subsequent block must be formatted as follows (with no additional lines between the blocks): The first line MUST be formatted as a b c d where a=0 means oriented anchoring (such as homeotropic or oriented planar), a=1 means degenerate Planar b is a scalar setting the anchoring strength c is the preferred value of S0 d is an integer specifying the number of sites, N_B.

Subsequently, there MUST be N_b lines of the form x y z C1 C2 C3 C4 C5, where x, y, and z are the integer lattice sites, and C1, C2, C3, C4, C5 are real numbers corresponding to the desired anchoring conditions: For oriented anchoring, C1, C2, C3, C4, C5 correspond directly to the surface-preferred Q-tensor: C1 = Qxx, C2 = Qxy, C3=Qxz, C4 = Qyy, C5=Qyz, where one often will set the Q-tensor by choosing a locally preferred director, \nu^s, and setting Q^B = 3 S_0/2*(\nu^s \nu^s - \delta_{ab}/3).

For degenerate planar anchoring the five constants should be specified as, C1 = \hat{nu}_x C2 = \hat{nu}_y C3 = \hat{nu}_z C4 = 0.0 C5 = 0.0, where \nu^s = {Cos[[Phi]] Sin[[theta]], Sin[[Phi]] Sin[[theta]], Cos[[theta]]} is the direction to which the LC should try to be orthogonal

References createMultirankBoundaryObject(), degeneratePlanar, homeotropic, and scalar.

◆ setDipolarField() [1/2]

void multirankSimulation::setDipolarField ( scalar3  center,
scalar  ThetaD,
scalar  radius,
scalar  range,
scalar  S0 
)

a function of convenience... make a dipolar field a la Lubensky et al.

for lattice sites within range of center, create a z-aligned dipolar field with opening angle ThetaD about an object of radius "radius" Functional form from Lubensky, Pettey, Currier, and Stark. PRE 57, 610, 1998

References Cos, latticeMinPosition, mConfiguration, norm(), PI, qTensorFromDirector(), rankParity, rankTopology, scalar, scalar3, and Sin.

◆ setDipolarField() [2/2]

void multirankSimulation::setDipolarField ( scalar3  center,
scalar3  direction,
scalar  radius,
scalar  range,
scalar  S0 
)

a function of convenience... make a dipolar field a la Lubensky et al. Simpler, uses the ravnik expression

for lattice sites within range of center, create a dipolar field in the direction specified about an object of radius "radius"

References latticeMinPosition, mConfiguration, norm(), qTensorFromDirector(), rankParity, rankTopology, scalar, and scalar3.

◆ communicateHaloSitesRoutine()

void multirankSimulation::communicateHaloSitesRoutine ( )
virtual

◆ synchronizeAndTransferBuffers()

void multirankSimulation::synchronizeAndTransferBuffers ( )
virtual

◆ getPointer()

shared_ptr<multirankSimulation> multirankSimulation::getPointer ( )
inline

return a shared pointer to this Simulation

Referenced by addUpdater().

◆ setConfiguration()

void multirankSimulation::setConfiguration ( MConfigPtr  _config)

Pass in a reference to the configuration.

Set a pointer to the configuration

References basicSimulation::Box, communicateHaloSitesRoutine(), latticeMinPosition, mConfiguration, and rankParity.

◆ computeForces()

void multirankSimulation::computeForces ( )
virtual

Call the force computer to compute the forces.

Calls all force computers, and evaluate the self force calculation if the model demands it

Implements basicSimulation.

References forceComputers, mConfiguration, synchronizeAndTransferBuffers(), and basicSimulation::useGPU.

◆ addUpdater() [1/2]

void multirankSimulation::addUpdater ( UpdaterPtr  _upd)
inline

Add an updater.

References updaters.

◆ addUpdater() [2/2]

void multirankSimulation::addUpdater ( UpdaterPtr  _upd,
MConfigPtr  _config 
)

Add an updater with a reference to a configuration.

Add a pointer to the list of updaters, and give that updater a reference to the model...

References getPointer(), and updaters.

◆ addForce() [1/2]

virtual void multirankSimulation::addForce ( ForcePtr  _force)
inlinevirtual

Add a force computer configuration.

References forceComputers.

◆ addForce() [2/2]

void multirankSimulation::addForce ( ForcePtr  _force,
MConfigPtr  _config 
)
virtual

Add a force computer configuration.

Add a pointer to the list of force computers, and give that FC a reference to the model...

References forceComputers.

◆ clearForceComputers()

void multirankSimulation::clearForceComputers ( )
inline

Clear out the vector of forceComputes.

References forceComputers.

◆ clearUpdaters()

void multirankSimulation::clearUpdaters ( )
inline

Clear out the vector of updaters.

References updaters.

◆ getMaxForce()

scalar multirankSimulation::getMaxForce ( )
inline

A utility function that just checks the first updater for a max force.

References updaters.

◆ performTimestep()

void multirankSimulation::performTimestep ( )

◆ sumUpdaterData()

void multirankSimulation::sumUpdaterData ( vector< scalar > &  data)
virtual

manipulate data from updaters

Reimplemented from basicSimulation.

References data, dataBuffer, profiler::end(), MPI_SCALAR, basicSimulation::nRanks, p1, and profiler::start().

Referenced by computePotentialEnergy().

◆ computePotentialEnergy()

scalar multirankSimulation::computePotentialEnergy ( bool  verbose = false)
virtual

compute the potential energy associated with all of the forces

Reimplemented from basicSimulation.

References forceComputers, basicSimulation::NActive, scalar, and sumUpdaterData().

Referenced by computeEnergy().

◆ computeKineticEnergy()

scalar multirankSimulation::computeKineticEnergy ( bool  verbose = false)
virtual

compute the kinetic energy

Reimplemented from basicSimulation.

References mConfiguration.

Referenced by computeEnergy().

◆ computeEnergy()

virtual scalar multirankSimulation::computeEnergy ( bool  verbose = false)
inlinevirtual

compute the total energy

References computeKineticEnergy(), and computePotentialEnergy().

◆ setIntegrationTimestep()

void multirankSimulation::setIntegrationTimestep ( scalar  dt)

Set the simulation timestep.

Postcondition
the cell configuration and e.o.m. timestep is set to the input value

References basicSimulation::integrationTimestep, and updaters.

◆ setCPUOperation()

void multirankSimulation::setCPUOperation ( bool  setcpu)

turn on CPU-only mode for all components

Postcondition
the cell configuration and e.o.m. timestep is set to the input value

References forceComputers, mConfiguration, updaters, and basicSimulation::useGPU.

◆ setReproducible()

void multirankSimulation::setReproducible ( bool  reproducible)

Enforce reproducible dynamics.

Precondition
the updaters already know if the GPU will be used
Postcondition
the updaters are set to be reproducible if the boolean is true, otherwise the RNG is initialized

References updaters.

◆ saveState()

void multirankSimulation::saveState ( string  fname,
int  latticeSkip = 1,
int  defectType = 0 
)

save a file for each rank recording the expanded lattice; lattice skip controls the sparsity of saved sites

Saves a file for each rank recording the current state of the system (each rank will produce an independent file, named StringJoin[fname,"_x(ToString[X])_y(ToString[Y])_z(ToString[Z]).txt"], where X, Y, and Z denote the position in the partitioning that the rank controls.

The latticeSkip parameter thins out the number of lattice sites that are saved (crucial for storing the results of large simulations). In the default setting of latticeSkip=1, every site is saved. Otherwise, only sites (x_i,y_i,j_i) will be saved, where (i mod(latticeSkip) == 0)

The format of the saved state is a simple text file where every line is the following: i j k Qxx Qxy Qxz Qyy Qyz type defectValue where i, j, and k are the coordinates of the lattice site (in the global frame), followed by the five independent components of the Q tensor, the "type" of the site, and the value computed by computeDefectMeasures(defectType). By default, this is just the maximum eigenvalue of the Q-tensor at the indicated site; see the documentation of computeDefectMeasures for more information.

References ArrayHandle< T >::data, access_location::host, idx, mConfiguration, rankParity, and access_mode::read.

◆ setRankTopology()

void multirankSimulation::setRankTopology ( int  x,
int  y,
int  z 
)
protected

◆ determineCommunicationPattern()

void multirankSimulation::determineCommunicationPattern ( bool  _edges,
bool  _corners 
)
protected

◆ setBox()

void basicSimulation::setBox ( BoxPtr  _box)
inherited

This changes the contents of the Box pointed to by Box to match that of _box.

Set a new Box for the simulation...This is the function that should be called to propagate a change in the box dimensions throughout the simulation...

References basicSimulation::Box, and basicSimulation::configuration.

◆ setSortPeriod()

void basicSimulation::setSortPeriod ( int  sp)
inlineinherited

Set the time between spatial sorting operations.

References basicSimulation::sortPeriod.

◆ setCurrentTime()

void basicSimulation::setCurrentTime ( scalar  _cTime)
virtualinherited

reset the simulation clock

Postcondition
the cell configuration and e.o.m. timestep is set to the input value

References basicSimulation::Time.

◆ setCurrentTimestep()

virtual void basicSimulation::setCurrentTimestep ( int  _cTime)
inlinevirtualinherited

reset the simulation clock counter

References basicSimulation::integerTimestep.

Member Data Documentation

◆ p1

profiler multirankSimulation::p1 = profiler("total communication time")

Referenced by moveParticles(), and sumUpdaterData().

◆ p4

profiler multirankSimulation::p4 = profiler("MPI recv")

◆ p3

profiler multirankSimulation::p3 = profiler("MPI send")

◆ p2

profiler multirankSimulation::p2 = profiler("GPU data buffering kernels time")

◆ mConfiguration

WeakMConfigPtr multirankSimulation::mConfiguration

◆ updaters

vector<WeakUpdaterPtr> multirankSimulation::updaters

A vector of updaters that the simulation will loop through.

Referenced by addUpdater(), clearUpdaters(), finalizeObjects(), getMaxForce(), performTimestep(), setCPUOperation(), setIntegrationTimestep(), and setReproducible().

◆ forceComputers

vector<WeakForcePtr> multirankSimulation::forceComputers

A vector of force computes the simulation will loop through.

Referenced by addForce(), clearForceComputers(), computeForces(), computePotentialEnergy(), and setCPUOperation().

◆ latticeMinPosition

int3 multirankSimulation::latticeMinPosition

in multi-rank simulations, this stores the lowest (x,y,z) coordinate controlled by the current rank

Referenced by createMultirankBoundaryObject(), setConfiguration(), and setDipolarField().

◆ communicationDirections

vector<int2> multirankSimulation::communicationDirections
protected

◆ communicationDirectionParity

vector<bool> multirankSimulation::communicationDirectionParity
protected

◆ communicationTargets

vector<int> multirankSimulation::communicationTargets
protected

◆ rankTopology

int3 multirankSimulation::rankTopology
protected

◆ rankParity

int3 multirankSimulation::rankParity
protected

◆ parityTest

Index3D multirankSimulation::parityTest
protected

◆ edges

bool multirankSimulation::edges
protected

do edges need to be communicated?

Referenced by determineCommunicationPattern().

◆ corners

bool multirankSimulation::corners
protected

do corner sites?

Referenced by determineCommunicationPattern().

◆ transfersUpToDate

bool multirankSimulation::transfersUpToDate
protected

have the halo sites been communicated?

Referenced by moveParticles(), performTimestep(), and synchronizeAndTransferBuffers().

◆ mpiStatus

MPI_Status multirankSimulation::mpiStatus
protected

◆ mpiStatuses

vector<MPI_Status> multirankSimulation::mpiStatuses
protected

◆ mpiRequests

vector<MPI_Request> multirankSimulation::mpiRequests
protected

◆ dataBuffer

vector<scalar> multirankSimulation::dataBuffer
protected

Referenced by sumUpdaterData().

◆ configuration

WeakConfigPtr basicSimulation::configuration
inherited

◆ Box

BoxPtr basicSimulation::Box
inherited

◆ integerTimestep

int basicSimulation::integerTimestep
inherited

An integer that keeps track of how often performTimestep has been called.

Referenced by Simulation::performTimestep(), performTimestep(), and basicSimulation::setCurrentTimestep().

◆ Time

scalar basicSimulation::Time
inherited

The current simulation time.

Referenced by Simulation::performTimestep(), performTimestep(), and basicSimulation::setCurrentTime().

◆ integrationTimestep

scalar basicSimulation::integrationTimestep
inherited

◆ useGPU

bool basicSimulation::useGPU
inherited

◆ myRank

int basicSimulation::myRank
inherited

integer for this rank

Referenced by multirankSimulation(), and setRankTopology().

◆ nRanks

int basicSimulation::nRanks
inherited

◆ NActive

int basicSimulation::NActive = 0
inherited

some measure of the number of active degrees of freedom

Referenced by computePotentialEnergy(), and finalizeObjects().

◆ sortPeriod

int basicSimulation::sortPeriod
protectedinherited

Determines how frequently the spatial sorter be called...once per sortPeriod Timesteps. When sortPeriod < 0 no sorting occurs.

Referenced by basicSimulation::setSortPeriod().

◆ spatialSortThisStep

bool basicSimulation::spatialSortThisStep
protectedinherited

A flag that determins if a spatial sorting is due to occur this Timestep.


The documentation for this class was generated from the following files: