From: David Wells Date: Sat, 12 Nov 2016 19:21:35 +0000 (-0500) Subject: Improve the output of PETSc exceptions. X-Git-Tag: v8.5.0-rc1~382^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=99dd20f7a32197d17d71ff326d7af7867bc0bf59;p=dealii.git Improve the output of PETSc exceptions. --- diff --git a/include/deal.II/lac/exceptions.h b/include/deal.II/lac/exceptions.h index e40044bacf..8bb77eedad 100644 --- a/include/deal.II/lac/exceptions.h +++ b/include/deal.II/lac/exceptions.h @@ -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 + * PetscErrorMessage 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 diff --git a/include/deal.II/lac/petsc_solver.h b/include/deal.II/lac/petsc_solver.h index 61adf24b45..3ae0852178 100644 --- a/include/deal.II/lac/petsc_solver.h +++ b/include/deal.II/lac/petsc_solver.h @@ -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: diff --git a/source/lac/CMakeLists.txt b/source/lac/CMakeLists.txt index ca536009dc..5161250dee 100644 --- a/source/lac/CMakeLists.txt +++ b/source/lac/CMakeLists.txt @@ -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 index 0000000000..77c61a88df --- /dev/null +++ b/source/lac/exceptions.cc @@ -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 + +#include + +#ifdef DEAL_II_WITH_PETSC +#include +#include +#include +#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 index 0000000000..3dd689108c --- /dev/null +++ b/tests/petsc/exception_messages.cc @@ -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 + +#include +#include +#include + +#include +#include + +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 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 index 0000000000..ac40a7b1f8 --- /dev/null +++ b/tests/petsc/exception_messages.mpirun=1.output @@ -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.