]> https://gitweb.dealii.org/ - dealii.git/commitdiff
remove compress without args from all PETSc sparse matrix types
authorTimo Heister <timo.heister@gmail.com>
Tue, 24 Feb 2015 14:45:01 +0000 (09:45 -0500)
committerTimo Heister <timo.heister@gmail.com>
Tue, 24 Feb 2015 22:27:24 +0000 (17:27 -0500)
- compress() without arguments was deprecated and now got remove
- replace redundant LastAction enum in by VectorOperation
- add additional compress() check.

include/deal.II/lac/petsc_matrix_base.h
source/lac/petsc_matrix_base.cc
source/lac/petsc_parallel_sparse_matrix.cc
source/lac/petsc_sparse_matrix.cc
source/numerics/matrix_tools.cc

index bce1aa14f64c12d68cd44a273bd728af32c77146..8a22c13cf7f0d32e75fb75cbd631f027d05418ad 100644 (file)
@@ -527,12 +527,7 @@ namespace PETScWrappers
      * @ref GlossCompress "Compressing distributed objects"
      * for more information.
      */
-    void compress (::dealii::VectorOperation::values operation);
-
-    /**
-     * @deprecated: use compress() with VectorOperation instead.
-     */
-    void compress () DEAL_II_DEPRECATED;
+    void compress (const VectorOperation::values operation);
 
     /**
      * Return the value of the entry (<i>i,j</i>).  This may be an expensive
@@ -912,34 +907,24 @@ namespace PETScWrappers
      */
     Mat matrix;
 
-    /**
-     * PETSc doesn't allow to mix additions to matrix entries and overwriting
-     * them (to make synchronisation of parallel computations simpler). Since
-     * the interface of the existing classes don't support the notion of not
-     * interleaving things, we have to emulate this ourselves. The way we do
-     * it is to, for each access operation, store whether it is an insertion
-     * or an addition. If the previous one was of different type, then we
-     * first have to flush the PETSc buffers; otherwise, we can simply go on.
-     *
-     * The following structure and variable declare and store the previous
-     * state.
-     */
-    struct LastAction
-    {
-      enum Values { none, insert, add };
-    };
-
     /**
      * Store whether the last action was a write or add operation.
      */
-    LastAction::Values last_action;
+    VectorOperation::values last_action;
 
     /**
      * Ensure that the add/set mode that is required for actions following
      * this call is compatible with the current mode. Should be called from
      * all internal functions accessing matrix elements.
      */
