Skip to content

Commit

Permalink
fix copies to rank zero for correctness tests (#531)
Browse files Browse the repository at this point in the history
### Description
This fixes parallel copies of particles and/or MultiFabs to rank zero so
that they actually work. In this case, constructing a
`amrex::DistributionMapping` must be done by explicitly setting it equal
to `{0}`.

This is used when testing the correctness of the results for the
`HydroRichtmyerMeshkov` and `BinaryOrbitCIC` tests.

Thanks to Andrew Myers for the fix.

### Related issues
N/A

### Checklist
_Before this pull request can be reviewed, all of these tasks should be
completed. Denote completed tasks with an `x` inside the square brackets
`[ ]` in the Markdown source below:_
- [x] I have added a description (see above).
- [x] I have added a link to any related issues see (see above).
- [x] I have read the [Contributing
Guide](https://github.com/quokka-astro/quokka/blob/development/CONTRIBUTING.md).
- [ ] I have added tests for any new physics that this PR adds to the
code.
- [x] I have tested this PR on my local computer and all tests pass.
- [x] I have manually triggered the GPU tests with the magic comment
`/azp run`.
- [x] I have requested a reviewer for this PR.
  • Loading branch information
BenWibking authored Feb 7, 2024
1 parent edc1dbc commit a16b88f
Show file tree
Hide file tree
Showing 2 changed files with 16 additions and 18 deletions.
6 changes: 3 additions & 3 deletions src/BinaryOrbitCIC/binary_orbit.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -101,14 +101,14 @@ template <> void RadhydroSimulation<BinaryOrbit>::computeAfterTimestep()
amrex::Box const box(amrex::IntVect{AMREX_D_DECL(0, 0, 0)}, amrex::IntVect{AMREX_D_DECL(1, 1, 1)});
amrex::Geometry const geom(box);
amrex::BoxArray const boxArray(box);
amrex::DistributionMapping const dmap(boxArray, 1);
amrex::Vector<int> const ranks({0}); // workaround nvcc bug
amrex::DistributionMapping const dmap(ranks);
analysisPC.Define(geom, dmap, boxArray);
analysisPC.copyParticles(*CICParticles);
// do we need to redistribute??

if (amrex::ParallelDescriptor::IOProcessor()) {
quokka::CICParticleIterator const pIter(analysisPC, 0);
if (pIter.isValid()) { // this returns false when there is more than 1 MPI rank (?)
if (pIter.isValid()) {
amrex::Print() << "Computing particle statistics...\n";
const amrex::Long np = pIter.numParticles();
auto &particles = pIter.GetArrayOfStructs();
Expand Down
28 changes: 13 additions & 15 deletions src/HydroRichtmeyerMeshkov/test_hydro2d_rm.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -42,34 +42,33 @@ template <> struct Physics_Traits<RichtmeyerMeshkovProblem> {
static constexpr int nGroups = 1; // number of radiation groups
};

// #define DEBUG_SYMMETRY
template <> void RadhydroSimulation<RichtmeyerMeshkovProblem>::computeAfterTimestep()
{
#ifdef DEBUG_SYMMETRY
// this code does not actually work with Nranks > 1 ...
const int ncomp_cc = Physics_Indices<RichtmeyerMeshkovProblem>::nvarTotal_cc;

// copy all FABs to a local FAB across the entire domain
amrex::BoxArray localBoxes(domain_);
amrex::DistributionMapping localDistribution(localBoxes, 1);
// copy all FAB data to a single FAB on rank zero
amrex::Box const domainBox = geom[0].Domain();
amrex::BoxArray const localBoxes(domainBox);
amrex::Vector<int> const ranks({0}); // workaround nvcc bug
amrex::DistributionMapping const localDistribution(ranks);
amrex::MultiFab state_mf(localBoxes, localDistribution, ncomp_cc, 0);
state_mf.ParallelCopy(state_new_cc_);
state_mf.ParallelCopy(state_new_cc_[0]);

if (amrex::ParallelDescriptor::IOProcessor()) {
auto const &state = state_mf.array(0);
auto const prob_lo = simGeometry_.ProbLoArray();
auto dx = dx_;
auto const prob_lo = geom[0].ProbLoArray();
auto const dx = geom[0].CellSizeArray();

amrex::Long asymmetry = 0;
auto nx = nx_;
auto ny = ny_;
auto nz = nz_;
auto nx = domainBox.length(0);
auto ny = domainBox.length(1);
auto nz = domainBox.length(2);
auto ncomp = ncomp_cc;
for (int i = 0; i < nx; ++i) {
for (int j = 0; j < ny; ++j) {
for (int k = 0; k < nz; ++k) {
amrex::Real const x = prob_lo[0] + (i + amrex::Real(0.5)) * dx[0];
amrex::Real const y = prob_lo[1] + (j + amrex::Real(0.5)) * dx[1];
amrex::Real const x = prob_lo[0] + (i + static_cast<amrex::Real>(0.5)) * dx[0];
amrex::Real const y = prob_lo[1] + (j + static_cast<amrex::Real>(0.5)) * dx[1];
for (int n = 0; n < ncomp; ++n) {
const amrex::Real comp_upper = state(i, j, k, n);

Expand Down Expand Up @@ -100,7 +99,6 @@ template <> void RadhydroSimulation<RichtmeyerMeshkovProblem>::computeAfterTimes
}
AMREX_ALWAYS_ASSERT_WITH_MESSAGE(asymmetry == 0, "x/y not symmetric!");
}
#endif
}

template <> void RadhydroSimulation<RichtmeyerMeshkovProblem>::setInitialConditionsOnGrid(quokka::grid grid_elem)
Expand Down

0 comments on commit a16b88f

Please sign in to comment.