Path Integral Quantum Monte Carlo
lookuptable.h
Go to the documentation of this file.
1 
9 #ifndef LOOKUPTABLE_H
10 #define LOOKUPTABLE_H
11 
12 #include "common.h"
13 #include "constants.h"
14 #include "container.h"
15 
16 class Path;
17 
18 // ========================================================================
19 // LookupTable Class
20 // ========================================================================
29 class LookupTable {
30 
31  public:
32  LookupTable(const Container *, const int, const int);
33  ~LookupTable();
34 
35  const Container *boxPtr;
36 
38  int numNN;
39 
40  Array <beadLocator,1> beadList;
41  Array <beadLocator,1> fullBeadList;
42  Array <dVec,1> beadSep;
43  int numBeads;
45 
47  iVec getNumNNGrid() {return numNNGrid;}
49  int getTotNumGridBoxes() {return totNumGridBoxes;}
50 
51  /* Update the NN table interaction list */
52  void updateInteractionList(const Path &, const beadLocator &);
53  void updateFullInteractionList(const beadLocator &, const int);
54  void updateFullInteractionList(const int, const int);
55  void updateGrid(const Path &);
56  void updateGrid(const Array <dVec,1>&);
57 
58  void printGrid();
59 
60  /* Remove a bead from the grid */
61  void delBead(const beadLocator &);
62  void addBead(const beadLocator &, const dVec &);
63  void updateBead(const beadLocator &, const dVec &);
64 
65  /* Returns the grid index where the suplied position resides */
66  inline iVec gridIndex(const dVec &);
67 
68  /* Returns the grid number where the suplied position resides */
69  inline int gridNumber(const dVec &);
70 
71  /* Converts a grid number to a grid index */
72  inline iVec gridIndex(int);
73 
74  /* Converts a grid index to a grid number */
75  inline int gridNumber(const iVec &);
76 
78  void resizeList(int _numParticles) {
79  beadList.resize(_numParticles);
80  fullBeadList.resize(_numParticles);
81  beadSep.resize(_numParticles);
82  grid.resizeAndPreserve(numLookupTimeSlices,_numParticles);
83  beadLabel.resizeAndPreserve(numLookupTimeSlices,_numParticles);
84  }
85 
86  /* Determine if two positions are in neighboring grid boxes */
87  bool gridNeighbors(const beadLocator&, const beadLocator&);
88 
90  bool gridShare(const beadLocator &bead1, const beadLocator &bead2) {
91  return all(grid(bead1)==grid(bead2));
92  }
93 
94  private:
95  int numLookupTimeSlices; // The number of time slices that
96  int totNumGridBoxes; // The total number of nn grid boxes
97 
98  iVec numNNGrid; // The number of nn grid boxes in each direction
99  iVec gIndex; // A commonly used grid box index vector
100 
101  double rc2; // A local copy of the potential cutoff squared
102 
103  TinyVector<int,NDIM+1> nnIndex; // Comonly used nn index vector
104  TinyVector<int,NDIM+1> nI; // Used for indexing numLabels
105  TinyVector<int,NDIM+2> hI; // Used for indexing hash
106 
107  Array <iVec,NDIM+1> gridNN; // The nearest neighbors of each grid box
108  Array <iVec,NDIM+1> gridNNReduced; // The nn reduced to contain no back links
109 
110  dVec sizeNNGrid; // The size of the nn grid boxes in each direction
111 
112  Array <int,NDIM+2> hash; // The main worldline lookup array
113  Array <iVec,2> grid; // The grid index of a bead
114  Array <int,2> beadLabel; // The label of a bead in a cell
115  Array <int,NDIM+1> numLabels; // The number of beads in a cell
116 
117  TinyVector <int,NDIM+2> hashSize; // The current size of the hash array
118 
119  beadLocator swap; // Used when deleting/updating beads
120 
121  void setupNNGrid(); // We initialize the nearest neighbor grid
122 
123  /* Return tiny vectors suitable for indexing the numLabel and hash array */
124  inline TinyVector <int,NDIM+1> numLabelIndex(const beadLocator&);
125  inline TinyVector <int,NDIM+1> numLabelIndex(const iVec&, const int);
126  inline TinyVector <int,NDIM+2> hashIndex(const beadLocator&, const int);
127  inline TinyVector <int,NDIM+2> hashIndex(const iVec&, const int, const int);
128 };
129 
130 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
131 // INLINE FUNCTION DEFINITIONS
132 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
133 
134 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
139 inline iVec LookupTable::gridIndex(const dVec &pos) {
140  iVec index; // The grid index in each dimension
141  for (int i = 0; i < NDIM; i++) {
142  index[i] = static_cast<int>( abs( pos[i] + 0.5 * boxPtr->side[i] - EPS )
143  / (sizeNNGrid[i] + EPS) );
144  PIMC_ASSERT(index[i]<numNNGrid[i]);
145  }
146  return index;
147 }
148 
149 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
154 inline int LookupTable::gridNumber(const dVec &pos) {
155  int gNumber = 0;
156  for (int i = 0; i < NDIM; i++) {
157  int scale = 1;
158  for (int j = i+1; j < NDIM; j++)
159  scale *= numNNGrid[j];
160  gNumber += scale *
161  static_cast<int>(abs( pos[i] + 0.5 * boxPtr->side[i] - EPS ) / (sizeNNGrid[i] + EPS));
162  }
163  PIMC_ASSERT(gNumber<totNumGridBoxes);
164  return gNumber;
165 }
166 
167 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
172  iVec _gridIndex;
173  for (int i = 0; i < NDIM; i++) {
174  int scale = 1;
175  for (int j = i+1; j < NDIM; j++)
176  scale *= numNNGrid[j];
177  _gridIndex[i] = (n/scale) % numNNGrid[i];
178  PIMC_ASSERT(_gridIndex[i]<numNNGrid[i]);
179  }
180  return _gridIndex;
181 }
182 
183 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
187 inline int LookupTable::gridNumber(const iVec &index) {
188  int gNumber = 0;
189  for (int i = 0; i < NDIM; i++) {
190  int scale = 1;
191  for (int j = i+1; j < NDIM; j++)
192  scale *= numNNGrid[j];
193  gNumber += index[i]*scale;
194  }
195  PIMC_ASSERT(gNumber<totNumGridBoxes);
196  return gNumber;
197 }
198 
199 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
204 inline TinyVector <int,NDIM+1> LookupTable::numLabelIndex(const beadLocator &beadIndex) {
205  TinyVector<int,NDIM+1> index;
206  for (int i = 0; i < NDIM; i++)
207  index[i] = grid(beadIndex)[i];
208  index[NDIM] = beadIndex[0];
209  return index;
210 }
211 
212 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
217 inline TinyVector <int,NDIM+1> LookupTable::numLabelIndex(const iVec &gI, const int slice) {
218  TinyVector<int,NDIM+1> index;
219  for (int i = 0; i < NDIM; i++)
220  index[i] = gI[i];
221  index[NDIM] = slice;
222  return index;
223 }
224 
225 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
230 inline TinyVector <int,NDIM+2> LookupTable::hashIndex(const beadLocator &beadIndex,
231  const int label) {
232  TinyVector<int,NDIM+2> index;
233  for (int i = 0; i < NDIM; i++)
234  index[i] = grid(beadIndex)[i];
235  index[NDIM] = beadIndex[0];
236  index[NDIM+1] = label;
237  return index;
238 }
239 
240 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
245 inline TinyVector <int,NDIM+2> LookupTable::hashIndex(const iVec &gI, const int slice,
246  const int label) {
247  TinyVector<int,NDIM+2> index;
248  for (int i = 0; i < NDIM; i++)
249  index[i] = gI[i];
250  index[NDIM] = slice;
251  index[NDIM+1] = label;
252  return index;
253 }
254 
255 #endif
The base class which holds details on the generalized box that our system will be simulated inside of...
Definition: container.h:24
dVec side
The linear dimensions of the box.
Definition: container.h:31
The particle (bead) lookup table.
Definition: lookuptable.h:29
int fullNumBeads
The full number of active beads in beadList;.
Definition: lookuptable.h:44
int numUniqueNN
The number of unique nearest neighbors of each box.
Definition: lookuptable.h:37
int getTotNumGridBoxes()
Return the total number of grid boxes.
Definition: lookuptable.h:49
int numNN
The total number of nearest neighbors of each box.
Definition: lookuptable.h:38
iVec gridIndex(const dVec &)
Given a particle position, return the grid index for the nearest neighbor lookup table.
Definition: lookuptable.h:139
void updateGrid(const Path &)
We update the full nearest neighbor grid by filling up all data arrays at all time slices.
Array< beadLocator, 1 > beadList
The cutoff dynamic list of interacting beads.
Definition: lookuptable.h:40
void updateFullInteractionList(const beadLocator &, const int)
Fill up the fullBeadList array with a list of beads in the same grid box as the supplied beadIndex an...
void updateInteractionList(const Path &, const beadLocator &)
Update the NN lookup table and the array of beadLocators containing all beads which 'interact' with t...
void resizeList(int _numParticles)
Resize the bead and grid lists.
Definition: lookuptable.h:78
iVec getNumNNGrid()
Return the number of NN grid boxes.
Definition: lookuptable.h:47
LookupTable(const Container *, const int, const int)
Initilialize the nearest neighbor lookup table.
Definition: lookuptable.cpp:26
void addBead(const beadLocator &, const dVec &)
Add a single bead to the NN grid and position pos.
Array< dVec, 1 > beadSep
The separation between beads.
Definition: lookuptable.h:42
bool gridNeighbors(const beadLocator &, const beadLocator &)
Given two beadIndices, determine if the beads lie in neighboring grid boxes.
int numBeads
The cutoff number of active beads in beadList;.
Definition: lookuptable.h:43
~LookupTable()
Free all blitz arrays.
Definition: lookuptable.cpp:60
Array< beadLocator, 1 > fullBeadList
The full dynamic list of interacting beads.
Definition: lookuptable.h:41
const Container * boxPtr
The simulation cell.
Definition: lookuptable.h:35
bool gridShare(const beadLocator &bead1, const beadLocator &bead2)
Determine if two beads are in the same grid box.
Definition: lookuptable.h:90
void updateBead(const beadLocator &, const dVec &)
Update the NN grid with a new bead position.
void delBead(const beadLocator &)
Remove a single bead from the NN grid.
int gridNumber(const dVec &)
Given a particle position, return the grid number for the nearest neighbor lookup table.
Definition: lookuptable.h:154
void printGrid()
Print the NN Lookup table.
The space-time trajectories.
Definition: path.h:29
Global common header with shared dependencies and methods.
#define NDIM
Number of spatial dimnsions.
Definition: common.h:71
TinyVector< int, 2 > beadLocator
time-slice,bead-number world line index
Definition: common.h:117
#define EPS
A small number.
Definition: common.h:94
TinyVector< double, NDIM > dVec
A NDIM-vector of type double.
Definition: common.h:111
#define PIMC_ASSERT(X)
Rename assert method.
Definition: common.h:64
TinyVector< int, NDIM > iVec
A NDIM-vector of type integer.
Definition: common.h:114
ConstantParameters class definition.
The simulation cell.