Fix for particles in empty local triangulations

diff --git a/doc/news/changes/minor/20210505Gassmoeller b/doc/news/changes/minor/20210505Gassmoeller
new file mode 100644 (file)
index 0000000..ec4cf86
--- /dev/null
+++ b/doc/news/changes/minor/20210505Gassmoeller
@@ -0,0 +1,4 @@
+Fixed: The function Particles::Generators::probabilistic_locations could crash
+if some MPI ranks had no active cells. This is fixed now.
+<br>
+(Rene Gassmoeller, 2021/05/05)

diff --git a/source/particles/generators.cc b/source/particles/generators.cc
index 60c78be8c075c5590c8b733d014933bf748d6d9a..d013e4dc8669ab1c9d8b0b250719ff4083742e9f 100644 (file)
--- a/source/particles/generators.cc
+++ b/source/particles/generators.cc
@@ -261,7 +261,10 @@ namespace Particles
                                                 probability_density_function);
 
         // Sum the local integrals over all nodes
-        double local_weight_integral = cumulative_cell_weights.back();
+        double local_weight_integral = (cumulative_cell_weights.size() > 0) ?
+                                         cumulative_cell_weights.back() :
+                                         0.0;
+
         double global_weight_integral;
 
         if (const auto tria =

diff --git a/tests/particles/generators_13.cc b/tests/particles/generators_13.cc
new file mode 100644 (file)
index 0000000..9718103
--- /dev/null
+++ b/tests/particles/generators_13.cc
@@ -0,0 +1,107 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 - 2019 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// like generators_05, but check that empty triangulations do not crash
+// (by using more cores than cells).
+
+#include <deal.II/base/function_lib.h>
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_nothing.h>
+#include <deal.II/fe/mapping_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+
+#include <deal.II/particles/generators.h>
+#include <deal.II/particles/particle_handler.h>
+
+#include "../tests.h"
+
+template <int dim, int spacedim>
+void
+test(Function<spacedim> &probability_density_function)
+{
+  {
+    parallel::distributed::Triangulation<dim, spacedim> tr(MPI_COMM_WORLD);
+
+    GridGenerator::hyper_cube(tr);
+    MappingQ<dim, spacedim> mapping(1);
+
+    Particles::ParticleHandler<dim, spacedim> particle_handler(tr, mapping);
+
+    Particles::Generators::probabilistic_locations(
+      tr, probability_density_function, true, 16, particle_handler, mapping);
+
+    deallog << "Rank: " << Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)
+            << ". Locally owned active cells: "
+            << tr.n_locally_owned_active_cells() << ". Local particles: "
+            << particle_handler.n_locally_owned_particles() << std::endl;
+
+    for (const auto &cell : tr.active_cell_iterators())
+      {
+        if (cell->is_locally_owned())
+          {
+            deallog << "Cell " << cell << " has "
+                    << particle_handler.n_particles_in_cell(cell)
+                    << " particles." << std::endl;
+          }
+      }
+
+    for (const auto &particle : particle_handler)
+      {
+        deallog << "Particle index " << particle.get_id() << " is in cell "
+                << particle.get_surrounding_cell(tr) << std::endl;
+        deallog << "Particle location: " << particle.get_location()
+                << std::endl;
+      }
+  }
+
+  deallog << "OK" << std::endl;
+}
+
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  MPILogInitAll init;
+
+  {
+    deallog.push("2d/2d");
+    Functions::ConstantFunction<2> probability_density_function(1.0);
+    test<2, 2>(probability_density_function);
+    deallog.pop();
+  }
+  {
+    deallog.push("2d/3d");
+    Functions::ConstantFunction<3> probability_density_function(1.0);
+    test<2, 3>(probability_density_function);
+  }
+  {
+    deallog.pop();
+    deallog.push("3d/3d");
+    Functions::ConstantFunction<3> probability_density_function(1.0);
+    test<3, 3>(probability_density_function);
+    deallog.pop();
+  }
+}

