]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Extend model for compatible partitioners.
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Sat, 12 Sep 2015 19:05:31 +0000 (21:05 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 22 Sep 2015 14:56:47 +0000 (16:56 +0200)
Fix bug in operator=.

include/deal.II/base/partitioner.h
include/deal.II/lac/parallel_vector.h
include/deal.II/lac/parallel_vector.templates.h
source/base/partitioner.cc

index 25239333d85ed10f5b182a0bc22afb8b69d98c26..5d856e5e034be3208ef4e83853f8cd560fe3ca6a 100644 (file)
@@ -212,7 +212,7 @@ namespace Utilities
       /**
        * Checks whether the given partitioner is compatible with the
        * partitioner used for this vector. Two partitioners are compatible if
-       * the have the same local size and the same ghost indices. They do not
+       * they have the same local size and the same ghost indices. They do not
        * necessarily need to be the same data field. This is a local operation
        * only, i.e., if only some processors decide that the partitioning is
        * not compatible, only these processors will return @p false, whereas
@@ -220,6 +220,21 @@ namespace Utilities
        */
       bool is_compatible (const Partitioner &part) const;
 
+      /**
+       * Checks whether the given partitioner is compatible with the
+       * partitioner used for this vector. Two partitioners are compatible if
+       * they have the same local size and the same ghost indices. They do not
+       * necessarily need to be the same data field. As opposed to
+       * is_compatible(), this method checks for compatibility among all
+       * processors and the method only returns @p true if the partitioner is
+       * the same on all processors.
+       *
+       * This method performs global communication, so make sure to use it
+       * only in a context where all processors call it the same number of
+       * times.
+       */
+      bool is_globally_compatible (const Partitioner &part) const;
+
       /**
        * Returns the MPI ID of the calling processor. Cached to have simple
        * access.
@@ -485,21 +500,6 @@ namespace Utilities
 
 
 
-    inline
-    bool
-    Partitioner::is_compatible (const Partitioner &part) const
-    {
-      // is the partitioner points to the same memory location as the calling
-      // processor
-      if (&part == this)
-        return true;
-      else
-        return (global_size == part.global_size &&
-                local_range_data == part.local_range_data &&
-                ghost_indices_data == part.ghost_indices_data);
-    }
-
-
     inline
     unsigned int
     Partitioner::this_mpi_process() const
index 5398967de3d0afb0a44bc8f393a91f14b877e491..7c9b50aa7e4cffa822a12aa1024125e25cf3254c 100644 (file)
@@ -962,7 +962,7 @@ namespace parallel
       /**
        * Checks whether the given partitioner is compatible with the
        * partitioner used for this vector. Two partitioners are compatible if
-       * the have the same local size and the same ghost indices. They do not
+       * they have the same local size and the same ghost indices. They do not
        * necessarily need to be the same data field. This is a local operation
        * only, i.e., if only some processors decide that the partitioning is
        * not compatible, only these processors will return @p false, whereas
@@ -971,6 +971,21 @@ namespace parallel
       bool
       partitioners_are_compatible (const Utilities::MPI::Partitioner &part) const;
 
+      /**
+       * Checks whether the given partitioner is compatible with the
+       * partitioner used for this vector. Two partitioners are compatible if
+       * they have the same local size and the same ghost indices. They do not
+       * necessarily need to be the same data field. As opposed to
+       * partitioners_are_compatible(), this method checks for compatibility
+       * among all processors and the method only returns @p true if the
+       * partitioner is the same on all processors.
+       *
+       * This method performs global communication, so make sure to use it
+       * only in a context where all processors call it the same number of
+       * times.
+       */
+      bool
+      partitioners_are_globally_compatible (const Utilities::MPI::Partitioner &part) const;
 
       /**
        * Prints the vector to the output stream @p out.
@@ -1309,13 +1324,18 @@ namespace parallel
         reinit (c, true);
       else if (partitioner.get() != c.partitioner.get())
         {
-          size_type local_ranges_different_loc = (local_range() !=
-                                                  c.local_range());
-          if ((partitioner->n_mpi_processes() > 1 &&
-               Utilities::MPI::max(local_ranges_different_loc,
-                                   partitioner->get_communicator()) != 0)
+          // local ranges are also the same if both partitioners are empty
+          // (even if they happen to define the empty range as [0,0) or [c,c)
+          // for some c!=0 in a different way).
+          int local_ranges_are_identical =
+            (local_range() == c.local_range() ||
+             (local_range().second == local_range().first &&
+              c.local_range().second == c.local_range().first));
+          if ((c.partitioner->n_mpi_processes() > 1 &&
+               Utilities::MPI::min(local_ranges_are_identical,
+                                   c.partitioner->get_communicator()) == 0)
               ||
-              local_ranges_different_loc)
+              !local_ranges_are_identical)
             reinit (c, true);
           else
             must_update_ghost_values |= vector_is_ghosted;
@@ -1334,6 +1354,8 @@ namespace parallel
 
       if (must_update_ghost_values)
         update_ghost_values();
+      else
+        zero_out_ghosts();
       return *this;
     }
 
@@ -2341,17 +2363,6 @@ namespace parallel
       return partitioner->get_communicator();
     }
 
-
-
-    template <typename Number>
-    inline
-    bool
-    Vector<Number>::partitioners_are_compatible
-    (const Utilities::MPI::Partitioner &part) const
-    {
-      return partitioner->is_compatible (part);
-    }
-
 #endif  // ifndef DOXYGEN
 
   } // end of namespace distributed
index 89b280dfc9cf2ce895f8fb230a45eb13d6a013e5..cce0c112a3544575540cadc4799b1f46efe2bdf2 100644 (file)
@@ -673,6 +673,28 @@ namespace parallel
 
 
 
+    template <typename Number>
+    inline
+    bool
+    Vector<Number>::partitioners_are_compatible
+    (const Utilities::MPI::Partitioner &part) const
+    {
+      return partitioner->is_compatible (part);
+    }
+
+
+
+    template <typename Number>
+    inline
+    bool
+    Vector<Number>::partitioners_are_globally_compatible
+    (const Utilities::MPI::Partitioner &part) const
+    {
+      return partitioner->is_globally_compatible (part);
+    }
+
+
+
     template <typename Number>
     std::size_t
     Vector<Number>::memory_consumption () const
index 1b216a9d4ef1de6eafe59c98ef30b411513a0293..b39efd0207a52f32844e1026c99c6fbcd631c50f 100644 (file)
@@ -325,6 +325,40 @@ namespace Utilities
 
 
 
+    bool
+    Partitioner::is_compatible (const Partitioner &part) const
+    {
+      // if the partitioner points to the same memory location as the calling
+      // processor
+      if (&part == this)
+        return true;
+#ifdef DEAL_II_WITH_MPI
+      if (Utilities::MPI::job_supports_mpi())
+        {
+          int communicators_same = 0;
+          MPI_Comm_compare (part.communicator, communicator,
+                            &communicators_same);
+          if (!(communicators_same == MPI_IDENT ||
+                communicators_same == MPI_CONGRUENT))
+            return false;
+        }
+#endif
+      return (global_size == part.global_size &&
+              local_range_data == part.local_range_data &&
+              ghost_indices_data == part.ghost_indices_data);
+    }
+
+
+
+    bool
+    Partitioner::is_globally_compatible (const Partitioner &part) const
+    {
+      return Utilities::MPI::min(static_cast<int>(is_compatible(part)),
+                                 communicator) == 1;
+    }
+
+
+
     std::size_t
     Partitioner::memory_consumption() const
     {

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.