-    void prepare_action(const LastAction::Values new_action);
+    void prepare_action(const VectorOperation::values new_action);
+
+    /**
+     * Internal function that checks that there are no pending insert/add
+     * operations. Throws an exception otherwise. Useful before calling
+     * any PETSc internal functions modifying the matrix.
+     */
+    void assert_is_compressed();
 
     /**
      * For some matrix storage formats, in particular for the PETSc
@@ -1241,7 +1226,7 @@ namespace PETScWrappers
                    const PetscScalar  *values,
                    const bool         elide_zero_values)
   {
-    prepare_action(LastAction::insert);
+    prepare_action(VectorOperation::insert);
 
     const PetscInt petsc_i = row;
     PetscInt *col_index_ptr;
@@ -1316,7 +1301,7 @@ namespace PETScWrappers
         // zero. However, these actions are done
         // in case we pass on to the other
         // function.
-        prepare_action(LastAction::add);
+        prepare_action(VectorOperation::add);
 
         return;
       }
@@ -1387,7 +1372,7 @@ namespace PETScWrappers
                    const bool         elide_zero_values,
                    const bool          /*col_indices_are_sorted*/)
   {
-    prepare_action(LastAction::add);
+    prepare_action(VectorOperation::add);
 
     const PetscInt petsc_i = row;
     PetscInt *col_index_ptr;
@@ -1520,14 +1505,24 @@ namespace PETScWrappers
 
   inline
   void
-  MatrixBase::prepare_action(const LastAction::Values new_action)
+  MatrixBase::prepare_action(const VectorOperation::values new_action)
   {
-    if (last_action == new_action)
-      ;
-    else if (last_action == LastAction::none)
+    if (last_action == VectorOperation::unknown)
       last_action = new_action;
-    else
-      Assert (false, ExcWrongMode (last_action, new_action));
+
+    Assert (last_action == new_action, ExcWrongMode (last_action, new_action));
+  }
+
+
+
+  inline
+  void
+  MatrixBase::assert_is_compressed ()
+  {
+    // compress() sets the last action to none, which allows us to check if there
+    // are pending add/insert operations:
+    AssertThrow (last_action == VectorOperation::unknown,
+                 ExcMessage("Error: missing compress() call."));
   }
 
 
@@ -1536,7 +1531,7 @@ namespace PETScWrappers
   void
   MatrixBase::prepare_add()
   {
-    prepare_action(LastAction::add);
+    prepare_action(VectorOperation::add);
   }
 
 
@@ -1545,7 +1540,7 @@ namespace PETScWrappers
   void
   MatrixBase::prepare_set()
   {
-    prepare_action(LastAction::insert);
+    prepare_action(VectorOperation::insert);
   }
 
 #endif // DOXYGEN
index 7e97c93f8b95db70e0629990d9212fe7210e2223..ae82ff13947b609a4b81f9fad17cb55afeb31c7c 100644 (file)
@@ -44,11 +44,7 @@ namespace PETScWrappers
           return;
         }
 
-      // otherwise first flush PETSc caches
-      matrix->compress ();
-
-      // get a representation of the present
-      // row
+      // get a representation of the present row
       PetscInt           ncols;
       const PetscInt    *colnums;
       const PetscScalar *values;
@@ -77,7 +73,7 @@ namespace PETScWrappers
 
   MatrixBase::MatrixBase ()
     :
-    last_action (LastAction::none)
+    last_action (VectorOperation::unknown)
   {}
 
 
@@ -119,11 +115,7 @@ namespace PETScWrappers
   {
     Assert (d==value_type(), ExcScalarAssignmentOnlyForZeroValue());
 
-    // flush previously cached elements. this
-    // seems to be necessary since petsc
-    // 2.2.1, at least for parallel vectors
-    // (see test bits/petsc_64)
-    compress ();
+    assert_is_compressed ();
 
     const int ierr = MatZeroEntries (matrix);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
@@ -137,10 +129,9 @@ namespace PETScWrappers
   MatrixBase::clear_row (const size_type   row,
                          const PetscScalar new_diag_value)
   {
-    compress ();
+    assert_is_compressed ();
 
-    // now set all the entries of this row to
-    // zero
+    // now set all the entries of this row to zero
     const PetscInt petsc_row = row;
 
     IS index_set;
@@ -164,8 +155,6 @@ namespace PETScWrappers
 #else
     ISDestroy (&index_set);
 #endif
-
-    compress ();
   }
 
 
@@ -174,7 +163,7 @@ namespace PETScWrappers
   MatrixBase::clear_rows (const std::vector<size_type> &rows,
                           const PetscScalar             new_diag_value)
   {
-    compress ();
+    assert_is_compressed ();
 
     // now set all the entries of these rows
     // to zero
@@ -207,8 +196,6 @@ namespace PETScWrappers
 #else
     ISDestroy (&index_set);
 #endif
-
-    compress ();
   }
 
 
@@ -244,8 +231,27 @@ namespace PETScWrappers
 
 
   void
-  MatrixBase::compress (::dealii::VectorOperation::values)
+  MatrixBase::compress (const VectorOperation::values operation)
   {
+#ifdef DEBUG
+#ifdef DEAL_II_WITH_MPI
+    // Check that all processors agree that last_action is the same (or none!)
+
+    int my_int_last_action = last_action;
+    int all_int_last_action;
+
+    MPI_Allreduce(&my_int_last_action, &all_int_last_action, 1, MPI_INT,
+                  MPI_BOR, get_mpi_communicator());
+
+    AssertThrow(all_int_last_action != (VectorOperation::add | VectorOperation::insert),
+                ExcMessage("Error: not all processors agree on the last VectorOperation before this compress() call."));
+#endif
+#endif
+
+    AssertThrow(last_action == VectorOperation::unknown
+                || last_action == operation,
+                ExcMessage("Missing compress() or calling with wrong VectorOperation argument."));
+
     // flush buffers
     int ierr;
     ierr = MatAssemblyBegin (matrix,MAT_FINAL_ASSEMBLY);
@@ -254,15 +260,7 @@ namespace PETScWrappers
     ierr = MatAssemblyEnd (matrix,MAT_FINAL_ASSEMBLY);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
 
-    last_action = LastAction::none;
-  }
-
-
-
-  void
-  MatrixBase::compress ()
-  {
-    compress(::dealii::VectorOperation::unknown);
+    last_action = VectorOperation::unknown;
   }
 
 
@@ -571,8 +569,7 @@ namespace PETScWrappers
     PetscBool
 #endif
     truth;
-    // First flush PETSc caches
-    compress ();
+    assert_is_compressed ();
     MatIsSymmetric (matrix, tolerance, &truth);
     return truth;
   }
