]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Improve the output of PETSc exceptions.
authorDavid Wells <wellsd2@rpi.edu>
Sat, 12 Nov 2016 19:21:35 +0000 (14:21 -0500)
committerDavid Wells <wellsd2@rpi.edu>
Thu, 17 Nov 2016 15:35:04 +0000 (10:35 -0500)
include/deal.II/lac/exceptions.h
include/deal.II/lac/petsc_solver.h
source/lac/CMakeLists.txt
source/lac/exceptions.cc [new file with mode: 0644]
tests/petsc/exception_messages.cc [new file with mode: 0644]
tests/petsc/exception_messages.mpirun=1.output [new file with mode: 0644]

index e40044bacfd0bbf54ca040024498a8a285682c29..8bb77eedad79d36b7f2d7ce3777d3052f8c97a4c 100644 (file)
@@ -43,13 +43,22 @@ namespace LACExceptions
   DeclException0 (ExcDifferentBlockIndices);
 
   /**
-   * An error of a PETSc function was encountered. Check the PETSc
-   * documentation for details.
+   * Exception thrown when a PETSc function reports an error. If possible,
+   * this exception uses the message provided by
+   * <code>PetscErrorMessage</code> to print a description of the error.
+   *
+   * @note For backwards compatibility this is defined whether or not deal.II
+   * is compiled with PETSc.
    */
-  DeclException1 (ExcPETScError,
-                  int,
-                  << "An error with error number " << arg1
-                  << " occurred while calling a PETSc function");
+  class ExcPETScError : public dealii::ExceptionBase
+  {
+  public:
+    ExcPETScError (const int error_code);
+
+    virtual void print_info (std::ostream &out) const;
+
+    const int error_code;
+  };
 
   /**
    * An error of a Trilinos function was encountered. Check the Trilinos
index 61adf24b4560d0e138ec16940109434e7620345a..3ae08521784b6288e907ac355c5b765cf536c47e 100644 (file)
@@ -163,12 +163,15 @@ namespace PETScWrappers
     void initialize(const PreconditionerBase &preconditioner);
 
     /**
-     * Exception
+     * Exception.
+     *
+     * @deprecated This function has been deprecated in favor of the more
+     * general LACExceptions::ExcPETScError exception class.
      */
     DeclException1 (ExcPETScError,
                     int,
                     << "An error with error number " << arg1
-                    << " occurred while calling a PETSc function");
+                    << " occurred while calling a PETSc function") DEAL_II_DEPRECATED;
 
   protected:
 
index ca536009dc84f17faffb222fbbbe2ec906aa17de..5161250deedb7a6951555c91c995bc01c9926abe 100644 (file)
@@ -23,8 +23,9 @@ SET(_src
   block_vector.cc
   chunk_sparse_matrix.cc
   chunk_sparsity_pattern.cc
-  dynamic_sparsity_pattern.cc
   constraint_matrix.cc
+  dynamic_sparsity_pattern.cc
+  exceptions.cc
   full_matrix.cc
   lapack_full_matrix.cc
   la_vector.cc
diff --git a/source/lac/exceptions.cc b/source/lac/exceptions.cc
new file mode 100644 (file)
index 0000000..77c61a8
--- /dev/null
@@ -0,0 +1,64 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/lac/exceptions.h>
+
+#ifdef DEAL_II_WITH_PETSC
+#include <petscconf.h>
+#include <petscsys.h>
+#include <petscerror.h>
+#endif // DEAL_II_WITH_PETSC
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace LACExceptions
+{
+  ExcPETScError::ExcPETScError (const int error_code)
+    :
+    error_code (error_code)
+  {}
+
+  void ExcPETScError::print_info (std::ostream &out) const
+  {
+    out << "deal.II encountered an error while calling a PETSc function."
+        << std::endl;
+#ifdef DEAL_II_WITH_PETSC
+    // PetscErrorMessage changes the value in a pointer to refer to a
+    // statically allocated description of the current error message.
+    const char *petsc_message;
+    const int ierr = PetscErrorMessage (error_code, &petsc_message, /*specific=*/NULL);
+    if (ierr == 0 && petsc_message != NULL)
+      {
+        out << "The description of the error provided by PETSc is \""
+            << petsc_message
+            << "\"."
+            << std::endl;
+      }
+    else
+      {
+        out << "PETSc was not able to determine a description for this particular error code."
+            << std::endl;
+      }
+#endif // DEAL_II_WITH_PETSC
+    out << "The numerical value of the original error code is "
+        << error_code
+        << "."
+        << std::endl;
+  }
+} // namespace LACExceptions
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/tests/petsc/exception_messages.cc b/tests/petsc/exception_messages.cc
new file mode 100644 (file)
index 0000000..3dd6891
--- /dev/null
@@ -0,0 +1,59 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+// Check that we can catch some standard PETSc error codes and print
+// appropriate exception messages to the screen.
+
+#include "../tests.h"
+#include <deal.II/lac/exceptions.h>
+
+#include <petscconf.h>
+#include <petscsys.h>
+#include <petscerror.h>
+
+#include <fstream>
+#include <vector>
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  std::ofstream logfile("output");
+  deallog.attach(logfile);
+  deallog.threshold_double(1.e-10);
+
+  std::vector<int> petsc_error_codes;
+  petsc_error_codes.push_back(PETSC_ERR_MIN_VALUE);
+  petsc_error_codes.push_back(PETSC_ERR_MEM);
+  petsc_error_codes.push_back(PETSC_ERR_LIB);
+  petsc_error_codes.push_back(PETSC_ERR_ARG_WRONGSTATE);
+  petsc_error_codes.push_back(PETSC_ERR_ARG_NOTSAMETYPE);
+
+  for (unsigned int i = 0; i < petsc_error_codes.size(); ++i)
+    {
+      try
+        {
+          throw ExcPETScError(petsc_error_codes[i]);
+        }
+      catch (const ExceptionBase &exc)
+        {
+          deallog << exc.get_exc_name() << std::endl;
+          exc.print_info(deallog.get_file_stream());
+        }
+    }
+
+  return 0;
+}
diff --git a/tests/petsc/exception_messages.mpirun=1.output b/tests/petsc/exception_messages.mpirun=1.output
new file mode 100644 (file)
index 0000000..ac40a7b
--- /dev/null
@@ -0,0 +1,21 @@
+
+DEAL::
+deal.II encountered an error while calling a PETSc function.
+PETSc was not able to determine a description for this particular error code.
+The numerical value of the original error code is 54.
+DEAL::
+deal.II encountered an error while calling a PETSc function.
+The description of the error provided by PETSc is "Out of memory".
+The numerical value of the original error code is 55.
+DEAL::
+deal.II encountered an error while calling a PETSc function.
+The description of the error provided by PETSc is "Error in external library".
+The numerical value of the original error code is 76.
+DEAL::
+deal.II encountered an error while calling a PETSc function.
+The description of the error provided by PETSc is "Object is in wrong state".
+The numerical value of the original error code is 73.
+DEAL::
+deal.II encountered an error while calling a PETSc function.
+The description of the error provided by PETSc is "Arguments must have same type".
+The numerical value of the original error code is 69.

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.