diff --git a/src/core/cells.cpp b/src/core/cells.cpp index ca21de4072b..9fe1ff2ba0c 100644 --- a/src/core/cells.cpp +++ b/src/core/cells.cpp @@ -306,33 +306,6 @@ void fold_and_reset(Particle &p) { p.l.p_old = p.r.p; } - -/** - * @brief Extract an indexed particle from a list. - * - * Removes a particle from a particle list and - * from the particle index. - */ -Particle extract_indexed_particle(ParticleList *sl, int i) { - Particle *src = &sl->part[i]; - Particle *end = &sl->part[sl->n - 1]; - - Particle p = std::move(*src); - - assert(p.p.identity <= max_seen_particle); - local_particles[p.p.identity] = nullptr; - - if (src != end) { - new (src) Particle(std::move(*end)); - } - - if (realloc_particlelist(sl, --sl->n)) { - update_local_particles(sl); - } else if (src != end) { - local_particles[src->p.identity] = src; - } - return p; -} } // namespace /** @@ -491,13 +464,11 @@ void cells_update_ghosts() { } Cell *find_current_cell(const Particle &p) { - auto c = cell_structure.position_to_cell(p.r.p); - if (c) { - return c; - } else if (!p.l.ghost) { - // Old pos must lie within the cell system - return cell_structure.position_to_cell(p.l.p_old); - } else { + assert(not resort_particles); + + if (p.l.ghost) { return nullptr; } + + return cell_structure.position_to_cell(p.l.p_old); } diff --git a/src/core/cells.hpp b/src/core/cells.hpp index 3487d57264e..7ed2a331216 100644 --- a/src/core/cells.hpp +++ b/src/core/cells.hpp @@ -298,13 +298,12 @@ void check_resort_particles(); /*@}*/ -/* @brief Finds the cell in which a particle is stored - - Uses position_to_cell on p.r.p. If this is not on the node's domain, - uses position at last Verlet list rebuild (p.l.p_old). - - @return pointer to the cell or nullptr if the particle is not on the node -*/ +/** + * @brief Finds the cell in which a particle is stored + * + * + * @return pointer to the cell or nullptr if the particle is not on the node + */ Cell *find_current_cell(const Particle &p); #endif diff --git a/src/core/particle_data.cpp b/src/core/particle_data.cpp index c0016b96d1c..a7b705499d2 100644 --- a/src/core/particle_data.cpp +++ b/src/core/particle_data.cpp @@ -220,6 +220,7 @@ void init_particlelist(ParticleList *pList) { } int realloc_particlelist(ParticleList *l, int size) { + assert(size >= 0); int old_max = l->max; Particle *old_start = l->part; @@ -280,6 +281,9 @@ Particle *append_indexed_particle(ParticleList *l, Particle &&part) { } Particle *move_unindexed_particle(ParticleList *dl, ParticleList *sl, int i) { + assert(sl->n > 0); + assert(i < sl->n); + realloc_particlelist(dl, ++dl->n); auto dst = &dl->part[dl->n - 1]; auto src = &sl->part[i]; @@ -290,12 +294,13 @@ Particle *move_unindexed_particle(ParticleList *dl, ParticleList *sl, int i) { new (src) Particle(std::move(*end)); } - sl->n -= 1; - realloc_particlelist(sl, sl->n); + realloc_particlelist(sl, --sl->n); return dst; } Particle *move_indexed_particle(ParticleList *dl, ParticleList *sl, int i) { + assert(sl->n > 0); + assert(i < sl->n); int re = realloc_particlelist(dl, ++dl->n); Particle *dst = &dl->part[dl->n - 1]; Particle *src = &sl->part[i]; @@ -322,6 +327,41 @@ Particle *move_indexed_particle(ParticleList *dl, ParticleList *sl, int i) { return dst; } +/** + * @brief Extract an indexed particle from a list. + * + * Removes a particle from a particle list and + * from the particle index. + * + * @param i Index of particle to remove, + * needs to be valid. + * @param sl List to remove the particle from, + * needs to be non-empty. + * @return The extracted particle. + */ +Particle extract_indexed_particle(ParticleList *sl, int i) { + assert(sl->n > 0); + assert(i < sl->n); + Particle *src = &sl->part[i]; + Particle *end = &sl->part[sl->n - 1]; + + Particle p = std::move(*src); + + assert(p.p.identity <= max_seen_particle); + local_particles[p.p.identity] = nullptr; + + if (src != end) { + new (src) Particle(std::move(*end)); + } + + if (realloc_particlelist(sl, --sl->n)) { + update_local_particles(sl); + } else if (src != end) { + local_particles[src->p.identity] = src; + } + return p; +} + namespace { /* Limit cache to 100 MiB */ std::size_t const max_cache_size = (100ul * 1048576ul) / sizeof(Particle); @@ -803,42 +843,50 @@ int remove_particle(int p_id) { return ES_OK; } -void local_remove_particle(int part) { - int ind, c; - Particle *p = local_particles[part]; - ParticleList *pl = nullptr, *tmp; - - /* the tricky - say ugly - part: determine - the cell the particle is located in by checking - whether the particle address is inside the array */ - for (c = 0; c < local_cells.n; c++) { - tmp = local_cells.cell[c]; - ind = p - tmp->part; - if (ind >= 0 && ind < tmp->n) { - pl = tmp; - break; +namespace { +std::pair find_particle(Particle *p, Cell *c) { + for (int i = 0; i < c->n; ++i) { + if ((c->part + i) == p) { + return {c, i}; } } - if (!pl) { - fprintf(stderr, - "%d: INTERNAL ERROR: could not find cell of particle %d, exiting\n", - this_node, part); - errexit(); - } + return {nullptr, 0}; +} - free_particle(p); +std::pair find_particle(Particle *p, CellPList cells) { + for (auto &c : cells) { + auto res = find_particle(p, c); + if (res.first) { + return res; + } + } - /* remove local_particles entry */ - local_particles[p->p.identity] = nullptr; + return {nullptr, 0}; +} +} // namespace - if (&pl->part[pl->n - 1] != p) { - /* move last particle to free position */ - *p = pl->part[pl->n - 1]; +void local_remove_particle(int part) { + Particle *p = local_particles[part]; + assert(p); + assert(not p->l.ghost); + + /* If the particles are sorted we can use the + * cell system to find the cell containing the + * particle. Otherwise we do a brute force search + * of the cells. */ + Cell *cell = nullptr; + size_t n = 0; + if (Cells::RESORT_NONE == get_resort_particles()) { + std::tie(cell, n) = find_particle(p, find_current_cell(*p)); + } - /* update the local_particles array for the moved particle */ - local_particles[p->p.identity] = p; + if (not cell) { + std::tie(cell, n) = find_particle(p, local_cells); } - pl->n--; + + assert(cell && cell->part && (n < cell->n) && ((cell->part + n) == p)); + + Particle p_destroy = extract_indexed_particle(cell, n); } void local_place_particle(int part, const double p[3], int _new) { diff --git a/src/core/particle_data.hpp b/src/core/particle_data.hpp index fcf5b9b1bc1..3cc090c688a 100644 --- a/src/core/particle_data.hpp +++ b/src/core/particle_data.hpp @@ -562,6 +562,8 @@ Particle *move_unindexed_particle(ParticleList *destList, Particle *move_indexed_particle(ParticleList *destList, ParticleList *sourceList, int ind); +Particle extract_indexed_particle(ParticleList *sl, int i); + /* Other Functions */ /************************************************/