Open Qmin
0.8.0
GPU-accelerated Q-tensor-based liquid crystal simulations
|
#include <multirankSimulation.h>
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< multirankSimulation > | getPointer () |
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< WeakUpdaterPtr > | updaters |
A vector of updaters that the simulation will loop through. More... | |
vector< WeakForcePtr > | forceComputers |
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< scalar > | dataBuffer |
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... | |
|
inline |
References determineCommunicationPattern(), basicSimulation::myRank, and setRankTopology().
|
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.
void multirankSimulation::finalizeObjects | ( | ) |
transfer buffers and make sure sites on the skin of each rank have correct type
References communicateHaloSitesRoutine(), basicSimulation::NActive, and updaters.
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().
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.
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.
void multirankSimulation::createSphericalCavity | ( | scalar3 | center, |
scalar | radius, | ||
boundaryObject & | bObj | ||
) |
make a simple sphere, setting all points farther than radius of center to be the object
References boundaryObject::boundary, createMultirankBoundaryObject(), degeneratePlanar, homeotropic, mConfiguration, norm(), boundaryObject::P1, boundaryObject::P2, qTensorFromDirector(), rankTopology, scalar, scalar3, and UNWRITTENCODE.
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.
colloidOrCapillary | if 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.
void multirankSimulation::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
References boundaryObject::boundary, createMultirankBoundaryObject(), degeneratePlanar, homeotropic, mConfiguration, norm(), boundaryObject::P1, boundaryObject::P2, pointSegmentDistance(), qTensorFromDirector(), rankTopology, scalar, scalar3, and UNWRITTENCODE.
|
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.
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.
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.
|
virtual |
handles calls to all necessary halo site transfer
References communicationDirectionParity, communicationDirections, communicationTargets, ArrayHandle< T >::data, access_location::device, access_location::host, mConfiguration, MPI_SCALAR, mpiRequests, basicSimulation::nRanks, access_mode::overwrite, access_mode::read, synchronizeAndTransferBuffers(), and basicSimulation::useGPU.
Referenced by finalizeObjects(), moveParticles(), and setConfiguration().
|
virtual |
synchronize mpi and make transfer buffers
References communicationDirections, mConfiguration, mpiRequests, mpiStatuses, basicSimulation::nRanks, transfersUpToDate, and basicSimulation::useGPU.
Referenced by communicateHaloSitesRoutine(), and computeForces().
|
inline |
return a shared pointer to this Simulation
Referenced by addUpdater().
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.
|
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.
|
inline |
Add an updater.
References updaters.
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.
|
inlinevirtual |
Add a force computer configuration.
References forceComputers.
|
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.
|
inline |
Clear out the vector of forceComputes.
References forceComputers.
|
inline |
Clear out the vector of updaters.
References updaters.
|
inline |
A utility function that just checks the first updater for a max force.
References updaters.
void multirankSimulation::performTimestep | ( | ) |
Call every updater to advance one time step.
References basicSimulation::integerTimestep, basicSimulation::integrationTimestep, basicSimulation::Time, transfersUpToDate, and updaters.
|
virtual |
manipulate data from updaters
Reimplemented from basicSimulation.
References data, dataBuffer, profiler::end(), MPI_SCALAR, basicSimulation::nRanks, p1, and profiler::start().
Referenced by computePotentialEnergy().
|
virtual |
compute the potential energy associated with all of the forces
Reimplemented from basicSimulation.
References forceComputers, basicSimulation::NActive, scalar, and sumUpdaterData().
Referenced by computeEnergy().
|
virtual |
compute the kinetic energy
Reimplemented from basicSimulation.
References mConfiguration.
Referenced by computeEnergy().
|
inlinevirtual |
compute the total energy
References computeKineticEnergy(), and computePotentialEnergy().
void multirankSimulation::setIntegrationTimestep | ( | scalar | dt | ) |
Set the simulation timestep.
References basicSimulation::integrationTimestep, and updaters.
void multirankSimulation::setCPUOperation | ( | bool | setcpu | ) |
turn on CPU-only mode for all components
References forceComputers, mConfiguration, updaters, and basicSimulation::useGPU.
void multirankSimulation::setReproducible | ( | bool | reproducible | ) |
Enforce reproducible dynamics.
References updaters.
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.
|
protected |
References Index3D::inverseIndex(), basicSimulation::myRank, basicSimulation::nRanks, parityTest, rankParity, and rankTopology.
Referenced by multirankSimulation().
|
protected |
References communicationDirectionParity, communicationDirections, communicationTargets, corners, edges, mpiRequests, mpiStatuses, parityTest, rankParity, rankTopology, and Index3D::sizes.
Referenced by multirankSimulation().
|
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.
|
inlineinherited |
Set the time between spatial sorting operations.
References basicSimulation::sortPeriod.
|
virtualinherited |
reset the simulation clock
References basicSimulation::Time.
|
inlinevirtualinherited |
reset the simulation clock counter
References basicSimulation::integerTimestep.
Referenced by moveParticles(), and sumUpdaterData().
WeakMConfigPtr multirankSimulation::mConfiguration |
The configuration of latticeSites.
Referenced by communicateHaloSitesRoutine(), computeForces(), computeKineticEnergy(), createCylindricalObject(), createMultirankBoundaryObject(), createSphericalCavity(), createSpherocylinder(), createWall(), moveParticles(), saveState(), setConfiguration(), setCPUOperation(), setDipolarField(), and synchronizeAndTransferBuffers().
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().
vector<WeakForcePtr> multirankSimulation::forceComputers |
A vector of force computes the simulation will loop through.
Referenced by addForce(), clearForceComputers(), computeForces(), computePotentialEnergy(), and setCPUOperation().
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().
|
protected |
Referenced by communicateHaloSitesRoutine(), determineCommunicationPattern(), and synchronizeAndTransferBuffers().
|
protected |
Referenced by communicateHaloSitesRoutine(), and determineCommunicationPattern().
|
protected |
Referenced by communicateHaloSitesRoutine(), and determineCommunicationPattern().
|
protected |
the number of ranks per {x,y,z} axis
Referenced by createCylindricalObject(), createMultirankBoundaryObject(), createSphericalCavity(), createSpherocylinder(), createWall(), determineCommunicationPattern(), setDipolarField(), and setRankTopology().
|
protected |
the local {x,y,z} rank coordinate
Referenced by createMultirankBoundaryObject(), determineCommunicationPattern(), saveState(), setConfiguration(), setDipolarField(), and setRankTopology().
|
protected |
Referenced by determineCommunicationPattern(), and setRankTopology().
|
protected |
do edges need to be communicated?
Referenced by determineCommunicationPattern().
|
protected |
do corner sites?
Referenced by determineCommunicationPattern().
|
protected |
have the halo sites been communicated?
Referenced by moveParticles(), performTimestep(), and synchronizeAndTransferBuffers().
|
protected |
|
protected |
Referenced by determineCommunicationPattern(), and synchronizeAndTransferBuffers().
|
protected |
Referenced by communicateHaloSitesRoutine(), determineCommunicationPattern(), and synchronizeAndTransferBuffers().
|
protected |
Referenced by sumUpdaterData().
|
inherited |
The configuration of particles.
Referenced by Simulation::computeForces(), Simulation::computeKineticEnergy(), Simulation::moveParticles(), basicSimulation::setBox(), Simulation::setConfiguration(), Simulation::setCPUOperation(), and Simulation::setNThreads().
|
inherited |
The domain of the simulation.
Referenced by basicSimulation::basicSimulation(), basicSimulation::setBox(), Simulation::setConfiguration(), and setConfiguration().
|
inherited |
An integer that keeps track of how often performTimestep has been called.
Referenced by Simulation::performTimestep(), performTimestep(), and basicSimulation::setCurrentTimestep().
|
inherited |
The current simulation time.
Referenced by Simulation::performTimestep(), performTimestep(), and basicSimulation::setCurrentTime().
|
inherited |
The dt of a time step.
Referenced by Simulation::performTimestep(), performTimestep(), Simulation::setIntegrationTimestep(), and setIntegrationTimestep().
|
inherited |
A flag controlling whether to use the GPU.
Referenced by communicateHaloSitesRoutine(), computeForces(), Simulation::setCPUOperation(), setCPUOperation(), and synchronizeAndTransferBuffers().
|
inherited |
integer for this rank
Referenced by multirankSimulation(), and setRankTopology().
|
inherited |
total number of ranks
Referenced by communicateHaloSitesRoutine(), setRankTopology(), sumUpdaterData(), and synchronizeAndTransferBuffers().
|
inherited |
some measure of the number of active degrees of freedom
Referenced by computePotentialEnergy(), and finalizeObjects().
|
protectedinherited |
Determines how frequently the spatial sorter be called...once per sortPeriod Timesteps. When sortPeriod < 0 no sorting occurs.
Referenced by basicSimulation::setSortPeriod().
|
protectedinherited |
A flag that determins if a spatial sorting is due to occur this Timestep.