From 0a8592912a6c13704177234719e7d97c1ede6e2c Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Fri, 15 May 2020 00:25:31 +0200 Subject: [PATCH] Bugfix: Cells may not receive any variable size data. --- source/distributed/tria.cc | 8 +- tests/particles/particle_handler_14.cc | 96 +++++++++++++++++++ ...handler_14.with_p4est=true.mpirun=2.output | 4 + 3 files changed, 104 insertions(+), 4 deletions(-) create mode 100644 tests/particles/particle_handler_14.cc create mode 100644 tests/particles/particle_handler_14.with_p4est=true.mpirun=2.output diff --git a/source/distributed/tria.cc b/source/distributed/tria.cc index ff7c033397..f0f395dde4 100644 --- a/source/distributed/tria.cc +++ b/source/distributed/tria.cc @@ -1611,16 +1611,16 @@ namespace parallel const bool callback_variable_transfer = (handle % 2 == 0); const unsigned int callback_index = handle / 2; + // Cells will always receive fixed size data (i.e., CellStatus + // information), but not necessarily variable size data (e.g., with a + // ParticleHandler a cell might not contain any particle at all). + // Thus it is sufficient to check if fixed size data has been received. Assert(sizes_fixed_cumulative.size() > 0, ExcMessage("No data has been packed!")); if (quad_cell_relations.size() > 0) { Assert(dest_data_fixed.size() > 0, ExcMessage("No data has been received!")); - - if (callback_variable_transfer) - Assert(dest_data_variable.size() > 0, - ExcMessage("No data has been received!")); } std::vector::const_iterator dest_data_it; diff --git a/tests/particles/particle_handler_14.cc b/tests/particles/particle_handler_14.cc new file mode 100644 index 0000000000..ae31bf7713 --- /dev/null +++ b/tests/particles/particle_handler_14.cc @@ -0,0 +1,96 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// check variable size data transfer for parallel distributed applications +// if one cell does not contain any particle at all + + +#include + +#include + +#include + +#include + +#include "../tests.h" + +template +void +test() +{ + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + + MappingQ mapping(1); + + // setup triangulation + const unsigned int n_cells = 2; + + std::vector rep(dim, 1); + rep[0] = n_cells; + Point p1, p2; + for (unsigned int d = 0; d < dim; ++d) + { + p1[d] = 0; + p2[d] = (d == 0) ? n_cells : 1; + } + GridGenerator::subdivided_hyper_rectangle(tr, rep, p1, p2); + + // setup one particle in first cell + Particles::ParticleHandler particle_handler(tr, mapping); + + Particles::Particle particle(Point(), + Point(), + 0); + + particle_handler.insert_particle(particle, tr.begin_active()); + particle_handler.update_cached_numbers(); + + // initiate data transfer + tr.signals.pre_distributed_repartition.connect([&particle_handler]() { + particle_handler.register_store_callback_function(); + }); + + tr.signals.post_distributed_repartition.connect([&particle_handler]() { + particle_handler.register_load_callback_function(false); + }); + + tr.repartition(); + + deallog << "OK" << std::endl; +} + + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0) + initlog(); + + deallog.push("2d/2d"); + test<2, 2>(); + deallog.pop(); + deallog.push("2d/3d"); + test<2, 3>(); + deallog.pop(); + deallog.push("3d/3d"); + test<3, 3>(); + deallog.pop(); +} diff --git a/tests/particles/particle_handler_14.with_p4est=true.mpirun=2.output b/tests/particles/particle_handler_14.with_p4est=true.mpirun=2.output new file mode 100644 index 0000000000..988eda8574 --- /dev/null +++ b/tests/particles/particle_handler_14.with_p4est=true.mpirun=2.output @@ -0,0 +1,4 @@ + +DEAL:2d/2d::OK +DEAL:2d/3d::OK +DEAL:3d/3d::OK -- 2.39.5