]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Wrap multiple versions of PETSc's MatDestroy.
authorDavid Wells <wellsd2@rpi.edu>
Sat, 9 Jul 2016 15:09:22 +0000 (11:09 -0400)
committerDavid Wells <wellsd2@rpi.edu>
Tue, 12 Jul 2016 00:08:58 +0000 (20:08 -0400)
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.

include/deal.II/lac/petsc_compatibility.h
source/lac/petsc_full_matrix.cc
source/lac/petsc_matrix_base.cc
source/lac/petsc_matrix_free.cc
source/lac/petsc_parallel_sparse_matrix.cc
source/lac/petsc_sparse_matrix.cc

index 939b9a01c872e2f10476c0b98f9ee218c3556faf..ae5e02945f09635ebe1e18ac2a42f4920fda78ed 100644 (file)
@@ -25,6 +25,7 @@
 #ifdef DEAL_II_WITH_PETSC
 
 #include <petscconf.h>
+#include <petscmat.h>
 #include <petscpc.h>
 
 #include <string>
@@ -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
index 46e8b05ddea8ce40f4a906e37a024b33a2ef570d..9d2201c157a1cfbac945ed3d45f270c27d664cbd 100644 (file)
@@ -17,6 +17,7 @@
 
 #ifdef DEAL_II_WITH_PETSC
 
+#include <deal.II/lac/petsc_compatibility.h>
 #include <deal.II/lac/exceptions.h>
 
 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);
   }
index 0d412a84ad551f3f042434b5b33eb119d59622fc..80cbcb4196bdfffdb93472988acfe5a35cd9f05b 100644 (file)
@@ -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));
   }
 
index f3027673df76885ca6dbce6c76016b84356fc5ef..40e4034839dbccb7e204b44455a87b512c67e665 100644 (file)
@@ -19,6 +19,7 @@
 #ifdef DEAL_II_WITH_PETSC
 
 #include <deal.II/lac/exceptions.h>
+#include <deal.II/lac/petsc_compatibility.h>
 
 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);
index 686346a918d9edd017493f5a541614a4113a1b83..5b706537d591b1aa7d3007f06157a8366d732227 100644 (file)
@@ -19,6 +19,7 @@
 
 #  include <deal.II/base/mpi.h>
 #  include <deal.II/lac/exceptions.h>
+#  include <deal.II/lac/petsc_compatibility.h>
 #  include <deal.II/lac/petsc_vector.h>
 #  include <deal.II/lac/sparsity_pattern.h>
 #  include <deal.II/lac/dynamic_sparsity_pattern.h>
@@ -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);
     }
index 9132bf123de136e978c0c11534cf494bacf5682d..5e4e12e1bc163c0f746304d97252431a233043ae 100644 (file)
@@ -18,6 +18,7 @@
 #ifdef DEAL_II_WITH_PETSC
 
 #  include <deal.II/lac/exceptions.h>
+#  include <deal.II/lac/petsc_compatibility.h>
 #  include <deal.II/lac/petsc_vector.h>
 #  include <deal.II/lac/sparsity_pattern.h>
 #  include <deal.II/lac/dynamic_sparsity_pattern.h>
@@ -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);
   }

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.