diff --git a/tests/particles/generators_13.with_p4est=true.mpirun=2.output b/tests/particles/generators_13.with_p4est=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..42fb023
--- /dev/null
+++ b/tests/particles/generators_13.with_p4est=true.mpirun=2.output
@@ -0,0 +1,114 @@
+
+DEAL:0:2d/2d::Rank: 0. Locally owned active cells: 0. Local particles: 0
+DEAL:0:2d/2d::OK
+DEAL:0:2d/3d::Rank: 0. Locally owned active cells: 0. Local particles: 0
+DEAL:0:2d/3d::OK
+DEAL:0:3d/3d::Rank: 0. Locally owned active cells: 0. Local particles: 0
+DEAL:0:3d/3d::OK
+
+DEAL:1:2d/2d::Rank: 1. Locally owned active cells: 1. Local particles: 16
+DEAL:1:2d/2d::Cell 0.0 has 16 particles.
+DEAL:1:2d/2d::Particle index 0 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.246806 0.602738
+DEAL:1:2d/2d::Particle index 1 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.756380 0.916857
+DEAL:1:2d/2d::Particle index 2 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.595957 0.950968
+DEAL:1:2d/2d::Particle index 3 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.793715 0.387718
+DEAL:1:2d/2d::Particle index 4 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.731510 0.439705
+DEAL:1:2d/2d::Particle index 5 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.135658 0.326733
+DEAL:1:2d/2d::Particle index 6 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.251525 0.339426
+DEAL:1:2d/2d::Particle index 7 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.595620 0.859503
+DEAL:1:2d/2d::Particle index 8 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.362813 0.242154
+DEAL:1:2d/2d::Particle index 9 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.299327 0.941832
+DEAL:1:2d/2d::Particle index 10 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.910383 0.872167
+DEAL:1:2d/2d::Particle index 11 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.375099 0.884627
+DEAL:1:2d/2d::Particle index 12 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.138433 0.419489
+DEAL:1:2d/2d::Particle index 13 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.692605 0.474755
+DEAL:1:2d/2d::Particle index 14 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.382099 0.844530
+DEAL:1:2d/2d::Particle index 15 is in cell 0.0
+DEAL:1:2d/2d::Particle location: 0.319646 0.518143
+DEAL:1:2d/2d::OK
+DEAL:1:2d/3d::Rank: 1. Locally owned active cells: 1. Local particles: 16
+DEAL:1:2d/3d::Cell 0.0 has 16 particles.
+DEAL:1:2d/3d::Particle index 0 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.246806 0.602738 0.00000
+DEAL:1:2d/3d::Particle index 1 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.916857 0.595957 0.00000
+DEAL:1:2d/3d::Particle index 2 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.793715 0.387718 0.00000
+DEAL:1:2d/3d::Particle index 3 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.439705 0.135658 0.00000
+DEAL:1:2d/3d::Particle index 4 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.251525 0.339426 0.00000
+DEAL:1:2d/3d::Particle index 5 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.859503 0.362813 0.00000
+DEAL:1:2d/3d::Particle index 6 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.299327 0.941832 0.00000
+DEAL:1:2d/3d::Particle index 7 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.872167 0.375099 0.00000
+DEAL:1:2d/3d::Particle index 8 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.138433 0.419489 0.00000
+DEAL:1:2d/3d::Particle index 9 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.474755 0.382099 0.00000
+DEAL:1:2d/3d::Particle index 10 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.319646 0.518143 0.00000
+DEAL:1:2d/3d::Particle index 11 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.456149 0.635687 0.00000
+DEAL:1:2d/3d::Particle index 12 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.588260 0.522915 0.00000
+DEAL:1:2d/3d::Particle index 13 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.256491 0.534039 0.00000
+DEAL:1:2d/3d::Particle index 14 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.957838 0.850786 0.00000
+DEAL:1:2d/3d::Particle index 15 is in cell 0.0
+DEAL:1:2d/3d::Particle location: 0.635976 0.937905 0.00000
+DEAL:1:2d/3d::OK
+DEAL:1:3d/3d::Rank: 1. Locally owned active cells: 1. Local particles: 16
+DEAL:1:3d/3d::Cell 0.0 has 16 particles.
+DEAL:1:3d/3d::Particle index 0 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.246806 0.602738 0.756380
+DEAL:1:3d/3d::Particle index 1 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.916857 0.595957 0.950968
+DEAL:1:3d/3d::Particle index 2 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.793715 0.387718 0.731510
+DEAL:1:3d/3d::Particle index 3 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.439705 0.135658 0.326733
+DEAL:1:3d/3d::Particle index 4 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.251525 0.339426 0.595620
+DEAL:1:3d/3d::Particle index 5 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.859503 0.362813 0.242154
+DEAL:1:3d/3d::Particle index 6 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.299327 0.941832 0.910383
+DEAL:1:3d/3d::Particle index 7 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.872167 0.375099 0.884627
+DEAL:1:3d/3d::Particle index 8 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.138433 0.419489 0.692605
+DEAL:1:3d/3d::Particle index 9 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.474755 0.382099 0.844530
+DEAL:1:3d/3d::Particle index 10 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.319646 0.518143 0.716901
+DEAL:1:3d/3d::Particle index 11 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.456149 0.635687 0.609423
+DEAL:1:3d/3d::Particle index 12 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.588260 0.522915 0.510393
+DEAL:1:3d/3d::Particle index 13 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.256491 0.534039 0.682149
+DEAL:1:3d/3d::Particle index 14 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.957838 0.850786 0.290195
+DEAL:1:3d/3d::Particle index 15 is in cell 0.0
+DEAL:1:3d/3d::Particle location: 0.635976 0.937905 0.412123
+DEAL:1:3d/3d::OK
+