]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add tests.
authorWolfgang Bangerth <bangerth@colostate.edu>
Tue, 6 Apr 2021 21:48:36 +0000 (15:48 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 22 Apr 2021 13:41:24 +0000 (07:41 -0600)
include/deal.II/base/aligned_vector.h
tests/base/aligned_vector_replicate_01.cc [new file with mode: 0644]
tests/base/aligned_vector_replicate_01.mpirun=3.output [new file with mode: 0644]
tests/base/aligned_vector_replicate_02.cc [new file with mode: 0644]
tests/base/aligned_vector_replicate_02.mpirun=3.output [new file with mode: 0644]
tests/base/aligned_vector_replicate_04.cc [new file with mode: 0644]
tests/base/aligned_vector_replicate_04.mpirun=3.output [new file with mode: 0644]

index 46c5478e51650f25ac01514351ec123778fb3670..73c60c042851b54c194f8dd16f71c583f9dfd89f 100644 (file)
@@ -325,6 +325,23 @@ public:
    *   typically in memory areas not accessible to the other processes.
    *   As a consequence, the usual use case for this function is to share
    *   arrays of simple objects such as `double`s or `int`s.
+   *
+   * @note After calling this function, objects on different MPI processes
+   *   share a common state. That means that certain operations become
+   *   "collective", i.e., they must be called on all participating
+   *   processors at the same time. In particular, you can no longer call
+   *   resize(), reserve(), or clear() on one MPI process -- you have to do
+   *   so on all processes at the same time, because they have to communicate
+   *   for these operations. If you do not do so, you will likely get
+   *   a deadlock that may be difficult to debug. By extension, this rule of
+   *   only collectively resizing extends to this function itself: You can
+   *   not call it twice in a row because that implies that first all but the
+   *   `root_process` throw away their data, which is not a collective
+   *   operation. Generally, these restrictions on what can and can not be
+   *   done hint at the correctness of the comments above: You should treat
+   *   an AlignedVector on which the current function has been called as
+   *   `const`, on which no further operations can be performed until
+   *   the destructor is called.
    */
   void
   replicate_across_communicator(const MPI_Comm &   communicator,
@@ -1114,6 +1131,17 @@ AlignedVector<T>::replicate_across_communicator(const MPI_Comm &   communicator,
 {
 #  ifdef DEAL_II_WITH_MPI
 #    if DEAL_II_MPI_VERSION_GTE(3, 0)
+
+  // **** Step 0 ****
+  // All but the root process no longer need their data, so release the memory
+  // used to store the previous elements.
+  if (Utilities::MPI::this_mpi_process(communicator) != root_process)
+    {
+      elements.reset();
+      used_elements_end      = nullptr;
+      allocated_elements_end = nullptr;
+    }
+
   // **** Step 1 ****
   // Create communicators for each group of processes that can use
   // shared memory areas. Within each of these groups, we don't care about
diff --git a/tests/base/aligned_vector_replicate_01.cc b/tests/base/aligned_vector_replicate_01.cc
new file mode 100644 (file)
index 0000000..7b64858
--- /dev/null
@@ -0,0 +1,68 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2012 - 2018 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.
+//
+// ---------------------------------------------------------------------
+
+
+// Test AlignedVector::replicate_across_communicator()
+
+#include <deal.II/base/aligned_vector.h>
+
+#include "../tests.h"
+
+
+void
+test()
+{
+  const MPI_Comm     communicator = MPI_COMM_WORLD;
+  const unsigned int root         = 1;
+  Assert(root < Utilities::MPI::n_mpi_processes(communicator),
+         ExcInternalError());
+
+  // Create a vector of differing sizes on all processes, but only on
+  // the root process do we put something useful into it
+  AlignedVector<int> avec(
+    Utilities::MPI::this_mpi_process(communicator) == root ? 10 : 5);
+  if (Utilities::MPI::this_mpi_process(communicator) == root)
+    {
+      avec[2] = 2;
+      avec[4] = 4;
+      avec[6] = 6;
+    }
+  else
+    {
+      for (unsigned int i = 0; i < avec.size(); ++i)
+        avec[i] = 100 * Utilities::MPI::this_mpi_process(communicator) + i;
+    }
+
+  // Now replicate this information across all processes
+  avec.replicate_across_communicator(communicator, root);
+
+  // Final step, let every process output what it now has
+  deallog << "On process " << Utilities::MPI::this_mpi_process(communicator)
+          << ": " << std::endl;
+  for (const auto i : avec)
+    deallog << i << ' ';
+  deallog << std::endl;
+}
+
+
+
+int
+main(int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll                    all;
+
+  test();
+}
diff --git a/tests/base/aligned_vector_replicate_01.mpirun=3.output b/tests/base/aligned_vector_replicate_01.mpirun=3.output
new file mode 100644 (file)
index 0000000..513345c
--- /dev/null
@@ -0,0 +1,11 @@
+
+DEAL:0::On process 0: 
+DEAL:0::0 0 2 0 4 0 6 0 0 0 
+
+DEAL:1::On process 1: 
+DEAL:1::0 0 2 0 4 0 6 0 0 0 
+
+
+DEAL:2::On process 2: 
+DEAL:2::0 0 2 0 4 0 6 0 0 0 
+
diff --git a/tests/base/aligned_vector_replicate_02.cc b/tests/base/aligned_vector_replicate_02.cc
new file mode 100644 (file)
index 0000000..ce7dc8b
--- /dev/null
@@ -0,0 +1,84 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2012 - 2018 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.
+//
+// ---------------------------------------------------------------------
+
+
+// Test AlignedVector::replicate_across_communicator(). This function
+// puts data into a shared memory area for all MPI processes that are
+// on the same machine. This is a read-write area, and so writing into
+// it becomes visible also to those processes on the same
+// machine. Since we run the test suite only in setups where *all*
+// processes are on the same machine, we know that a change made on
+// one process should be visible on *all*.
+
+#include <deal.II/base/aligned_vector.h>
+
+#include "../tests.h"
+
+
+void
+test()
+{
+  const MPI_Comm     communicator = MPI_COMM_WORLD;
+  const unsigned int root         = 1;
+  Assert(root < Utilities::MPI::n_mpi_processes(communicator),
+         ExcInternalError());
+
+  // Same as for the _01 test:
+  AlignedVector<int> avec(
+    Utilities::MPI::this_mpi_process(communicator) == root ? 10 : 5);
+  if (Utilities::MPI::this_mpi_process(communicator) == root)
+    {
+      avec[2] = 2;
+      avec[4] = 4;
+      avec[6] = 6;
+    }
+  else
+    {
+      for (unsigned int i = 0; i < avec.size(); ++i)
+        avec[i] = 100 * Utilities::MPI::this_mpi_process(communicator) + i;
+    }
+
+  // Now replicate this information across all processes
+  avec.replicate_across_communicator(communicator, root);
+
+
+  // Now let each process change one element:
+  if (Utilities::MPI::this_mpi_process(communicator) == 0)
+    avec[0] = 42;
+  else if (Utilities::MPI::this_mpi_process(communicator) == 1)
+    avec[1] = 43;
+  else if (Utilities::MPI::this_mpi_process(communicator) == 2)
+    avec[3] = 44;
+
+  // Final step, let every process output what it now has. Every
+  // process should be able to see the changed element, even though
+  // only one changed it.
+  deallog << "On process " << Utilities::MPI::this_mpi_process(communicator)
+          << ": " << std::endl;
+  for (const auto i : avec)
+    deallog << i << ' ';
+  deallog << std::endl;
+}
+
+
+
+int
+main(int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll                    all;
+
+  test();
+}
diff --git a/tests/base/aligned_vector_replicate_02.mpirun=3.output b/tests/base/aligned_vector_replicate_02.mpirun=3.output
new file mode 100644 (file)
index 0000000..688e26f
--- /dev/null
@@ -0,0 +1,11 @@
+
+DEAL:0::On process 0: 
+DEAL:0::42 43 2 44 4 0 6 0 0 0 
+
+DEAL:1::On process 1: 
+DEAL:1::42 43 2 44 4 0 6 0 0 0 
+
+
+DEAL:2::On process 2: 
+DEAL:2::42 43 2 44 4 0 6 0 0 0 
+
diff --git a/tests/base/aligned_vector_replicate_04.cc b/tests/base/aligned_vector_replicate_04.cc
new file mode 100644 (file)
index 0000000..0fac6cc
--- /dev/null
@@ -0,0 +1,51 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2012 - 2018 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.
+//
+// ---------------------------------------------------------------------
+
+
+// Test AlignedVector::replicate_across_communicator().
+//
+// Check what happens if we call this function on an empty object.
+
+#include <deal.II/base/aligned_vector.h>
+
+#include "../tests.h"
+
+
+void
+test()
+{
+  const MPI_Comm     communicator = MPI_COMM_WORLD;
+  const unsigned int root         = 1;
+  Assert(root < Utilities::MPI::n_mpi_processes(communicator),
+         ExcInternalError());
+
+  // Create an empty object and replicate it.
+  AlignedVector<int> avec;
+  avec.replicate_across_communicator(communicator, root);
+
+  deallog << "On process " << Utilities::MPI::this_mpi_process(communicator)
+          << ": " << avec.size() << std::endl;
+}
+
+
+
+int
+main(int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll                    all;
+
+  test();
+}
diff --git a/tests/base/aligned_vector_replicate_04.mpirun=3.output b/tests/base/aligned_vector_replicate_04.mpirun=3.output
new file mode 100644 (file)
index 0000000..e4fe2b2
--- /dev/null
@@ -0,0 +1,8 @@
+
+DEAL:0::On process 0: 0
+
+DEAL:1::On process 1: 0
+
+
+DEAL:2::On process 2: 0
+

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.