From 087d85546b94bdde19c00d93723da2e1afd55279 Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Tue, 24 Feb 2015 09:45:01 -0500 Subject: [PATCH] remove compress without args from all PETSc sparse matrix types - 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 | 67 ++++++++++------------ source/lac/petsc_matrix_base.cc | 65 ++++++++++----------- source/lac/petsc_parallel_sparse_matrix.cc | 5 +- source/lac/petsc_sparse_matrix.cc | 2 +- source/numerics/matrix_tools.cc | 4 -- 5 files changed, 65 insertions(+), 78 deletions(-) diff --git a/include/deal.II/lac/petsc_matrix_base.h b/include/deal.II/lac/petsc_matrix_base.h index bce1aa14f6..8a22c13cf7 100644 --- a/include/deal.II/lac/petsc_matrix_base.h +++ b/include/deal.II/lac/petsc_matrix_base.h @@ -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,j). 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 diff --git a/source/lac/petsc_matrix_base.cc b/source/lac/petsc_matrix_base.cc index 7e97c93f8b..ae82ff1394 100644 --- a/source/lac/petsc_matrix_base.cc +++ b/source/lac/petsc_matrix_base.cc @@ -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 &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, diff --git a/source/lac/petsc_parallel_sparse_matrix.cc b/source/lac/petsc_parallel_sparse_matrix.cc index ae8ae2a9a8..0232e5a143 100644 --- a/source/lac/petsc_parallel_sparse_matrix.cc +++ b/source/lac/petsc_parallel_sparse_matrix.cc @@ -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 diff --git a/source/lac/petsc_sparse_matrix.cc b/source/lac/petsc_sparse_matrix.cc index 4b3396a994..2835a27b46 100644 --- a/source/lac/petsc_sparse_matrix.cc +++ b/source/lac/petsc_sparse_matrix.cc @@ -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 diff --git a/source/numerics/matrix_tools.cc b/source/numerics/matrix_tools.cc index 87ca717ef0..e90e551af9 100644 --- a/source/numerics/matrix_tools.cc +++ b/source/numerics/matrix_tools.cc @@ -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 (); } -- 2.39.5