diff --git a/CHANGELOG.md b/CHANGELOG.md index 8a69281fe45f..39e4edf1ce00 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -3,6 +3,7 @@ ## Current develop ### Added (new features/APIs/variables/...) +- [[PR 1184]](https://github.com/parthenon-hpc-lab/parthenon/pull/1184) Fix swarm block neighbor indexing in 1D, 2D - [[PR 1183]](https://github.com/parthenon-hpc-lab/parthenon/pull/1183) Fix particle leapfrog example initialization data - [[PR 1179]](https://github.com/parthenon-hpc-lab/parthenon/pull/1179) Make a global variable for whether simulation is a restart - [[PR 1171]](https://github.com/parthenon-hpc-lab/parthenon/pull/1171) Add PARTHENON_USE_SYSTEM_PACKAGES build option diff --git a/example/particles/particles.cpp b/example/particles/particles.cpp index 3f99407ded54..77c0c8ca585b 100644 --- a/example/particles/particles.cpp +++ b/example/particles/particles.cpp @@ -307,6 +307,14 @@ TaskStatus CreateSomeParticles(MeshBlock *pmb, const double t0) { weight(n) = 1.0; + // Check that we are on current meshblock + bool is_on_current_mesh_block = false; + PARTHENON_REQUIRE(swarm_d.GetNeighborBlockIndex(n, x(n), y(n), z(n), + is_on_current_mesh_block) == -1, + "Particle must be on current meshblock!"); + PARTHENON_REQUIRE(is_on_current_mesh_block == true, + "Particle must be on current meshblock!"); + rng_pool.free_state(rng_gen); }); } else { @@ -333,6 +341,14 @@ TaskStatus CreateSomeParticles(MeshBlock *pmb, const double t0) { weight(n) = 1.0; + // Check that we are on current meshblock + bool is_on_current_mesh_block = false; + PARTHENON_REQUIRE(swarm_d.GetNeighborBlockIndex(n, x(n), y(n), z(n), + is_on_current_mesh_block) == -1, + "Particle must be on current meshblock!"); + PARTHENON_REQUIRE(is_on_current_mesh_block == true, + "Particle must be on current meshblock!"); + rng_pool.free_state(rng_gen); }); } diff --git a/src/interface/swarm_comms.cpp b/src/interface/swarm_comms.cpp index 650bcef09304..3ee4e2af7512 100644 --- a/src/interface/swarm_comms.cpp +++ b/src/interface/swarm_comms.cpp @@ -44,21 +44,6 @@ void Swarm::SetNeighborIndices_() { } } - // Indicate which neighbor regions correspond to this meshblock - const int kmin = ndim < 3 ? 0 : 1; - const int kmax = ndim < 3 ? 3 : 2; - const int jmin = ndim < 2 ? 0 : 1; - const int jmax = ndim < 2 ? 3 : 2; - const int imin = 1; - const int imax = 2; - for (int k = kmin; k <= kmax; k++) { - for (int j = jmin; j <= jmax; j++) { - for (int i = imin; i <= imax; i++) { - neighbor_indices_h(k, j, i) = this_block_; - } - } - } - // Create a point in the center of each ghost halo region at maximum refinement level // and then test whether each neighbor block includes that point. const auto &bsize = pmb->block_size; @@ -130,6 +115,21 @@ void Swarm::SetNeighborIndices_() { } } + // Indicate which neighbor regions correspond to this meshblock + const int kmin = ndim < 3 ? 0 : 1; + const int kmax = ndim < 3 ? 3 : 2; + const int jmin = ndim < 2 ? 0 : 1; + const int jmax = ndim < 2 ? 3 : 2; + const int imin = 1; + const int imax = 2; + for (int k = kmin; k <= kmax; k++) { + for (int j = jmin; j <= jmax; j++) { + for (int i = imin; i <= imax; i++) { + neighbor_indices_h(k, j, i) = this_block_; + } + } + } + neighbor_indices_.DeepCopy(neighbor_indices_h); } @@ -437,6 +437,9 @@ void Swarm::AllocateComms(std::weak_ptr wpmb) { // Enroll SwarmVariable object vbswarm->bswarm_index = pmb->pbswarm->bswarms.size(); pmb->pbswarm->bswarms.push_back(vbswarm); + + // Evaluate neigbor indices + SetNeighborIndices_(); } } // namespace parthenon diff --git a/src/interface/swarm_device_context.hpp b/src/interface/swarm_device_context.hpp index 1dbcf3383901..936d2d56ad35 100644 --- a/src/interface/swarm_device_context.hpp +++ b/src/interface/swarm_device_context.hpp @@ -75,9 +75,9 @@ class SwarmDeviceContext { // Ignore k,j indices as necessary based on problem dimension if (ndim_ == 1) { - block_index_(n) = neighbor_indices_(0, 0, i); + block_index_(n) = neighbor_indices_(1, 1, i); } else if (ndim_ == 2) { - block_index_(n) = neighbor_indices_(0, j, i); + block_index_(n) = neighbor_indices_(1, j, i); } else { block_index_(n) = neighbor_indices_(k, j, i); }