]> https://gitweb.dealii.org/ - dealii.git/commitdiff
All tests passing
authorDaniel Arndt <d.arndt@math.uni-goettingen.de>
Sun, 11 Sep 2016 12:37:15 +0000 (14:37 +0200)
committerDaniel Arndt <d.arndt@math.uni-goettingen.de>
Mon, 12 Sep 2016 10:56:22 +0000 (12:56 +0200)
include/deal.II/base/index_set.h
source/base/index_set.cc
tests/trilinos/precondition_muelu_q_iso_q1.cc
tests/trilinos/precondition_muelu_q_iso_q1.with_trilinos.geq.11.14.with_64bit_indices=off.output

index 2ec7a302709786fe9383798272e3b0c16361a619..34400387a1057d6478e73fb3686ef703889efe8d 100644 (file)
@@ -219,13 +219,14 @@ public:
 
   /**
    * Return whether the IndexSets are ascending with respect to MPI process
-   * number, i.e., each index is contained in exactly one IndexSet,
-   * the first indices are contained in the IndexSet of the first MPI process,
-   * the second indices are contained in the IndexSet of the second MPI process
-   * and so on.
-   * In case there is only one MPI process, this is always true.
+   * number and 1:1, i.e., each index is contained in exactly one IndexSet
+   * (among those stored on the different processes), each process stores
+   * contiguous subset of indices, and the index set on process $p+1$ starts
+   * at the index one larger than the last one stored on process $p$.
+   * In case there is only one MPI process, this just means that the IndexSet
+   * is complete.
    */
-  bool is_globally_ascending() const;
+  bool is_ascending_and_one_to_one(const MPI_Comm &communicator) const;
 
   /**
    * Return the number of elements stored in this index set.
index 24d5c8ca9340886360391a08e8a5a58a2d4c6b0e..b174c0a33c2d0329d35ac1a9672ee4f8e85515b9 100644 (file)
@@ -512,11 +512,11 @@ IndexSet::make_trilinos_map (const MPI_Comm &communicator,
     }
 #endif
 
-  // Find out if the IndexSet is ascending.
-  // Overlapping IndexSets are never ascending.
-  const bool ascending = overlapping ? false : is_globally_ascending();
+  // Find out if the IndexSet is ascending and 1:1. This correspinds to a
+  // linear EpetraMap. Overlapping IndexSets are never 1:1.
+  const bool linear = overlapping ? false : is_ascending_and_one_to_one(communicator);
 
-  if (ascending)
+  if (linear)
     return Epetra_Map (TrilinosWrappers::types::int_type(size()),
                        TrilinosWrappers::types::int_type(n_elements()),
                        0,
@@ -551,36 +551,38 @@ IndexSet::make_trilinos_map (const MPI_Comm &communicator,
 
 
 bool
-IndexSet::is_globally_ascending () const
+IndexSet::is_ascending_and_one_to_one (const MPI_Comm &communicator) const
 {
-  // If there is only one process, the IndexSet is always ascending.
-  const unsigned int n_ranks = Utilities::MPI::n_mpi_processes(communicator);
-  if (n_ranks==1)
-    return true;
-  // / Non-contiguous IndexSets can't be linear.
+  // If the sum of local elements does not add up to the total size,
+  // the IndexSet can't be complete.
+  const size_type n_global_elements
+    = Utilities::MPI::sum (n_elements(), communicator);
+  if (n_global_elements != size())
+    return false;
+
+  // Non-contiguous IndexSets can't be linear.
   const bool all_contiguous = (Utilities::MPI::min (is_contiguous() ? 1 : 0, communicator) == 1);
   if (!all_contiguous)
     return false;
 
   bool is_globally_ascending = true;
   // we know that there is only one interval
-  types::global_dof_index local_dofs[2];
-  local_dofs[0] = (n_elements()>0) ? *(begin_intervals()->begin())
-                  : numbers::invalid_dof_index ;
-  local_dofs[1] = (n_elements()>0) ? begin_intervals()->last()
-                  : numbers::invalid_dof_index;
+  types::global_dof_index first_local_dof
+    = (n_elements()>0) ? *(begin_intervals()->begin())
+      : numbers::invalid_dof_index ;
 
   const unsigned int my_rank = Utilities::MPI::this_mpi_process(communicator);
+  const unsigned int n_ranks = Utilities::MPI::n_mpi_processes(communicator);
   // first gather all information on process 0
-  const unsigned int gather_size = (my_rank==0)?2*n_ranks:1;
+  const unsigned int gather_size = (my_rank==0)?n_ranks:1;
   std::vector<types::global_dof_index> global_dofs(gather_size);
 
-  MPI_Gather(local_dofs, 2, DEAL_II_DOF_INDEX_MPI_TYPE,
-             &(global_dofs[0]), 2, DEAL_II_DOF_INDEX_MPI_TYPE, 0,
+  MPI_Gather(&first_local_dof, 1, DEAL_II_DOF_INDEX_MPI_TYPE,
+             &(global_dofs[0]), 1, DEAL_II_DOF_INDEX_MPI_TYPE, 0,
              communicator);
   if (my_rank == 0)
     {
-      // find out if the received std::vector is linear
+      // find out if the received std::vector is ascending
       types::global_dof_index old_dof = global_dofs[0], new_dof = 0;
       types::global_dof_index index = 0;
       while (global_dofs[index] == numbers::invalid_dof_index)
@@ -591,7 +593,7 @@ IndexSet::is_globally_ascending () const
           new_dof = global_dofs[index];
           if (new_dof == numbers::invalid_dof_index)
             new_dof = old_dof;
-          else if (new_dof<old_dof)
+          else if (new_dof<=old_dof)
             {
               is_globally_ascending = false;
               break;
index de8d33e387e3e76cc2ddf9a0a03cf69dc209b8c9..3bb317ef2932e34f0b4da66bcefc66869e38efea 100644 (file)
@@ -276,7 +276,7 @@ void Step4<dim>::solve ()
     check_solver_within_range(
       solver.solve (system_matrix, solution, system_rhs,
                     preconditioner),
-      solver_control.last_step(), 25, 34);
+      solver_control.last_step(), 20, 34);
   }
   deallog.pop();
 
@@ -290,7 +290,7 @@ void Step4<dim>::solve ()
     check_solver_within_range(
       solver.solve (system_matrix, solution, system_rhs,
                     preconditioner),
-      solver_control.last_step(), 24, 40);
+      solver_control.last_step(), 20, 40);
   }
   deallog.pop();
   deallog << std::endl;
index ac0864f637890fc57387ffb245e936e84156d6b8..1e03d7a759d6bf3274c99ddab2db162611f4fe9c 100644 (file)
@@ -1,7 +1,7 @@
 
-DEAL:02401:MueLu_Q::Solver stopped within 25 - 34 iterations
-DEAL:02401:MueLu_Q_iso_Q1::Solver stopped within 24 - 40 iterations
+DEAL:02401:MueLu_Q::Solver stopped within 20 - 34 iterations
+DEAL:02401:MueLu_Q_iso_Q1::Solver stopped within 20 - 40 iterations
 DEAL:02401::
-DEAL:09409:MueLu_Q::Solver stopped within 25 - 34 iterations
-DEAL:09409:MueLu_Q_iso_Q1::Solver stopped within 24 - 40 iterations
+DEAL:09409:MueLu_Q::Solver stopped within 20 - 34 iterations
+DEAL:09409:MueLu_Q_iso_Q1::Solver stopped within 20 - 40 iterations
 DEAL:09409::

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.