]> https://gitweb.dealii.org/ - dealii.git/commitdiff
improve MPI::Partitioner to support empty ArrayView for owned DoFs 6099/head
authorDenis Davydov <davydden@gmail.com>
Sat, 24 Mar 2018 19:26:25 +0000 (20:26 +0100)
committerDenis Davydov <davydden@gmail.com>
Sat, 24 Mar 2018 21:57:48 +0000 (22:57 +0100)
include/deal.II/base/partitioner.templates.h
tests/mpi/parallel_partitioner_08.cc [new file with mode: 0644]
tests/mpi/parallel_partitioner_08.mpirun=2.output [new file with mode: 0644]

index 5ba72ce87dfba922a5462d5f88c823318742d856..66d24d5f1f355f5abb2be70484709995b8c5f412 100644 (file)
@@ -42,7 +42,6 @@ namespace Utilities
                                                const ArrayView<Number>       &ghost_array,
                                                std::vector<MPI_Request>      &requests) const
     {
-      AssertDimension(locally_owned_array.size(), local_size());
       AssertDimension(temporary_storage.size(), n_import_indices());
       Assert(ghost_array.size() == n_ghost_indices() ||
              ghost_array.size() == n_ghost_indices_in_larger_set,
@@ -52,6 +51,9 @@ namespace Utilities
       const unsigned int n_import_targets = import_targets_data.size();
       const unsigned int n_ghost_targets = ghost_targets_data.size();
 
+      if (n_import_targets>0)
+        AssertDimension(locally_owned_array.size(), local_size());
+
       Assert(requests.size() == 0,
              ExcMessage("Another operation seems to still be running. "
                         "Call update_ghost_values_finish() first."));
@@ -316,7 +318,6 @@ namespace Utilities
                                                   const ArrayView<Number>       &ghost_array,
                                                   std::vector<MPI_Request>      &requests) const
     {
-      AssertDimension(locally_owned_array.size(), local_size());
       AssertDimension(temporary_storage.size(), n_import_indices());
       Assert(ghost_array.size() == n_ghost_indices() ||
              ghost_array.size() == n_ghost_indices_in_larger_set,
@@ -353,6 +354,7 @@ namespace Utilities
       // first wait for the receive to complete
       if (requests.size() > 0 && n_import_targets > 0)
         {
+          AssertDimension(locally_owned_array.size(), local_size());
           const int ierr = MPI_Waitall (n_import_targets, requests.data(),
                                         MPI_STATUSES_IGNORE);
           AssertThrowMPI(ierr);
diff --git a/tests/mpi/parallel_partitioner_08.cc b/tests/mpi/parallel_partitioner_08.cc
new file mode 100644 (file)
index 0000000..cd0a8b3
--- /dev/null
@@ -0,0 +1,172 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 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 at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// test MPI::Partitioner update_ghosts() and compress() in case we have
+// empty owned DoFs
+
+#include "../tests.h"
+#include <deal.II/base/mpi.h>
+#include <deal.II/base/index_set.h>
+#include <deal.II/base/partitioner.h>
+#include <deal.II/base/aligned_vector.h>
+
+template <typename Number=double>
+void test()
+{
+  const unsigned int rank = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+
+  // setup index sets
+  //                            subset:                    is2
+  //                            ghost:       8 9  10 11    is3
+  //      rank 0 :  00 01 02 03 04 05 06 07 00 00 00 00
+  //      rank 1 :  00 00 00 00 00 00 00 00 00 00 00 00
+  //                            ghost:      0  1  2  3     is3
+  //                            subset:        1  2        is2
+  //
+  // expected result udpate ghosts()
+  //
+  //      rank 0 :  00 01 02 03 04 05 06 07 00 00 00 00
+  //      rank 1 :  00 00 00 00 00 00 00 00 00 01 02 00
+  //
+  // compress(insert) -- does not change anything but zero ghosts
+  //
+  // set rank1 ghosts to: 00 10 20 00
+  // compress(add)
+  //
+  //      rank 0 :  00 11 22 03 04 05 06 07 00 00 00 00
+  //      rank 1 :  00 00 00 00 00 00 00 00 00 10 20 00
+
+
+  IndexSet is1(16), is2(16), is3(16);
+
+  if (rank==0)
+    {
+      is1.add_range(0,8);
+      // note: empty is2
+      is3.add_range(8,12);
+    }
+  else if (rank==1)
+    {
+      is1.add_range(8,16);
+      is2.add_index(1);
+      is2.add_index(2);
+      is3.add_range(0,4);
+    }
+
+  // create partitioner
+  std::shared_ptr<Utilities::MPI::Partitioner> partitioner
+  (new Utilities::MPI::Partitioner(is1, MPI_COMM_WORLD));
+  partitioner->set_ghost_indices(is3);
+  std::shared_ptr<Utilities::MPI::Partitioner> tight_partitioner
+  (new Utilities::MPI::Partitioner(is1, MPI_COMM_WORLD));
+  tight_partitioner->set_ghost_indices(is2,is3);
+
+  // create vector
+  AlignedVector<Number> owned(rank==0?8:0);
+  AlignedVector<Number> ghost(4);
+
+  for (unsigned int i = 0; i < 4; ++i)
+    ghost[i] = 0.;
+
+  if (rank==0)
+    for (int i = 0; i < 8; i++)
+      owned[i] = i;
+
+  // update ghost values
+  // vector of requests
+  std::vector<MPI_Request> requests;
+  std::vector<MPI_Request> compress_requests;
+
+  // allocate temporal array
+  AlignedVector<Number> tmp_data;
+  tmp_data.resize_fast(tight_partitioner->n_import_indices());
+
+  // begin exchange, and ...
+  tight_partitioner->export_to_ghosted_array_start(
+    0,
+    ArrayView<const Number>(owned.begin(), owned.size()),
+    ArrayView<Number>(tmp_data.begin(), tight_partitioner->n_import_indices()),
+    ArrayView<Number>(ghost.begin(), ghost.size()),
+    requests);
+
+  // ... finish exchange
+  tight_partitioner->export_to_ghosted_array_finish(
+    ArrayView<Number>(ghost.begin(), ghost.size()),
+    requests);
+
+  auto print = [&] ()
+  {
+    deallog << "owned:" << std::endl;
+    for (auto el : owned)
+      deallog << el << " ";
+    deallog << std::endl
+            << "ghost:" << std::endl;
+    for (auto el : ghost)
+      deallog << el << " ";
+    deallog << std::endl;
+  };
+
+  deallog << "update ghosts()"<<std::endl;
+  print();
+
+  AlignedVector<Number> import_data;
+  import_data.resize_fast(tight_partitioner->n_import_indices());
+
+  // now do insert:
+  auto compress = [&] (VectorOperation::values operation)
+  {
+    const unsigned int counter = 0;
+    tight_partitioner->import_from_ghosted_array_start
+    (operation, counter,
+     ArrayView<Number>(ghost.begin(), ghost.size()),
+     ArrayView<Number>(import_data.begin(), tight_partitioner->n_import_indices()),
+     compress_requests);
+
+    tight_partitioner->import_from_ghosted_array_finish
+    (operation,
+     ArrayView<const Number>(import_data.begin(), tight_partitioner->n_import_indices()),
+     ArrayView<Number>(owned.begin(), owned.size()),
+     ArrayView<Number>(ghost.begin(), ghost.size()),
+     compress_requests);
+  };
+
+  deallog << "compress(insert)"<<std::endl;
+  compress(VectorOperation::insert);
+  print();
+
+  if (rank==1)
+    {
+      ghost[1]=10;
+      ghost[2]=20;
+    }
+
+  deallog << "compress(add)"<<std::endl;
+  compress(VectorOperation::add);
+  print();
+
+}
+
+int main(int argc, char **argv)
+{
+  using namespace dealii;
+
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  MPILogInitAll log;
+
+  test();
+
+  return 0;
+}
diff --git a/tests/mpi/parallel_partitioner_08.mpirun=2.output b/tests/mpi/parallel_partitioner_08.mpirun=2.output
new file mode 100644 (file)
index 0000000..d60bc99
--- /dev/null
@@ -0,0 +1,33 @@
+
+DEAL:0::update ghosts()
+DEAL:0::owned:
+DEAL:0::0.00000 1.00000 2.00000 3.00000 4.00000 5.00000 6.00000 7.00000 
+DEAL:0::ghost:
+DEAL:0::0.00000 0.00000 0.00000 0.00000 
+DEAL:0::compress(insert)
+DEAL:0::owned:
+DEAL:0::0.00000 1.00000 2.00000 3.00000 4.00000 5.00000 6.00000 7.00000 
+DEAL:0::ghost:
+DEAL:0::0.00000 0.00000 0.00000 0.00000 
+DEAL:0::compress(add)
+DEAL:0::owned:
+DEAL:0::0.00000 11.0000 22.0000 3.00000 4.00000 5.00000 6.00000 7.00000 
+DEAL:0::ghost:
+DEAL:0::0.00000 0.00000 0.00000 0.00000 
+
+DEAL:1::update ghosts()
+DEAL:1::owned:
+DEAL:1::
+DEAL:1::ghost:
+DEAL:1::0.00000 1.00000 2.00000 0.00000 
+DEAL:1::compress(insert)
+DEAL:1::owned:
+DEAL:1::
+DEAL:1::ghost:
+DEAL:1::0.00000 0.00000 0.00000 0.00000 
+DEAL:1::compress(add)
+DEAL:1::owned:
+DEAL:1::
+DEAL:1::ghost:
+DEAL:1::0.00000 0.00000 0.00000 0.00000 
+

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.