@@ -591,8 +588,7 @@ namespace PETScWrappers
 #endif
     truth;
 
-    // First flush PETSc caches
-    compress ();
+    assert_is_compressed ();
     MatIsHermitian (matrix, tolerance, &truth);
 
     return truth;
@@ -601,8 +597,7 @@ namespace PETScWrappers
   void
   MatrixBase::write_ascii (const PetscViewerFormat format)
   {
-    // First flush PETSc caches
-    compress ();
+    assert_is_compressed ();
 
     // Set options
     PetscViewerSetFormat (PETSC_VIEWER_STDOUT_WORLD,
index ae8ae2a9a830e7248fa995b552b306403a6d464f..0232e5a14362faed39ffb39e18e96c959474137f 100644 (file)
@@ -567,6 +567,8 @@ namespace PETScWrappers
                                     local_columns_per_process.size()));
       Assert (this_process < local_rows_per_process.size(),
               ExcInternalError());
+      assert_is_compressed ();
+
       // for each row that we own locally, we
       // have to count how many of the
       // entries in the sparsity pattern lie
@@ -753,13 +755,12 @@ namespace PETScWrappers
               }
           }
 
-          compress ();
+          compress (VectorOperation::insert);
 
           // set the dummy entries set above
           // back to zero
           *this = 0;
 #endif // version <=2.3.3
-          compress ();
 
 
           // Tell PETSc that we are not
index 4b3396a994f5a9181f5a8109b56167c89a2f59b3..2835a27b46b67dd586602002f1830db47f1202d6 100644 (file)
@@ -262,7 +262,7 @@ namespace PETScWrappers
                           row_lengths[i], &row_entries[0],
                           &row_values[0], INSERT_VALUES);
           }
-        compress ();
+        compress (VectorOperation::insert);
 
 
         // Tell PETSc that we are not
index 87ca717ef0c760f64ea5b6e8a4ee42d0c7acbc6b..e90e551af920dbc76d8ae496b8029bb3c332400c 100644 (file)
@@ -2410,7 +2410,6 @@ namespace MatrixTools
           }
 
         // clean up
-        matrix.compress ();
         solution.compress (VectorOperation::insert);
         right_hand_side.compress (VectorOperation::insert);
       }
@@ -2445,9 +2444,6 @@ namespace MatrixTools
     // used for both petsc matrix types
     internal::PETScWrappers::apply_boundary_values (boundary_values, matrix, solution,
                                                     right_hand_side, eliminate_columns);
-
-    // compress the matrix once we're done
-    matrix.compress ();
   }
 
 

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.