From: David Wells Date: Sat, 9 Jul 2016 15:09:22 +0000 (-0400) Subject: Wrap multiple versions of PETSc's MatDestroy. X-Git-Tag: v8.5.0-rc1~872^2~5 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=f6b0d34e127eddddd8e095b2c3268ae7e8980dbe;p=dealii.git Wrap multiple versions of PETSc's MatDestroy. This commit puts the version checks of PETSc, which dictate the correct usage of MatDestroy, all in one place. The wrapped version does not throw exceptions, so this fixes a problem (throwing exceptions in a destructor if PETSc fails) with some of the current wrapper classes. --- diff --git a/include/deal.II/lac/petsc_compatibility.h b/include/deal.II/lac/petsc_compatibility.h index 939b9a01c8..ae5e02945f 100644 --- a/include/deal.II/lac/petsc_compatibility.h +++ b/include/deal.II/lac/petsc_compatibility.h @@ -25,6 +25,7 @@ #ifdef DEAL_II_WITH_PETSC #include +#include #include #include @@ -53,6 +54,28 @@ namespace PETScWrappers PetscOptionsSetValue (NULL, name.c_str (), value.c_str ()); #endif } + + + + /** + * Destroy a PETSc matrix. This function wraps MatDestroy with a version + * check (the signature of this function changed in PETSc 3.2.0). + * + * @warning Since the primary intent of this function is to enable RAII + * semantics in the PETSc wrappers, this function will not throw an + * exception if an error occurs, but instead just returns the error code + * given by MatDestroy. + * + */ + inline PetscErrorCode destroy_matrix (Mat &matrix) + { + // PETSc will check whether or not matrix is NULL. +#if DEAL_II_PETSC_VERSION_LT(3, 2, 0) + return MatDestroy (matrix); +#else + return MatDestroy (&matrix); +#endif + } } DEAL_II_NAMESPACE_CLOSE diff --git a/source/lac/petsc_full_matrix.cc b/source/lac/petsc_full_matrix.cc index 46e8b05dde..9d2201c157 100644 --- a/source/lac/petsc_full_matrix.cc +++ b/source/lac/petsc_full_matrix.cc @@ -17,6 +17,7 @@ #ifdef DEAL_II_WITH_PETSC +#include #include DEAL_II_NAMESPACE_OPEN @@ -42,12 +43,8 @@ namespace PETScWrappers { // get rid of old matrix and generate a // new one -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - const int ierr = MatDestroy (matrix); -#else - const int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + const PetscErrorCode ierr = destroy_matrix(matrix); + AssertThrow(ierr == 0, ExcPETScError(ierr)); do_reinit (m, n); } diff --git a/source/lac/petsc_matrix_base.cc b/source/lac/petsc_matrix_base.cc index 0d412a84ad..80cbcb4196 100644 --- a/source/lac/petsc_matrix_base.cc +++ b/source/lac/petsc_matrix_base.cc @@ -83,12 +83,7 @@ namespace PETScWrappers MatrixBase::~MatrixBase () { -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - const int ierr = MatDestroy (matrix); -#else - const int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + destroy_matrix (matrix); } @@ -97,17 +92,16 @@ namespace PETScWrappers MatrixBase::clear () { // destroy the matrix... -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - int ierr = MatDestroy (matrix); -#else - int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + { + const PetscErrorCode ierr = destroy_matrix (matrix); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + } + // ...and replace it by an empty // sequential matrix const int m=0, n=0, n_nonzero_per_row=0; - ierr = MatCreateSeqAIJ(PETSC_COMM_SELF, m, n, n_nonzero_per_row, - 0, &matrix); + const int ierr = MatCreateSeqAIJ(PETSC_COMM_SELF, m, n, n_nonzero_per_row, + 0, &matrix); AssertThrow (ierr == 0, ExcPETScError(ierr)); } diff --git a/source/lac/petsc_matrix_free.cc b/source/lac/petsc_matrix_free.cc index f3027673df..40e4034839 100644 --- a/source/lac/petsc_matrix_free.cc +++ b/source/lac/petsc_matrix_free.cc @@ -19,6 +19,7 @@ #ifdef DEAL_II_WITH_PETSC #include +#include DEAL_II_NAMESPACE_OPEN @@ -105,14 +106,9 @@ namespace PETScWrappers { this->communicator = communicator; - // destroy the matrix and - // generate a new one -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - int ierr = MatDestroy (matrix); -#else - int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + // destroy the matrix and generate a new one + const PetscErrorCode ierr = destroy_matrix (matrix); + AssertThrow (ierr == 0, ExcPETScError (ierr)); do_reinit (m, n, local_rows, local_columns); } @@ -133,13 +129,8 @@ namespace PETScWrappers ExcInternalError()); this->communicator = communicator; - -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - int ierr = MatDestroy (matrix); -#else - int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + const PetscErrorCode ierr = destroy_matrix (matrix); + AssertThrow (ierr != 0, ExcPETScError (ierr)); do_reinit (m, n, local_rows_per_process[this_process], @@ -171,12 +162,8 @@ namespace PETScWrappers void MatrixFree::clear () { -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - int ierr = MatDestroy (matrix); -#else - int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + const PetscErrorCode ierr = destroy_matrix (matrix); + AssertThrow (ierr == 0, ExcPETScError (ierr)); const int m=0; do_reinit (m, m, m, m); diff --git a/source/lac/petsc_parallel_sparse_matrix.cc b/source/lac/petsc_parallel_sparse_matrix.cc index 686346a918..5b706537d5 100644 --- a/source/lac/petsc_parallel_sparse_matrix.cc +++ b/source/lac/petsc_parallel_sparse_matrix.cc @@ -19,6 +19,7 @@ # include # include +# include # include # include # include @@ -45,14 +46,7 @@ namespace PETScWrappers SparseMatrix::~SparseMatrix () { - int ierr; - -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - ierr = MatDestroy (matrix); -#else - ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + destroy_matrix (matrix); } SparseMatrix::SparseMatrix (const MPI_Comm &communicator, @@ -116,15 +110,10 @@ namespace PETScWrappers this->communicator = other.communicator; - int ierr; -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - ierr = MatDestroy (matrix); -#else - ierr = MatDestroy (&matrix); -#endif + PetscErrorCode ierr = destroy_matrix (matrix); AssertThrow (ierr == 0, ExcPETScError(ierr)); - ierr = MatDuplicate(other.matrix, MAT_DO_NOT_COPY_VALUES, &matrix); + ierr = MatDuplicate (other.matrix, MAT_DO_NOT_COPY_VALUES, &matrix); AssertThrow (ierr == 0, ExcPETScError(ierr)); } @@ -160,14 +149,9 @@ namespace PETScWrappers { this->communicator = communicator; - // get rid of old matrix and generate a - // new one -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - const int ierr = MatDestroy (matrix); -#else - const int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + // get rid of old matrix and generate a new one + const PetscErrorCode ierr = destroy_matrix (matrix); + AssertThrow (ierr == 0, ExcPETScError (ierr)); do_reinit (m, n, local_rows, local_columns, n_nonzero_per_row, is_symmetric, @@ -190,12 +174,8 @@ namespace PETScWrappers // get rid of old matrix and generate a // new one -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - const int ierr = MatDestroy (matrix); -#else - const int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + const PetscErrorCode ierr = destroy_matrix (matrix); + AssertThrow (ierr == 0, ExcPETScError (ierr)); do_reinit (m, n, local_rows, local_columns, row_lengths, is_symmetric, offdiag_row_lengths); @@ -215,14 +195,11 @@ namespace PETScWrappers { this->communicator = communicator; - // get rid of old matrix and generate a - // new one -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - const int ierr = MatDestroy (matrix); -#else - const int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + // get rid of old matrix and generate a new one + destroy_matrix (matrix); + const PetscErrorCode ierr = destroy_matrix (matrix); + AssertThrow (ierr == 0, ExcPETScError (ierr)); + do_reinit (sparsity_pattern, local_rows_per_process, local_columns_per_process, this_process, @@ -239,14 +216,9 @@ namespace PETScWrappers { this->communicator = communicator; - // get rid of old matrix and generate a - // new one -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - const int ierr = MatDestroy (matrix); -#else - const int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + // get rid of old matrix and generate a new one + const PetscErrorCode ierr = destroy_matrix (matrix); + AssertThrow(ierr == 0, ExcPETScError (ierr)); do_reinit (local_rows, local_columns, sparsity_pattern); } diff --git a/source/lac/petsc_sparse_matrix.cc b/source/lac/petsc_sparse_matrix.cc index 9132bf123d..5e4e12e1bc 100644 --- a/source/lac/petsc_sparse_matrix.cc +++ b/source/lac/petsc_sparse_matrix.cc @@ -18,6 +18,7 @@ #ifdef DEAL_II_WITH_PETSC # include +# include # include # include # include @@ -85,12 +86,8 @@ namespace PETScWrappers { // get rid of old matrix and generate a // new one -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - const int ierr = MatDestroy (matrix); -#else - const int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + const PetscErrorCode ierr = destroy_matrix (matrix); + AssertThrow (ierr == 0, ExcPETScError (ierr)); do_reinit (m, n, n_nonzero_per_row, is_symmetric); } @@ -105,12 +102,8 @@ namespace PETScWrappers { // get rid of old matrix and generate a // new one -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - const int ierr = MatDestroy (matrix); -#else - const int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + const PetscErrorCode ierr = destroy_matrix (matrix); + AssertThrow (ierr == 0, ExcPETScError (ierr)); do_reinit (m, n, row_lengths, is_symmetric); } @@ -125,12 +118,8 @@ namespace PETScWrappers { // get rid of old matrix and generate a // new one -#if DEAL_II_PETSC_VERSION_LT(3,2,0) - const int ierr = MatDestroy (matrix); -#else - const int ierr = MatDestroy (&matrix); -#endif - AssertThrow (ierr == 0, ExcPETScError(ierr)); + const PetscErrorCode ierr = destroy_matrix (matrix); + AssertThrow (ierr == 0, ExcPETScError (ierr)); do_reinit (sparsity_pattern, preset_nonzero_locations); }