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
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:
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
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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;
+}
--- /dev/null
+
+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.