| 98 | } |
| 99 | |
| 100 | CellListResult CellList::getNeighboursForPosition(const double x, const double y, const double z) const |
| 101 | { |
| 102 | // The indices of the neighbouring atoms |
| 103 | vector<int> neighbours; |
| 104 | vector<double> distances; |
| 105 | vector<double> distancesSquared; |
| 106 | |
| 107 | // Find bin for the given position |
| 108 | int i0 = (x - this->xmin)/this->dx; |
| 109 | int j0 = (y - this->ymin)/this->dy; |
| 110 | int k0 = (z - this->zmin)/this->dz; |
| 111 | |
| 112 | // Get the bin ranges to check for each dimension. |
| 113 | int istart = max(i0-1, 0); |
| 114 | int iend = min(i0+1, this->nx-1); |
| 115 | int jstart = max(j0-1, 0); |
| 116 | int jend = min(j0+1, this->ny-1); |
| 117 | int kstart = max(k0-1, 0); |
| 118 | int kend = min(k0+1, this->nz-1); |
| 119 | |
| 120 | // Loop over neighbouring bins |
| 121 | for (int i = istart; i <= iend; i++){ |
| 122 | for (int j = jstart; j <= jend; j++){ |
| 123 | for (int k = kstart; k <= kend; k++){ |
| 124 | |
| 125 | // For each atom in the current bin, calculate the actual distance |
| 126 | vector<int> binIndices = this->bins[i][j][k]; |
| 127 | for (auto &idx : binIndices) { |
| 128 | double ix = this->positions(idx, 0); |
| 129 | double iy = this->positions(idx, 1); |
| 130 | double iz = this->positions(idx, 2); |
| 131 | double deltax = x - ix; |
| 132 | double deltay = y - iy; |
| 133 | double deltaz = z - iz; |
| 134 | double distanceSquared = deltax*deltax + deltay*deltay + deltaz*deltaz; |
| 135 | if (distanceSquared <= this->cutoffSquared) { |
| 136 | neighbours.push_back(idx); |
| 137 | distancesSquared.push_back(distanceSquared); |
| 138 | distances.push_back(sqrt(distanceSquared)); |
| 139 | } |
| 140 | } |
| 141 | } |
| 142 | } |
| 143 | } |
| 144 | return CellListResult{neighbours, distances, distancesSquared}; |
| 145 | } |
| 146 | |
| 147 | CellListResult CellList::getNeighboursForIndex(const int idx) const |
| 148 | { |
no test coverage detected