CellGPU
0.8.0
GPU-accelerated simulations of cells
|
Implement a 2D Voronoi model, with and without some extra bells and whistles, using kernels in SPV Kernels. More...
#include <voronoiQuadraticEnergy.h>
Public Member Functions | |
VoronoiQuadraticEnergy (int n, bool reprod=false) | |
initialize with random positions in a square box More... | |
VoronoiQuadraticEnergy (int n, Dscalar A0, Dscalar P0, bool reprod=false) | |
initialize with random positions and set all cells to have uniform target A_0 and P_0 parameters More... | |
VoronoiQuadraticEnergy () | |
Blank constructor. | |
void | initializeVoronoiQuadraticEnergy (int n) |
Initialize voronoiQuadraticEnergy and call the initializer chain. More... | |
virtual void | computeForces () |
compute the geometry and get the forces More... | |
virtual Dscalar | computeEnergy () |
compute the quadratic energy functional More... | |
virtual void | ComputeForceSetsGPU () |
Compute force sets on the GPU. More... | |
void | SumForcesGPU () |
Add up the force sets to get the net force per particle on the GPU. More... | |
virtual void | computeVoronoiForceCPU (int i) |
Compute the net force on particle i on the CPU. More... | |
virtual void | computeVoronoiForceSetsGPU () |
call gpu_force_sets kernel caller More... | |
void | sumForceSets () |
call gpu_sum_force_sets kernel caller More... | |
void | sumForceSetsWithExclusions () |
call gpu_sum_force_sets_with_exclusions kernel caller More... | |
void | reportCellInfo () |
Report various cell infor for testing and debugging. More... | |
void | reportForces (bool verbose) |
Report information about net forces... | |
virtual void | getDynMatEntries (vector< int2 > &rcs, vector< Dscalar > &vals, Dscalar unstress=1.0, Dscalar stress=1.0) |
Save tuples for half of the dynamical matrix. More... | |
virtual Dscalar | getSigmaXY () |
calculate the current global off-diagonal stress More... | |
void | initializeVoronoiModelBase (int n) |
A default initialization scheme. More... | |
void | setCPU (bool global=true) |
Enforce CPU-only operation. More... | |
virtual void | setCPU () |
Enforce CPU-only operation. Derived classes might have to do more work when the CPU mode is invoked. | |
void | writeTriangulation (ofstream &outfile) |
write triangulation to text file | |
void | readTriangulation (ifstream &infile) |
read positions from text file...for debugging | |
virtual void | enforceTopology () |
update/enforce the topology More... | |
void | setExclusions (vector< int > &exes) |
Declare which particles are to be excluded (exes[i]!=0) More... | |
virtual int | getNumberOfDegreesOfFreedom () |
In voronoi models the number of degrees of freedom is the number of cells. | |
virtual void | moveDegreesOfFreedom (GPUArray< Dscalar2 > &displacements, Dscalar scale=1.) |
moveDegrees of Freedom calls either the move points or move points CPU routines More... | |
virtual void | getForces (GPUArray< Dscalar2 > &forces) |
return the forces | |
virtual GPUArray< Dscalar2 > & | returnForces () |
return a reference to the GPUArray of the current forces | |
virtual void | computeGeometryCPU () |
Compute cell geometry on the CPU. More... | |
virtual void | computeGeometryGPU () |
call gpu_compute_geometry kernel caller More... | |
virtual void | cellDivision (const vector< int > ¶meters, const vector< Dscalar > &dParams) |
allow for cell division, according to a vector of model-dependent parameters More... | |
virtual void | cellDeath (int cellIndex) |
Kill the indexed cell by simply removing it from the simulation. More... | |
void | movePoints (GPUArray< Dscalar2 > &displacements, Dscalar scale) |
move particles on the GPU More... | |
void | movePointsCPU (GPUArray< Dscalar2 > &displacements, Dscalar scale) |
move particles on the CPU More... | |
void | resetDelLocPoints () |
Transfer particle data from the GPU to the CPU for use by delLoc. More... | |
void | updateCellList () |
Update the cell list structure after particles have moved. More... | |
void | updateNeighIdxs () |
update the NieghIdxs data structure More... | |
void | initializeSimple2DActiveCell (int n) |
initialize class' data structures and set default values More... | |
void | setv0Dr (Dscalar v0new, Dscalar drnew) |
Set uniform motility. More... | |
void | setCellMotility (vector< Dscalar > &v0s, vector< Dscalar > &drs) |
Set non-uniform cell motilites. More... | |
void | setCellDirectorsRandomly () |
Set random cell directors (for active cell models) More... | |
Dscalar | vicsekOrderParameter (Dscalar2 &vParallel, Dscalar2 &vPerpendicular) |
measure the viscek order parameter N^-1 {v_i}{|v_i} | |
Dscalar | vicsekOrderParameterDirector (Dscalar2 &vParallel, Dscalar2 &vPerpendicular) |
measure the viscek order parameter N^-1 {v_i}{|v_i} from the director only | |
void | initializeSimple2DCell (int n) |
initialize class' data structures and set default values More... | |
virtual void | setGPU () |
Enforce GPU-only operation. This is the default mode, so this method need not be called most of the time. | |
virtual void | computeGeometry () |
call either the computeGeometryCPU or GPU routines for the current model More... | |
Dscalar | computeKineticEnergy () |
Call masses and velocities to get the total kinetic energy. More... | |
Dscalar4 | computeKineticPressure () |
Call masses and velocities to get the average kinetic contribution to the pressure tensor. More... | |
void | setCellPreferencesUniform (Dscalar A0, Dscalar P0) |
Set uniform cell area and perimeter preferences. More... | |
void | setCellPreferences (vector< Dscalar2 > &AreaPeriPreferences) |
Set cell area and perimeter preferences according to input vector. More... | |
void | setCellPositionsRandomly () |
Set random cell positions, and set the periodic box to a square with average cell area=1. More... | |
void | setCellPositions (vector< Dscalar2 > newCellPositions) |
Set cell positions according to a user-specified vector. More... | |
void | setVertexPositions (vector< Dscalar2 > newVertexPositions) |
Set vertex positions according to a user-specified vector. More... | |
Dscalar | setCellVelocitiesMaxwellBoltzmann (Dscalar T) |
Set velocities via a temperature. The return value is the total kinetic energy. More... | |
Dscalar | setVertexVelocitiesMaxwellBoltzmann (Dscalar T) |
Set velocities via a temperature for the vertex degrees of freedom. More... | |
void | setModuliUniform (Dscalar newKA, Dscalar newKP) |
set uniform moduli for all cells More... | |
void | setCellTypeUniform (int i) |
Set all cells to the same "type". More... | |
void | setCellType (vector< int > &types) |
Set cells to different "type". More... | |
void | setVertexTopologyFromCells (vector< vector< int > > cellVertexIndices) |
An uncomfortable function to allow the user to set vertex topology "by hand". More... | |
virtual gpubox & | returnBox () |
return the gpubox | |
void | setBox (BoxPtr _box) |
This can be used, but should not normally be. This re-assigns the pointer. | |
virtual vector< int > & | returnItt () |
return the base "itt" re-indexing vector | |
virtual GPUArray< Dscalar2 > & | returnModuli () |
Return a reference to moduli. | |
virtual GPUArray< Dscalar2 > & | returnAreaPeri () |
Return a reference to AreaPeri array. | |
virtual GPUArray< Dscalar2 > & | returnAreaPeriPreferences () |
Return a reference to AreaPeriPreferences. | |
virtual GPUArray< Dscalar2 > & | returnVelocities () |
Return a reference to velocities on cells. VertexModelBase will instead return vertexVelocities. | |
virtual GPUArray< Dscalar2 > & | returnPositions () |
Return a reference to Positions on cells. | |
virtual GPUArray< Dscalar > & | returnMasses () |
Return a reference to Masses on cells. | |
virtual GPUArray< Dscalar > & | returnOtherData () |
Return other data just returns the masses; in this class it's not needed. | |
void | setDeltaT (Dscalar dt) |
Set the simulation time stepsize. | |
void | getCellNeighs (int idx, int &nNeighs, vector< int > &neighs) |
Dscalar | getMaxForce () |
Get the maximum force on a cell. | |
void | reportMeanCellForce (bool verbose) |
Report the current average force on each cell. More... | |
void | reportMeanVertexForce (bool verbose=false) |
Report the current average force per vertex...should be close to zero. | |
void | reportAP (bool verbose=false) |
report the current total area, and optionally the area and perimeter for each cell | |
Dscalar | reportq () |
Report the average value of p/sqrt(A) for the cells in the system. More... | |
Dscalar | reportVarq () |
Report the variance of p/sqrt(A) for the cells in the system. More... | |
Dscalar2 | reportVarAP () |
Report the variance of A and P for the cells in the system. More... | |
Dscalar | reportMeanP () |
Report the mean value of the perimeter. More... | |
virtual void | setTime (Dscalar time) |
set the time | |
Public Attributes | |
DelaunayLoc | delLoc |
The class' local Delaunay tester/updater. | |
Dscalar | repPerFrame |
Collect statistics of how many triangulation repairs are done per frame, etc. | |
int | skippedFrames |
How often were all circumcenters empty (so that no data transfers and no repairs were necessary)? | |
int | GlobalFixes |
How often were global re-triangulations performed? | |
GPUArray< Dscalar2 > | external_forces |
"exclusions" zero out the force on a cell...the external force needed to do this is stored in external_forces | |
GPUArray< int > | exclusions |
An array containing the indices of excluded particles. | |
int | localTopologyUpdates |
The number of topology updates performed at the individual particle level. | |
GPUArray< Dscalar > | cellDirectors |
An array of angles (relative to the x-axis) that the cell directors point. | |
GPUArray< Dscalar > | cellDirectorForces |
An array of forces acting on the cell directors. | |
Dscalar | v0 |
velocity of cells in mono-motile systems | |
Dscalar | Dr |
rotational diffusion of cell directors in mono-motile systems | |
GPUArray< Dscalar2 > | Motility |
The motility parameters (v0 and Dr) for each cell. | |
int | Ncells |
Number of cells in the simulation. | |
int | Nvertices |
Number of vertices. | |
GPUArray< Dscalar2 > | cellPositions |
Cell positions... not used for computation, but can track, e.g., MSD of cell centers. | |
GPUArray< Dscalar2 > | vertexPositions |
Position of the vertices. | |
GPUArray< Dscalar2 > | cellVelocities |
The velocity vector of cells (only relevant if the equations of motion use it) | |
GPUArray< Dscalar > | cellMasses |
The masses of the cells. | |
GPUArray< Dscalar2 > | vertexVelocities |
The velocity vector of vertices (only relevant if the equations of motion use it) | |
GPUArray< Dscalar > | vertexMasses |
The masses of the vertices. | |
GPUArray< int > | vertexNeighbors |
VERTEX neighbors of every vertex. More... | |
GPUArray< int > | vertexCellNeighbors |
Cell neighbors of every vertex. More... | |
Index2D | n_idx |
A 2dIndexer for computing where in the GPUArray to look for a given cell's vertices. | |
GPUArray< int > | cellNeighborNum |
The number of CELL neighbors of each cell. For simple models this is the same as cellVertexNum, but does not have to be. | |
GPUArray< int > | cellNeighbors |
CELL neighbors of every cell. | |
GPUArray< int > | cellVertexNum |
The number of vertices defining each cell. More... | |
GPUArray< Dscalar2 > | vertexForces |
an array containing net force on each vertex | |
GPUArray< Dscalar2 > | cellForces |
an array containing net force on each cell | |
GPUArray< int > | cellType |
An array of integers labeling cell type...an easy way of determining if cells are different. More... | |
Index2D | cellTypeIndexer |
A indexer for turning a pair of cells into a 1-D index. | |
Dscalar | Energy |
The current potential energy of the system; only updated when an explicit energy calculation is called (i.e. not by default each timestep) | |
Dscalar | KineticEnergy |
The current kinetic energy of the system; only updated when an explicit calculation is called. | |
vector< int > | tagToIdx |
To write consistent files...the cell that started the simulation as index i has current index tagToIdx[i]. More... | |
vector< int > | tagToIdxVertex |
To write consistent files...the vertex that started the simulation as index i has current index tagToIdx[i]. | |
BoxPtr | Box |
the box defining the periodic domain | |
int | Timestep |
Count the number of times "performTimeStep" has been called. | |
Dscalar | deltaT |
The time stepsize of the simulation. | |
bool | forcesUpToDate |
Are the forces (and hence, the geometry) up-to-date? | |
Dscalar | currentTime |
a time variable for keeping track of the simulation variable (for databases) | |
Protected Member Functions | |
Matrix2x2 | d2Edridrj (int i, int j, neighborType neighbor, Dscalar unstress=1.0, Dscalar stress=1.0) |
Second derivative of the energy w/r/t cell positions...for getting dynMat info. More... | |
void | spatialSorting () |
sort points along a Hilbert curve for data locality More... | |
void | fullTriangulation (bool verbose=false) |
construct the global periodic triangulation point-by-point using non-CGAL methods More... | |
void | globalTriangulationCGAL (bool verbose=false) |
Globally construct the triangulation via CGAL. More... | |
void | getCircumcenterIndices (bool secondtime=false, bool verbose=false) |
build the auxiliary data structure containing the indices of the particle circumcenters from the neighbor list More... | |
void | testTriangulation () |
Test the current neighbor list to see if it is still a valid triangulation. GPU function. More... | |
void | testTriangulationCPU () |
Test the validity of the triangulation on the CPU. More... | |
void | repairTriangulation (vector< int > &fixlist) |
repair any problems with the triangulation on the CPU More... | |
void | testAndRepairTriangulation (bool verb=false) |
A workhorse function that calls the appropriate topology testing and repairing routines. More... | |
void | allDelSets () |
call getDelSets for all particles More... | |
bool | getDelSets (int i) |
void | resetLists () |
resize all neighMax-related arrays More... | |
void | resizeAndReset () |
do resize and resetting operations common to cellDivision and cellDeath More... | |
Matrix2x2 | dHdri (Dscalar2 ri, Dscalar2 rj, Dscalar2 rk) |
The derivative of a voronoi vertex position with respect to change in the first cells position. More... | |
Dscalar2 | dAidrj (int i, int j) |
Derivative of the area of cell i with respect to the position of cell j. More... | |
Dscalar2 | dPidrj (int i, int j) |
Derivative of the perimeter of cell i with respect to the position of cell j. More... | |
Matrix2x2 | d2Areadvdr (Matrix2x2 &dvpdr, Matrix2x2 &dvmdr) |
Second derivative of area w/r/t voronoi and cell position. More... | |
Matrix2x2 | d2Peridvdr (Matrix2x2 &dvdr, Matrix2x2 &dvmdr, Matrix2x2 &dvpdr, Dscalar2 vm, Dscalar2 v, Dscalar2 vp) |
Second derivative of perimeter w/r/t voronoi and cell position. More... | |
vector< Dscalar > | d2Hdridrj (Dscalar2 rj, Dscalar2 rk, int jj) |
second derivatives of voronoi vertex with respect to cell positions More... | |
void | spatiallySortVerticesAndCellActivity () |
call the Simple2DCell spatial vertex sorter, and re-index arrays of cell activity More... | |
void | spatiallySortCellsAndCellActivity () |
call the Simple2DCell spatial cell sorter, and re-index arrays of cell activity More... | |
void | initializeCellSorting () |
set the size of the cell-sorting structures, initialize lists simply More... | |
void | initializeVertexSorting () |
set the size of the vertex-sorting structures, initialize lists simply More... | |
void | reIndexCellArray (GPUArray< int > &array) |
Re-index cell arrays after a spatial sorting has occured. More... | |
void | reIndexCellArray (GPUArray< Dscalar > &array) |
why use templates when you can type more? More... | |
void | reIndexCellArray (GPUArray< Dscalar2 > &array) |
why use templates when you can type more? More... | |
void | reIndexVertexArray (GPUArray< int > &array) |
Re-index vertex after a spatial sorting has occured. | |
void | reIndexVertexArray (GPUArray< Dscalar > &array) |
why use templates when you can type more? | |
void | reIndexVertexArray (GPUArray< Dscalar2 > &array) |
why use templates when you can type more? More... | |
void | spatiallySortCells () |
Perform a spatial sorting of the cells to try to maintain data locality. More... | |
void | spatiallySortVertices () |
Perform a spatial sorting of the vertices to try to maintain data locality. More... | |
Protected Attributes | |
cellListGPU | celllist |
The associated cell list structure. | |
Dscalar | cellsize |
The size of the cell list's underlying grid. | |
int | neighMax |
An upper bound for the maximum number of neighbors that any cell has. | |
GPUArray< int2 > | NeighIdxs |
int | NeighIdxNum |
A utility integer to help with NeighIdxs. | |
GPUArray< int3 > | circumcenters |
A data structure that holds the indices of particles forming the circumcircles of the Delaunay Triangulation. | |
int | NumCircumCenters |
The number of circumcircles...for a periodic system, this should never change. This provides one check that local updates to the triangulation are globally consistent. | |
GPUArray< int > | anyCircumcenterTestFailed |
A flag that can be accessed by child classes... serves as notification that any change in the network topology has occured. | |
int | completeRetriangulationPerformed |
A flag that notifies that a global re-triangulation has been performed. | |
bool | neighMaxChange |
A flag that notifies that the maximum number of neighbors may have changed, necessitating resizing of some data arrays. | |
GPUArray< int > | repair |
A a vector of zeros (everything is fine) and ones (that index needs to be repaired) | |
vector< int > | NeedsFixing |
A smaller vector that, after testing the triangulation, contains the particle indices that need their local topology to be updated. | |
bool | globalOnly |
When true, the CPU branch will execute global retriangulations through CGAL on every time step. More... | |
int | timestep |
Count the number of times that testAndRepair has been called, separately from the derived class' time. | |
bool | particleExclusions |
A flag that notifies the existence of any particle exclusions (for which the net force is set to zero by fictitious external forces) | |
GPUArray< int2 > | delSets |
delSet.data[n_idx(nn,i)] are the previous and next consecutive delaunay neighbors More... | |
GPUArray< int > | delOther |
delOther.data[n_idx(nn,i)] contains the index of the "other" delaunay neighbor. More... | |
GPUArray< Dscalar2 > | forceSets |
In GPU mode, interactions are computed "per voronoi vertex"...forceSets are summed up to get total force on a particle. | |
bool | GPUcompute |
Compute aspects of the model on the GPU. | |
bool | Reproducible |
A flag that determines whether the GPU RNG is the same every time. | |
noiseSource | noise |
A source of noise for random cell initialization. | |
Dscalar | KA |
the area modulus | |
Dscalar | KP |
The perimeter modulus. | |
GPUArray< Dscalar2 > | Moduli |
The area and perimeter moduli of each cell. CURRENTLY NOT SUPPORTED, BUT EASY TO IMPLEMENT. | |
GPUArray< Dscalar2 > | AreaPeri |
The current area and perimeter of each cell. | |
GPUArray< Dscalar2 > | AreaPeriPreferences |
The area and perimeter preferences of each cell. | |
GPUArray< int > | cellVertices |
A structure that indexes the vertices defining each cell. More... | |
int | vertexMax |
An upper bound for the maximum number of neighbors that any cell has. | |
GPUArray< Dscalar2 > | voroCur |
3*Nvertices length array of the position of vertices around cells More... | |
GPUArray< Dscalar4 > | voroLastNext |
vector< int > | itt |
A map between cell index and the spatially sorted version. More... | |
vector< int > | tti |
A temporary structure that inverts itt. | |
vector< int > | idxToTag |
A temporary structure that inverse tagToIdx. | |
vector< int > | ittVertex |
A map between vertex index and the spatially sorted version. | |
vector< int > | ttiVertex |
A temporary structure that inverts itt. | |
vector< int > | idxToTagVertex |
A temporary structure that inverse tagToIdx. | |
GPUArray< Dscalar2 > | displacements |
An array of displacements used only for the equations of motion. | |
Friends | |
class | SPVDatabaseNetCDF |
Implement a 2D Voronoi model, with and without some extra bells and whistles, using kernels in SPV Kernels.
A child class of voronoiModelBase, this implements a Voronoi model in 2D. This involves mostly calculating the forces in the Voronoi model. Optimizing these procedures for hybrid CPU/GPU-based computation involves declaring and maintaining several related auxiliary data structures that capture different features of the local topology and local geoemetry for each cell.
VoronoiQuadraticEnergy::VoronoiQuadraticEnergy | ( | int | n, |
bool | reprod = false |
||
) |
initialize with random positions in a square box
n | number of cells to initialize |
reprod | should the simulation be reproducible (i.e. call a RNG with a fixed seed) |
VoronoiQuadraticEnergy::VoronoiQuadraticEnergy | ( | int | n, |
Dscalar | A0, | ||
Dscalar | P0, | ||
bool | reprod = false |
||
) |
initialize with random positions and set all cells to have uniform target A_0 and P_0 parameters
n | number of cells to initialize |
A0 | set uniform preferred area for all cells |
P0 | set uniform preferred perimeter for all cells |
reprod | should the simulation be reproducible (i.e. call a RNG with a fixed seed) |
void VoronoiQuadraticEnergy::initializeVoronoiQuadraticEnergy | ( | int | n | ) |
Initialize voronoiQuadraticEnergy and call the initializer chain.
n | Number of cells to initialized |
Referenced by VoronoiQuadraticEnergy().
|
virtual |
compute the geometry and get the forces
goes through the process of computing the forces on either the CPU or GPU, either with or without exclusions, as determined by the flags. Assumes the geometry has NOT yet been computed.
Reimplemented from Simple2DCell.
Reimplemented in VoronoiQuadraticEnergyWithTension.
References ComputeForceSetsGPU(), Simple2DCell::computeGeometry(), computeVoronoiForceCPU(), Simple2DCell::forcesUpToDate, Simple2DCell::GPUcompute, Simple2DCell::Ncells, and SumForcesGPU().
Referenced by VoronoiQuadraticEnergy().
|
virtual |
compute the quadratic energy functional
Returns the quadratic energy functional: E = {cells} K_A(A_i-A_i,0)^2 + K_P(P_i-P_i,0)^2
Reimplemented from Simple2DCell.
Reimplemented in VoronoiQuadraticEnergyWithTension.
|
virtual |
Compute force sets on the GPU.
Reimplemented in VoronoiQuadraticEnergyWithTension.
References computeVoronoiForceSetsGPU().
Referenced by computeForces().
void VoronoiQuadraticEnergy::SumForcesGPU | ( | ) |
Add up the force sets to get the net force per particle on the GPU.
References voronoiModelBase::particleExclusions, sumForceSets(), and sumForceSetsWithExclusions().
Referenced by VoronoiQuadraticEnergyWithTension::computeForces(), and computeForces().
|
virtual |
Compute the net force on particle i on the CPU.
i | The particle index for which to compute the net force, assuming addition tension terms between unlike particles |
Referenced by VoronoiQuadraticEnergyWithTension::computeForces(), and computeForces().
|
virtual |
call gpu_force_sets kernel caller
Calculate the contributions to the net force on particle "i" from each of particle i's voronoi vertices
References Simple2DCell::AreaPeri, Simple2DCell::AreaPeriPreferences, Simple2DCell::Box, Simple2DCell::cellPositions, ArrayHandle< T >::data, voronoiModelBase::delOther, voronoiModelBase::delSets, access_location::device, voronoiModelBase::forceSets, gpu_force_sets(), Simple2DCell::KA, Simple2DCell::KP, Simple2DCell::n_idx, voronoiModelBase::NeighIdxNum, voronoiModelBase::NeighIdxs, access_mode::overwrite, access_mode::read, Simple2DCell::voroCur, and Simple2DCell::voroLastNext.
Referenced by VoronoiQuadraticEnergyWithTension::ComputeForceSetsGPU(), and ComputeForceSetsGPU().
void VoronoiQuadraticEnergy::sumForceSets | ( | ) |
call gpu_sum_force_sets kernel caller
References Simple2DCell::cellForces, Simple2DCell::cellNeighborNum, ArrayHandle< T >::data, access_location::device, voronoiModelBase::forceSets, gpu_sum_force_sets(), Simple2DCell::n_idx, Simple2DCell::Ncells, access_mode::overwrite, and access_mode::read.
Referenced by SumForcesGPU().
void VoronoiQuadraticEnergy::sumForceSetsWithExclusions | ( | ) |
call gpu_sum_force_sets_with_exclusions kernel caller
References Simple2DCell::cellForces, Simple2DCell::cellNeighborNum, ArrayHandle< T >::data, access_location::device, voronoiModelBase::exclusions, voronoiModelBase::external_forces, voronoiModelBase::forceSets, gpu_sum_force_sets_with_exclusions(), Simple2DCell::n_idx, Simple2DCell::Ncells, access_mode::overwrite, and access_mode::read.
Referenced by SumForcesGPU().
void VoronoiQuadraticEnergy::reportCellInfo | ( | ) |
Report various cell infor for testing and debugging.
a utility function...output some information assuming the system is uniform
References Simple2DActiveCell::Dr, Simple2DCell::Ncells, and Simple2DActiveCell::v0.
|
virtual |
Save tuples for half of the dynamical matrix.
rcs | a vector of (row,col) locations |
vals | a vector of the corresponding value of the dynamical matrix |
Reimplemented from Simple2DModel.
|
virtual |
calculate the current global off-diagonal stress
This function calculates {xy} = 1/Area_{total}*(dE/d), the normalized change in energy when deforming the box with a strain tensor given by 0 0 1 . Notably, this is done by taking analytic derivatives, not by doing a finite-difference computed by actually deforming the box a bit and recomputing the geometry.
|
protected |
Second derivative of the energy w/r/t cell positions...for getting dynMat info.
i | The index of cell i |
j | The index of cell j |
|
inherited |
A default initialization scheme.
a function that takes care of the initialization of the class.
n | the number of cells to initialize |
|
inlinevirtualinherited |
Enforce CPU-only operation.
global | defaults to true. When global is set to true, the CPU branch will try the local repair scheme. This is generally slower, but if working in a regime where things change very infrequently, it may be faster. |
Implements Simple2DCell.
References voronoiModelBase::enforceTopology(), voronoiModelBase::globalOnly, Simple2DCell::GPUcompute, voronoiModelBase::readTriangulation(), voronoiModelBase::setExclusions(), and voronoiModelBase::writeTriangulation().
|
virtualinherited |
update/enforce the topology
goes through the process of testing and repairing the topology on either the CPU or GPU
Reimplemented from Simple2DCell.
References voronoiModelBase::allDelSets(), voronoiModelBase::anyCircumcenterTestFailed, voronoiModelBase::completeRetriangulationPerformed, ArrayHandle< T >::data, voronoiModelBase::delOther, voronoiModelBase::delSets, access_location::device, voronoiModelBase::getDelSets(), voronoiModelBase::globalTriangulationCGAL(), Simple2DCell::GPUcompute, access_location::host, voronoiModelBase::NeedsFixing, voronoiModelBase::NeighIdxs, voronoiModelBase::neighMaxChange, access_mode::read, voronoiModelBase::resetLists(), and voronoiModelBase::testAndRepairTriangulation().
Referenced by voronoiModelBase::setCPU().
|
inherited |
Declare which particles are to be excluded (exes[i]!=0)
exes | a list of per-particle indications of whether a particle should be excluded (exes[i] !=0) or not/ |
References ArrayHandle< T >::data, voronoiModelBase::exclusions, voronoiModelBase::external_forces, access_location::host, Simple2DActiveCell::Motility, Simple2DCell::Ncells, access_mode::overwrite, voronoiModelBase::particleExclusions, access_mode::readwrite, and GPUArray< T >::resize().
Referenced by voronoiModelBase::setCPU().
|
virtualinherited |
moveDegrees of Freedom calls either the move points or move points CPU routines
Displace cells on either the GPU or CPU, according to the flag
displacements | a vector of Dscalar2 specifying how much to move every cell |
scale | a scalar that multiples the value of every index of displacements before things are moved |
Reimplemented from Simple2DCell.
Referenced by voronoiModelBase::getNumberOfDegreesOfFreedom().
|
virtualinherited |
Compute cell geometry on the CPU.
Reimplemented from Simple2DCell.
References Simple2DCell::AreaPeri, Simple2DCell::cellNeighborNum, Simple2DCell::cellNeighbors, Simple2DCell::cellPositions, ArrayHandle< T >::data, access_location::host, Simple2DCell::n_idx, Simple2DCell::Ncells, access_mode::overwrite, access_mode::read, access_mode::readwrite, Simple2DCell::voroCur, and Simple2DCell::voroLastNext.
Referenced by voronoiModelBase::returnForces().
|
virtualinherited |
call gpu_compute_geometry kernel caller
Reimplemented from Simple2DCell.
References Simple2DCell::AreaPeri, Simple2DCell::Box, Simple2DCell::cellNeighborNum, Simple2DCell::cellNeighbors, Simple2DCell::cellPositions, ArrayHandle< T >::data, access_location::device, gpu_compute_voronoi_geometry(), Simple2DCell::n_idx, Simple2DCell::Ncells, access_mode::overwrite, access_mode::read, access_mode::readwrite, Simple2DCell::voroCur, and Simple2DCell::voroLastNext.
Referenced by voronoiModelBase::returnForces().
|
virtualinherited |
allow for cell division, according to a vector of model-dependent parameters
Trigger a cell division event, which involves some laborious re-indexing of various data structures. This version uses a heavy-handed approach, hitting the cell positions with a full, global retriangulation (rather than just updating the targeted cell and its neighbor with the local topology repair routines). If a simulation is needed where the cell division rate is very rapid, this should be improved. The idea of the division is that a targeted cell will divide normal to an axis specified by the angle, theta, passed to the function. The final state cell positions are placed along the axis at a distance away from the initial cell position set by a multiplicative factor (<1) of the in-routine determined maximum distance in the cell along that axis. parameters[0] = the index of the cell to undergo a division event dParams[0] = an angle, theta dParams[1] = a fraction of maximum separation of the new cell positions along the axis of the cell specified by theta
This function is meant to be called before the start of a new timestep. It should be immediately followed by a computeGeometry call.
Reimplemented from Simple2DActiveCell.
Referenced by voronoiModelBase::returnForces().
|
virtualinherited |
Kill the indexed cell by simply removing it from the simulation.
Trigger a cell death event. In the Voronoi model this simply removes the targeted cell and instantaneously computes the new tesselation. Very violent, if the cell isn't already small.
Reimplemented from Simple2DActiveCell.
References Simple2DActiveCell::cellDeath(), and voronoiModelBase::resizeAndReset().
Referenced by voronoiModelBase::returnForces().
|
inherited |
move particles on the GPU
Displace cells on the GPU
displacements | a vector of Dscalar2 specifying how much to move every cell |
Referenced by voronoiModelBase::returnForces().
|
inherited |
move particles on the CPU
Displace cells on the CPU
displacements | a vector of Dscalar2 specifying how much to move every cell |
|
inherited |
Transfer particle data from the GPU to the CPU for use by delLoc.
The GPU moves the location of points in the GPU memory... this gets a local copy that can be used by the DelaunayLoc class
References Simple2DCell::cellPositions, voronoiModelBase::cellsize, voronoiModelBase::delLoc, access_location::host, DelaunayLoc::initialize(), Simple2DCell::Ncells, access_mode::read, and DelaunayLoc::setPoints().
Referenced by voronoiModelBase::fullTriangulation(), voronoiModelBase::resizeAndReset(), and voronoiModelBase::testTriangulationCPU().
|
inherited |
Update the cell list structure after particles have moved.
References voronoiModelBase::celllist, Simple2DCell::cellPositions, cellListGPU::compute(), cellListGPU::computeGPU(), ArrayHandle< T >::data, Simple2DCell::GPUcompute, access_location::host, Simple2DCell::Ncells, access_mode::read, and cellListGPU::setParticles().
Referenced by voronoiModelBase::testTriangulation().
|
inherited |
update the NieghIdxs data structure
References Simple2DCell::cellNeighborNum, access_location::host, voronoiModelBase::NeighIdxs, access_mode::overwrite, and access_mode::read.
Referenced by voronoiModelBase::allDelSets(), and voronoiModelBase::fullTriangulation().
|
protectedvirtualinherited |
sort points along a Hilbert curve for data locality
When sortPeriod < 0, this routine does not get called
Reimplemented from Simple2DModel.
References voronoiModelBase::allDelSets(), voronoiModelBase::exclusions, voronoiModelBase::globalTriangulationCGAL(), Simple2DCell::reIndexCellArray(), voronoiModelBase::resetLists(), and Simple2DActiveCell::spatiallySortCellsAndCellActivity().
|
protectedinherited |
construct the global periodic triangulation point-by-point using non-CGAL methods
The DelaunayLoc and DelaunayNP classes are invoked to performed to determine the Delaunay triangulation of the entire periodic domain.
References Simple2DCell::cellNeighborNum, Simple2DCell::cellNeighbors, voronoiModelBase::completeRetriangulationPerformed, ArrayHandle< T >::data, voronoiModelBase::delLoc, voronoiModelBase::getCircumcenterIndices(), DelaunayLoc::getNeighbors(), voronoiModelBase::GlobalFixes, access_location::host, Simple2DCell::n_idx, Simple2DCell::Ncells, voronoiModelBase::neighMax, voronoiModelBase::neighMaxChange, access_mode::overwrite, voronoiModelBase::repair, voronoiModelBase::resetDelLocPoints(), GPUArray< T >::resize(), voronoiModelBase::updateNeighIdxs(), and voronoiModelBase::writeTriangulation().
|
protectedinherited |
Globally construct the triangulation via CGAL.
This function calls the DelaunayCGAL class to determine the Delaunay triangulation of the entire square periodic domain this method is, obviously, better than the version written by DMS, so should be the default option. In addition to performing a triangulation, the function also automatically calls updateNeighIdxs and getCircumcenterIndices/
References Simple2DCell::cellPositions, voronoiModelBase::completeRetriangulationPerformed, voronoiModelBase::GlobalFixes, access_location::host, Simple2DCell::Ncells, and access_mode::read.
Referenced by voronoiModelBase::enforceTopology(), voronoiModelBase::getCircumcenterIndices(), voronoiModelBase::resizeAndReset(), voronoiModelBase::spatialSorting(), voronoiModelBase::testAndRepairTriangulation(), and voronoiModelBase::testTriangulationCPU().
|
protectedinherited |
build the auxiliary data structure containing the indices of the particle circumcenters from the neighbor list
Converts the neighbor list data structure into a list of the three particle indices defining all of the circumcenters in the triangulation. Keeping this version of the topology on the GPU allows for fast testing of what points need to be retriangulated.
References Simple2DCell::cellNeighborNum, Simple2DCell::cellNeighbors, voronoiModelBase::circumcenters, ArrayHandle< T >::data, voronoiModelBase::globalTriangulationCGAL(), access_location::host, Simple2DCell::n_idx, Simple2DCell::Ncells, voronoiModelBase::NumCircumCenters, access_mode::overwrite, access_mode::read, voronoiModelBase::timestep, and voronoiModelBase::writeTriangulation().
Referenced by voronoiModelBase::fullTriangulation().
|
protectedinherited |
Test the current neighbor list to see if it is still a valid triangulation. GPU function.
Call the GPU to test each circumcenter to see if it is still empty (i.e., how much of the triangulation from the last time step is still valid?). Note that because gpu_test_circumcenters always* copies at least a single integer back and forth (to answer the question "did any circumcircle come back non-empty?" for the cpu) this function is always an implicit cuda synchronization event. At least until non-default streams are added to the code.
References voronoiModelBase::anyCircumcenterTestFailed, Simple2DCell::Box, cellListGPU::cell_indexer, cellListGPU::cell_list_indexer, cellListGPU::cell_sizes, voronoiModelBase::celllist, Simple2DCell::cellPositions, voronoiModelBase::circumcenters, ArrayHandle< T >::data, access_location::device, cellListGPU::getBoxsize(), cellListGPU::getXsize(), cellListGPU::getYsize(), gpu_test_circumcenters(), access_location::host, cellListGPU::idxs, Simple2DCell::Ncells, voronoiModelBase::NumCircumCenters, access_mode::overwrite, access_mode::read, access_mode::readwrite, voronoiModelBase::repair, and voronoiModelBase::updateCellList().
Referenced by voronoiModelBase::testAndRepairTriangulation().
|
protectedinherited |
Test the validity of the triangulation on the CPU.
perform the same check on the CPU... because of the cost of checking circumcircles and the relatively poor performance of the 1-ring calculation in DelaunayLoc, it is sometimes better to just re-triangulate the entire point set with CGAL. At the moment that is the default behavior of the cpu branch.
References voronoiModelBase::anyCircumcenterTestFailed, Simple2DCell::cellNeighborNum, Simple2DCell::cellNeighbors, ArrayHandle< T >::data, voronoiModelBase::delLoc, voronoiModelBase::globalOnly, voronoiModelBase::globalTriangulationCGAL(), access_location::host, voronoiModelBase::localTopologyUpdates, Simple2DCell::n_idx, Simple2DCell::Ncells, voronoiModelBase::neighMax, access_mode::readwrite, voronoiModelBase::repair, voronoiModelBase::resetDelLocPoints(), voronoiModelBase::skippedFrames, and DelaunayLoc::testPointTriangulation().
Referenced by voronoiModelBase::testAndRepairTriangulation().
|
protectedinherited |
repair any problems with the triangulation on the CPU
Given a list of particle indices that need to be repaired, call CGAL to figure out their neighbors and then update the relevant data structures.
References voronoiModelBase::repPerFrame.
Referenced by voronoiModelBase::testAndRepairTriangulation().
|
protectedinherited |
A workhorse function that calls the appropriate topology testing and repairing routines.
This function calls the relevant testing and repairing functions, and increments the "timestep" by one. Note that the call to testTriangulation will always synchronize the gpu (via a memcpy of the "anyCircumcenterTestFailed" variable)
References voronoiModelBase::anyCircumcenterTestFailed, Simple2DCell::cellNeighborNum, Simple2DCell::cellNeighbors, voronoiModelBase::completeRetriangulationPerformed, ArrayHandle< T >::data, voronoiModelBase::globalTriangulationCGAL(), Simple2DCell::GPUcompute, access_location::host, Simple2DCell::n_idx, Simple2DCell::Ncells, voronoiModelBase::NeedsFixing, access_mode::read, access_mode::readwrite, voronoiModelBase::repair, voronoiModelBase::repairTriangulation(), voronoiModelBase::skippedFrames, voronoiModelBase::testTriangulation(), voronoiModelBase::testTriangulationCPU(), and voronoiModelBase::timestep.
Referenced by voronoiModelBase::enforceTopology().
|
protectedinherited |
call getDelSets for all particles
Calls updateNeighIdxs and then getDelSets(i) for all cells i
References voronoiModelBase::getDelSets(), Simple2DCell::Ncells, and voronoiModelBase::updateNeighIdxs().
Referenced by voronoiModelBase::enforceTopology(), voronoiModelBase::resizeAndReset(), and voronoiModelBase::spatialSorting().
|
protectedinherited |
Maintain the delSets and delOther data structure for particle i If it returns false there was a problem and a global re-triangulation is triggered.
i | the cell in question |
References Simple2DCell::cellNeighborNum, Simple2DCell::cellNeighbors, ArrayHandle< T >::data, voronoiModelBase::delOther, voronoiModelBase::delSets, access_location::host, Simple2DCell::n_idx, access_mode::read, and access_mode::readwrite.
Referenced by voronoiModelBase::allDelSets(), and voronoiModelBase::enforceTopology().
|
protectedinherited |
resize all neighMax-related arrays
As the code is modified, all GPUArrays whose size depend on neighMax should be added to this function
References voronoiModelBase::delOther, voronoiModelBase::delSets, voronoiModelBase::forceSets, Simple2DCell::Ncells, voronoiModelBase::neighMax, GPUArray< T >::resize(), Simple2DCell::voroCur, and Simple2DCell::voroLastNext.
Referenced by voronoiModelBase::enforceTopology(), voronoiModelBase::resizeAndReset(), and voronoiModelBase::spatialSorting().
|
protectedinherited |
do resize and resetting operations common to cellDivision and cellDeath
A utility function for resizing data arrays... used by the cell division and cell death routines
References voronoiModelBase::allDelSets(), Simple2DCell::cellForces, voronoiModelBase::celllist, Simple2DCell::cellNeighborNum, voronoiModelBase::circumcenters, Simple2DCell::displacements, voronoiModelBase::exclusions, voronoiModelBase::external_forces, voronoiModelBase::globalTriangulationCGAL(), Simple2DCell::Ncells, voronoiModelBase::NeighIdxs, voronoiModelBase::repair, voronoiModelBase::resetDelLocPoints(), voronoiModelBase::resetLists(), GPUArray< T >::resize(), and cellListGPU::setNp().
Referenced by voronoiModelBase::cellDeath().
|
protectedinherited |
The derivative of a voronoi vertex position with respect to change in the first cells position.
ri | The position of cell i |
rj | The position of cell j |
rk | The position of cell k Returns the derivative of the voronoi vertex shared by cells i, j , and k with respect to changing the position of cell i the (row, column) format specifies dH_{row}/dr_{i,column} |
|
protectedinherited |
Derivative of the area of cell i with respect to the position of cell j.
i | The index of cell i |
j | The index of cell j |
|
protectedinherited |
Derivative of the perimeter of cell i with respect to the position of cell j.
i | The index of cell i |
j | The index of cell j Returns the derivative of the perimeter of cell i w/r/t the position of cell j |
Second derivative of area w/r/t voronoi and cell position.
dvpdr | derivative of v_{i+1} w/r/t a cell position |
dvmdr | derivative of v_{i-1} w/r/t a cell position |
References Matrix2x2::x11.
|
protectedinherited |
Second derivative of perimeter w/r/t voronoi and cell position.
dvdr | derivative of v_{i} w/r/t a cell position |
dvpdr | derivative of v_{i+1} w/r/t a cell position |
dvmdr | derivative of v_{i-1} w/r/t a cell position |
vm | position of v_{i-1} |
v | position of v_{i} |
vp | position of v_{i+1} |
|
protectedinherited |
second derivatives of voronoi vertex with respect to cell positions
ri | The position of cell i |
rj | The position of cell j |
rk | The position of cell k |
jj | the index EITHER 1 or 2 of the second derivative Returns an 8-component vector containing the derivatives of the voronoi vertex formed by cells i, j, and k with respect to r_i and r_{jj}... jj should be either 1 (to give d^2H/(d r_i)^2 or 2 (to give d^2H/dridrj) The vector is laid out as (H_x/r_{i,x}r_{j,x}, H_y/r_{i,x}r_{j,x} H_x/r_{i,y}r_{j,x}, H_y/r_{i,y}r_{j,x H_x/r_{i,x}r_{j,y}, H_y/r_{i,x}r_{j,y H_x/r_{i,y}r_{j,y}, H_y/r_{i,y}r_{j,y} ) NOTE: This function does not check that ri, rj, and rk actually share a voronoi vertex in the triangulation NOTE: This function assumes that rj and rk are w/r/t the position of ri, so ri = (0.,0.) |
|
inherited |
initialize class' data structures and set default values
Initialize the data structures to the size specified by n, and set default values, and call Simple2DCell's initilization routine.
|
virtualinherited |
Set uniform motility.
v0new | the new value of velocity for all cells |
drnew | the new value of the rotational diffusion of cell directors for all cells |
Implements Simple2DCell.
|
inherited |
Set non-uniform cell motilites.
v0s | the per-particle vector of what all velocities will be |
drs | the per-particle vector of what all rotational diffusions will be |
References Simple2DActiveCell::cellDirectors, Simple2DCell::cellVelocities, ArrayHandle< T >::data, access_location::host, Simple2DActiveCell::Motility, Simple2DCell::Ncells, access_mode::overwrite, access_mode::read, and GPUArray< T >::resize().
|
inherited |
Set random cell directors (for active cell models)
Assign cell directors via a simple, reproducible RNG
References Simple2DActiveCell::cellDirectorForces, Simple2DActiveCell::cellDirectors, Simple2DCell::cellVelocities, ArrayHandle< T >::data, noiseSource::getRealUniform(), access_location::host, Simple2DCell::Ncells, Simple2DCell::noise, access_mode::overwrite, noiseSource::Reproducible, Simple2DCell::Reproducible, and GPUArray< T >::resize().
|
protectedinherited |
call the Simple2DCell spatial vertex sorter, and re-index arrays of cell activity
Calls the spatial vertex sorting routine in Simple2DCell, and re-indexes the arrays for the cell RNGS, as well as the cell motility and cellDirector arrays
References Simple2DActiveCell::cellDirectors, Simple2DActiveCell::Motility, Simple2DCell::reIndexCellArray(), and Simple2DCell::spatiallySortVertices().
Referenced by vertexModelBase::spatialSorting().
|
protectedinherited |
call the Simple2DCell spatial cell sorter, and re-index arrays of cell activity
Calls the spatial vertex sorting routine in Simple2DCell, and re-indexes the arrays for the cell RNGS, as well as the cell motility and cellDirector arrays
References Simple2DActiveCell::cellDirectors, Simple2DActiveCell::Motility, Simple2DCell::reIndexCellArray(), and Simple2DCell::spatiallySortCells().
Referenced by voronoiModelBase::spatialSorting().
|
inherited |
initialize class' data structures and set default values
Initialize the data structures to the size specified by n, and set default values.
|
virtualinherited |
call either the computeGeometryCPU or GPU routines for the current model
Simply call either the CPU or GPU routine in the current or derived model
References Simple2DCell::computeGeometryCPU(), Simple2DCell::computeGeometryGPU(), and Simple2DCell::GPUcompute.
Referenced by VertexQuadraticEnergyWithTension::computeForces(), VoronoiQuadraticEnergyWithTension::computeForces(), VertexQuadraticEnergy::computeForces(), computeForces(), and Simple2DCell::computeForces().
|
inherited |
Call masses and velocities to get the total kinetic energy.
E = 0.5*m_i v_i^2
|
inherited |
Call masses and velocities to get the average kinetic contribution to the pressure tensor.
P_ab = m_i v_{ib}v_{ia}
|
inherited |
Set uniform cell area and perimeter preferences.
Generically believe that cells in 2D have a notion of a preferred area and perimeter
Referenced by Simple2DCell::enforceTopology().
|
inherited |
Set cell area and perimeter preferences according to input vector.
Set the Area and Perimeter preferences to the input vector
References Simple2DCell::AreaPeriPreferences, ArrayHandle< T >::data, access_location::host, Simple2DCell::Ncells, access_mode::overwrite, and GPUArray< T >::resize().
|
inherited |
Set random cell positions, and set the periodic box to a square with average cell area=1.
Resize the box so that every cell has, on average, area = 1, and place cells via either a simple, reproducible RNG or a non-reproducible RNG
References Simple2DCell::cellPositions, Simple2DCell::Ncells, and GPUArray< T >::resize().
|
inherited |
Set cell positions according to a user-specified vector.
Does not update any other lists – it is the user's responsibility to maintain topology, etc, when using this function.
References Simple2DCell::cellPositions, ArrayHandle< T >::data, GPUArray< T >::getNumElements(), access_location::host, Simple2DCell::Ncells, access_mode::overwrite, and GPUArray< T >::resize().
|
inherited |
Set vertex positions according to a user-specified vector.
Does not update any other lists – it is the user's responsibility to maintain topology, etc, when using this function.
References ArrayHandle< T >::data, GPUArray< T >::getNumElements(), access_location::host, Simple2DCell::Nvertices, access_mode::overwrite, GPUArray< T >::resize(), and Simple2DCell::vertexPositions.
|
inherited |
Set velocities via a temperature. The return value is the total kinetic energy.
Set the cell velocities by drawing from a Maxwell-Boltzmann distribution, and then make sure there is no net momentum. The return value is the total kinetic energy
|
inherited |
Set velocities via a temperature for the vertex degrees of freedom.
Set the vertex velocities by drawing from a Maxwell-Boltzmann distribution, and then make sure there is no net momentum. The return value is the total kinetic energy.
|
inherited |
set uniform moduli for all cells
set all cell K_A, K_P preferences to uniform values. PLEASE NOTE that as an optimization this data is not actually used at the moment, but the code could be trivially altered to use this
|
inherited |
Set all cells to the same "type".
set all cell types to i
References Simple2DCell::cellType, ArrayHandle< T >::data, access_location::host, Simple2DCell::Ncells, access_mode::overwrite, and GPUArray< T >::resize().
|
inherited |
Set cells to different "type".
types | a vector of integers that the cell types will be set to |
References Simple2DCell::cellType, ArrayHandle< T >::data, access_location::host, Simple2DCell::Ncells, access_mode::overwrite, and GPUArray< T >::resize().
|
inherited |
An uncomfortable function to allow the user to set vertex topology "by hand".
This function allows a user to set the vertex topology by hand. The user is responsible for making sure the input topology is sensible. DMS NOTE – this functionality has not been thoroughly tested
cellVertexIndices | a vector of vector of ints. Each vector of ints must correspond to the counter-clockwise ordering of vertices that make up the cell, and every vertex should appear at most three times in different cells |
References Simple2DCell::cellVertexNum, Simple2DCell::cellVertices, ArrayHandle< T >::data, access_location::host, Simple2DCell::n_idx, Simple2DCell::Ncells, Simple2DCell::Nvertices, access_mode::overwrite, GPUArray< T >::resize(), Simple2DCell::vertexCellNeighbors, Simple2DCell::vertexMax, and Simple2DCell::vertexNeighbors.
|
protectedinherited |
set the size of the cell-sorting structures, initialize lists simply
Sets the size of itt, tti, idxToTag, and tagToIdx, and sets all of them so that array[i] = i, i.e., unsorted
References Simple2DCell::idxToTag, Simple2DCell::itt, Simple2DCell::Ncells, Simple2DCell::tagToIdx, and Simple2DCell::tti.
|
protectedinherited |
set the size of the vertex-sorting structures, initialize lists simply
Sets the size of ittVertex, ttiVertex, idxToTagVertex, and tagToIdxVertex,and sets all of them so that array[i] = i, i.e., things are unsorted
References Simple2DCell::idxToTagVertex, Simple2DCell::ittVertex, Simple2DCell::Nvertices, Simple2DCell::tagToIdxVertex, and Simple2DCell::ttiVertex.
|
protectedinherited |
Re-index cell arrays after a spatial sorting has occured.
Re-indexes GPUarrays of ints
References ArrayHandle< T >::data, access_location::host, Simple2DCell::itt, Simple2DCell::Ncells, access_mode::read, and access_mode::readwrite.
Referenced by Simple2DCell::spatiallySortCells(), Simple2DActiveCell::spatiallySortCellsAndCellActivity(), Simple2DCell::spatiallySortVertices(), Simple2DActiveCell::spatiallySortVerticesAndCellActivity(), and voronoiModelBase::spatialSorting().
|
protectedinherited |
why use templates when you can type more?
Re-indexes GPUarrays of Dscalars
References ArrayHandle< T >::data, access_location::host, Simple2DCell::itt, Simple2DCell::Ncells, access_mode::read, and access_mode::readwrite.
|
protectedinherited |
why use templates when you can type more?
Always called after spatial sorting is performed, reIndexCellArray shuffles the order of an array based on the spatial sort order of the cells
References ArrayHandle< T >::data, access_location::host, Simple2DCell::itt, Simple2DCell::Ncells, access_mode::read, and access_mode::readwrite.
|
protectedinherited |
why use templates when you can type more?
Called if the vertices need to be spatially sorted, reIndexVertexArray shuffles the order of an array based on the spatial sort order of the vertices
References ArrayHandle< T >::data, access_location::host, Simple2DCell::ittVertex, Simple2DCell::Nvertices, access_mode::read, and access_mode::readwrite.
|
protectedinherited |
Perform a spatial sorting of the cells to try to maintain data locality.
take the current location of the cells and sort them according the their order along a 2D Hilbert curve
References Simple2DCell::AreaPeri, Simple2DCell::AreaPeriPreferences, Simple2DCell::Box, Simple2DCell::cellMasses, Simple2DCell::cellPositions, Simple2DCell::cellType, Simple2DCell::cellVelocities, ArrayHandle< T >::data, HilbertSorter::getIdx(), access_location::host, Simple2DCell::idxToTag, Simple2DCell::itt, Simple2DCell::Moduli, Simple2DCell::Ncells, access_mode::readwrite, Simple2DCell::reIndexCellArray(), Simple2DCell::tagToIdx, and Simple2DCell::tti.
Referenced by Simple2DActiveCell::spatiallySortCellsAndCellActivity().
|
protectedinherited |
Perform a spatial sorting of the vertices to try to maintain data locality.
take the current location of the vertices and sort them according the their order along a 2D Hilbert curve. This routine first sorts the vertices, and then uses the vertex sorting to derive a sorting of the cells
References Simple2DCell::AreaPeriPreferences, Simple2DCell::Box, Simple2DCell::cellPositions, Simple2DCell::cellType, Simple2DCell::cellVertexNum, Simple2DCell::cellVertices, ArrayHandle< T >::data, HilbertSorter::getIdx(), access_location::host, Simple2DCell::idxToTag, Simple2DCell::idxToTagVertex, Simple2DCell::itt, Simple2DCell::ittVertex, Simple2DCell::Moduli, Simple2DCell::n_idx, Simple2DCell::Ncells, Simple2DCell::Nvertices, access_mode::read, access_mode::readwrite, Simple2DCell::reIndexCellArray(), Simple2DCell::reIndexVertexArray(), Simple2DCell::tagToIdx, Simple2DCell::tagToIdxVertex, Simple2DCell::tti, Simple2DCell::ttiVertex, Simple2DCell::vertexCellNeighbors, Simple2DCell::vertexMasses, Simple2DCell::vertexNeighbors, Simple2DCell::vertexPositions, and Simple2DCell::vertexVelocities.
Referenced by Simple2DActiveCell::spatiallySortVerticesAndCellActivity().
|
virtualinherited |
Report the current average force on each cell.
a utility/testing function...output the currently computed mean net force to screen.
verbose | if true also print out the force on each cell |
Implements Simple2DModel.
References Simple2DCell::cellForces, Simple2DCell::cellPositions, access_location::host, and access_mode::read.
|
virtualinherited |
Report the average value of p/sqrt(A) for the cells in the system.
Returns the mean value of the shape parameter:
Implements Simple2DModel.
|
inherited |
Report the variance of p/sqrt(A) for the cells in the system.
Returns the variance of the shape parameter:
|
inherited |
Report the variance of A and P for the cells in the system.
Returns the variance of the A and P for the system:
|
inherited |
Report the mean value of the perimeter.
Returns the mean value of the perimeter
|
protectedinherited |
An array that holds (particle, neighbor_number) info to avoid intra-warp divergence in GPU -based force calculations that might be used by child classes
Referenced by computeVoronoiForceSetsGPU(), VoronoiQuadraticEnergyWithTension::computeVoronoiSimpleTensionForceSetsGPU(), VoronoiQuadraticEnergyWithTension::computeVoronoiTensionForceSetsGPU(), voronoiModelBase::enforceTopology(), voronoiModelBase::resizeAndReset(), and voronoiModelBase::updateNeighIdxs().
|
protectedinherited |
When true, the CPU branch will execute global retriangulations through CGAL on every time step.
When running on the CPU, should only global retriangulations be performed, or should local test-and-updates still be performed? Depending on parameters simulated, performance here can be quite difference, since the circumcircle test itself is CPU expensive
Referenced by voronoiModelBase::setCPU(), and voronoiModelBase::testTriangulationCPU().
|
protectedinherited |
delSet.data[n_idx(nn,i)] are the previous and next consecutive delaunay neighbors
These are orientationally ordered, of point i (for use in computing forces on GPU)
Referenced by computeVoronoiForceSetsGPU(), VoronoiQuadraticEnergyWithTension::computeVoronoiSimpleTensionForceSetsGPU(), VoronoiQuadraticEnergyWithTension::computeVoronoiTensionForceSetsGPU(), voronoiModelBase::enforceTopology(), voronoiModelBase::getDelSets(), and voronoiModelBase::resetLists().
|
protectedinherited |
delOther.data[n_idx(nn,i)] contains the index of the "other" delaunay neighbor.
i.e., the mutual neighbor of delSet.data[n_idx(nn,i)].y and delSet.data[n_idx(nn,i)].z that isn't point i
Referenced by computeVoronoiForceSetsGPU(), VoronoiQuadraticEnergyWithTension::computeVoronoiSimpleTensionForceSetsGPU(), VoronoiQuadraticEnergyWithTension::computeVoronoiTensionForceSetsGPU(), voronoiModelBase::enforceTopology(), voronoiModelBase::getDelSets(), and voronoiModelBase::resetLists().
|
inherited |
VERTEX neighbors of every vertex.
in general, we have: vertexNeighbors[3*i], vertexNeighbors[3*i+1], and vertexNeighbors[3*i+2] contain the indices of the three vertices that are connected to vertex i
Referenced by vertexModelBase::cellDeath(), vertexModelBase::flipEdgesGPU(), Simple2DCell::setVertexTopologyFromCells(), Simple2DCell::spatiallySortVertices(), vertexModelBase::testAndPerformT1TransitionsCPU(), and vertexModelBase::testEdgesForT1GPU().
|
inherited |
Cell neighbors of every vertex.
in general, we have: vertexCellNeighbors[3*i], vertexCellNeighbors[3*i+1], and vertexCellNeighbors[3*i+2] contain the indices of the three cells are neighbors of vertex i
Referenced by vertexModelBase::cellDeath(), VertexQuadraticEnergy::computeForcesCPU(), VertexQuadraticEnergy::computeForcesGPU(), vertexModelBase::computeGeometryCPU(), vertexModelBase::computeGeometryGPU(), VertexQuadraticEnergyWithTension::computeVertexTensionForcesCPU(), vertexModelBase::flipEdgesGPU(), vertexModelBase::getCellVertexSetForT1(), Simple2DCell::setVertexTopologyFromCells(), Simple2DCell::spatiallySortVertices(), vertexModelBase::testAndPerformT1TransitionsCPU(), and vertexModelBase::testEdgesForT1GPU().
|
inherited |
The number of vertices defining each cell.
cellVertexNum[c] is an integer storing the number of vertices that make up the boundary of cell c.
Referenced by vertexModelBase::cellDeath(), vertexModelBase::computeGeometryCPU(), vertexModelBase::computeGeometryGPU(), VertexQuadraticEnergyWithTension::computeVertexTensionForcesCPU(), vertexModelBase::flipEdgesGPU(), vertexModelBase::getCellCentroidsCPU(), vertexModelBase::getCellPositionsCPU(), vertexModelBase::getCellPositionsGPU(), vertexModelBase::getCellVertexSetForT1(), vertexModelBase::growCellVerticesList(), vertexModelBase::reportNeighborsCell(), Simple2DCell::setVertexTopologyFromCells(), Simple2DCell::spatiallySortVertices(), vertexModelBase::testAndPerformT1TransitionsCPU(), and vertexModelBase::testEdgesForT1GPU().
|
inherited |
An array of integers labeling cell type...an easy way of determining if cells are different.
Please note that "type" is not meaningful unless it is used by child classes. That is, things like area/perimeter preferences, or motility, or whatever are neither set nor accessed by cell type, but rather by cell index! Thus, this is just an additional data structure that can be useful. For instance, the VoronoiTension2D classes uses the integers of cellType to determine when to apply an additional line tension between cells.
Referenced by Simple2DCell::cellDeath(), Simple2DCell::cellDivision(), VertexQuadraticEnergyWithTension::computeVertexTensionForcesCPU(), VoronoiQuadraticEnergyWithTension::computeVoronoiSimpleTensionForceSetsGPU(), VoronoiQuadraticEnergyWithTension::computeVoronoiTensionForceSetsGPU(), Simple2DCell::setCellType(), Simple2DCell::setCellTypeUniform(), Simple2DCell::spatiallySortCells(), and Simple2DCell::spatiallySortVertices().
|
inherited |
To write consistent files...the cell that started the simulation as index i has current index tagToIdx[i].
The Hilbert sorting stuff makes keeping track of particles, and re-indexing things when particle number changes, a pain. Here's a description of the four relevant data structures. tagToIdx[i] = a. At the beginning of a simulation, a particle had index "i", meaning its current state was found in position "i" of all various data vectors and arrays. That same particle's data is now in position "a" of those data structures. Short version: "Where do I look to find info for what I orinally called partice i?" idxToTag[a] = i. That is, idxToTag just helps invert the tagToIdx list. idxToTag[tagToIdx[i]]=i The above two structures (and the vertex versions of them) tell you how to go back and forth between the current state of the system and the initial state of the system. What about going back and forth between the current sorted state and the previous sorted state? The "itt" and "tti" vectors give this information. The itt and tti vectors are completely overwritten each time a spatial sorting is called. By the way, I apologize if the nomenclature of "index" vs. "tag" is the opposite of what you, the reader of these code comments, might expect.
Referenced by Simple2DCell::cellDeath(), Simple2DCell::cellDivision(), Simple2DCell::initializeCellSorting(), Simple2DCell::spatiallySortCells(), and Simple2DCell::spatiallySortVertices().
|
protectedinherited |
A structure that indexes the vertices defining each cell.
cellVertices is a large, 1D array containing the vertices associated with each cell. It must be accessed with the help of the Index2D structure n_idx. the index of the kth vertex of cell c (where the ordering is counter-clockwise starting with a random vertex) is given by cellVertices[n_idx(k,c)];
Referenced by vertexModelBase::cellDeath(), vertexModelBase::computeGeometryCPU(), vertexModelBase::computeGeometryGPU(), VertexQuadraticEnergyWithTension::computeVertexTensionForcesCPU(), vertexModelBase::flipEdgesGPU(), vertexModelBase::getCellCentroidsCPU(), vertexModelBase::getCellPositionsCPU(), vertexModelBase::getCellPositionsGPU(), vertexModelBase::getCellVertexSetForT1(), vertexModelBase::growCellVerticesList(), vertexModelBase::reportNeighborsCell(), Simple2DCell::setVertexTopologyFromCells(), Simple2DCell::spatiallySortVertices(), vertexModelBase::testAndPerformT1TransitionsCPU(), and vertexModelBase::testEdgesForT1GPU().
|
protectedinherited |
3*Nvertices length array of the position of vertices around cells
For both vertex and Voronoi models, it may help to save the relative position of the vertices around a cell, either for easy force computation or in the geometry routine, etc. voroCur.data[n_idx(nn,i)] gives the nth vertex, in CCW order, of cell i
Referenced by VertexQuadraticEnergy::computeForcesCPU(), VertexQuadraticEnergy::computeForcesGPU(), vertexModelBase::computeGeometryCPU(), voronoiModelBase::computeGeometryCPU(), vertexModelBase::computeGeometryGPU(), voronoiModelBase::computeGeometryGPU(), VertexQuadraticEnergyWithTension::computeVertexTensionForcesCPU(), computeVoronoiForceSetsGPU(), VoronoiQuadraticEnergyWithTension::computeVoronoiSimpleTensionForceSetsGPU(), VoronoiQuadraticEnergyWithTension::computeVoronoiTensionForceSetsGPU(), and voronoiModelBase::resetLists().
|
protectedinherited |
3*Nvertices length array of the position of the last and next vertices along the cell Similarly, voroLastNext.data[n_idx(nn,i)] gives the previous and next vertex of the same
Referenced by VertexQuadraticEnergy::computeForcesCPU(), VertexQuadraticEnergy::computeForcesGPU(), vertexModelBase::computeGeometryCPU(), voronoiModelBase::computeGeometryCPU(), vertexModelBase::computeGeometryGPU(), voronoiModelBase::computeGeometryGPU(), VertexQuadraticEnergyWithTension::computeVertexTensionForcesCPU(), computeVoronoiForceSetsGPU(), VoronoiQuadraticEnergyWithTension::computeVoronoiSimpleTensionForceSetsGPU(), VoronoiQuadraticEnergyWithTension::computeVoronoiTensionForceSetsGPU(), and voronoiModelBase::resetLists().
|
protectedinherited |
A map between cell index and the spatially sorted version.
sortedArray[i] = unsortedArray[itt[i]] after a hilbert sort
Referenced by Simple2DCell::cellDeath(), Simple2DCell::cellDivision(), Simple2DCell::initializeCellSorting(), Simple2DCell::reIndexCellArray(), Simple2DCell::returnItt(), Simple2DCell::spatiallySortCells(), and Simple2DCell::spatiallySortVertices().