From: David Wells Date: Thu, 7 Jul 2016 09:50:59 +0000 (-0400) Subject: Fix compilation with PETSc 3.7.0. X-Git-Tag: v8.5.0-rc1~918^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F2754%2Fhead;p=dealii.git Fix compilation with PETSc 3.7.0. --- diff --git a/doc/news/changes.h b/doc/news/changes.h index 6e393a1905..4a537c0720 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -150,6 +150,12 @@ inconvenience this causes.

General

    +
  1. New: The library is now compatible with PETSc 3.7.0. Part of this change + included adding a new header, petsc_compatibility.h, which provides + some version-independent functions for using common PETSc functions. +
    + (David Wells, 2016/07/07) +
  2. New: Added TrilinosWrappers::SolveDirect::Initialize and TrilinosWrappers::SolverDirect::Solve to solve distributed linear systems diff --git a/include/deal.II/lac/petsc_compatibility.h b/include/deal.II/lac/petsc_compatibility.h new file mode 100644 index 0000000000..705ffa01a0 --- /dev/null +++ b/include/deal.II/lac/petsc_compatibility.h @@ -0,0 +1,55 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + +/* + * Rather than using ifdefs everywhere, try to wrap older versions of PETSc + * functions in one place. + */ +#ifndef dealii__petsc_compatibility_h +#define dealii__petsc_compatibility_h + +#include + +#ifdef DEAL_II_WITH_PETSC + +#include +#include + +#include + +DEAL_II_NAMESPACE_OPEN + +namespace PETScWrappers +{ + /** + * Set an option in the global PETSc database. This function just wraps + * PetscOptionsSetValue with a version check (the signature of this function + * changed in PETSc 3.7.0). + */ + inline void set_option_value (const std::string &name, + const std::string &value) + { +#if DEAL_II_PETSC_VERSION_LT(3, 7, 0) + PetscOptionsSetValue (name.c_str (), value.c_str ()); +#else + PetscOptionsSetValue (NULL, name.c_str (), value.c_str ()); +#endif + } +} + +DEAL_II_NAMESPACE_CLOSE + +#endif // DEAL_II_WITH_PETSC +#endif // dealii__petsc_compatibility_h diff --git a/source/lac/petsc_precondition.cc b/source/lac/petsc_precondition.cc index 7dd465f16e..96bd524b52 100644 --- a/source/lac/petsc_precondition.cc +++ b/source/lac/petsc_precondition.cc @@ -18,6 +18,7 @@ #ifdef DEAL_II_WITH_PETSC # include +# include # include # include # include @@ -488,32 +489,32 @@ namespace PETScWrappers AssertThrow (ierr == 0, ExcPETScError(ierr)); if (additional_data.output_details) - PetscOptionsSetValue("-pc_hypre_boomeramg_print_statistics","1"); + { + set_option_value("-pc_hypre_boomeramg_print_statistics", "1"); + } - PetscOptionsSetValue("-pc_hypre_boomeramg_agg_nl", - Utilities::int_to_string( - additional_data.aggressive_coarsening_num_levels - ).c_str()); + set_option_value("-pc_hypre_boomeramg_agg_nl", + Utilities::to_string(additional_data.aggressive_coarsening_num_levels)); std::stringstream ssStream; ssStream << additional_data.max_row_sum; - PetscOptionsSetValue("-pc_hypre_boomeramg_max_row_sum", ssStream.str().c_str()); + set_option_value("-pc_hypre_boomeramg_max_row_sum", ssStream.str()); ssStream.str(""); // empty the stringstream ssStream << additional_data.strong_threshold; - PetscOptionsSetValue("-pc_hypre_boomeramg_strong_threshold", ssStream.str().c_str()); + set_option_value("-pc_hypre_boomeramg_strong_threshold", ssStream.str()); if (additional_data.symmetric_operator) { - PetscOptionsSetValue("-pc_hypre_boomeramg_relax_type_up", "symmetric-SOR/Jacobi"); - PetscOptionsSetValue("-pc_hypre_boomeramg_relax_type_down", "symmetric-SOR/Jacobi"); - PetscOptionsSetValue("-pc_hypre_boomeramg_relax_type_coarse", "Gaussian-elimination"); + set_option_value("-pc_hypre_boomeramg_relax_type_up", "symmetric-SOR/Jacobi"); + set_option_value("-pc_hypre_boomeramg_relax_type_down", "symmetric-SOR/Jacobi"); + set_option_value("-pc_hypre_boomeramg_relax_type_coarse", "Gaussian-elimination"); } else { - PetscOptionsSetValue("-pc_hypre_boomeramg_relax_type_up", "SOR/Jacobi"); - PetscOptionsSetValue("-pc_hypre_boomeramg_relax_type_down", "SOR/Jacobi"); - PetscOptionsSetValue("-pc_hypre_boomeramg_relax_type_coarse", "Gaussian-elimination"); + set_option_value("-pc_hypre_boomeramg_relax_type_up", "SOR/Jacobi"); + set_option_value("-pc_hypre_boomeramg_relax_type_down", "SOR/Jacobi"); + set_option_value("-pc_hypre_boomeramg_relax_type_coarse", "Gaussian-elimination"); } ierr = PCSetFromOptions (pc); @@ -592,7 +593,9 @@ namespace PETScWrappers AssertThrow (ierr == 0, ExcPETScError(ierr)); if (additional_data.output_details) - PetscOptionsSetValue("-pc_hypre_parasails_logging","1"); + { + set_option_value("-pc_hypre_parasails_logging","1"); + } Assert ((additional_data.symmetric == 0 || additional_data.symmetric == 1 || @@ -626,20 +629,18 @@ namespace PETScWrappers ExcMessage("ParaSails parameter symmetric can only be equal to 0, 1, 2!")); }; - PetscOptionsSetValue("-pc_hypre_parasails_sym",ssStream.str().c_str()); + set_option_value("-pc_hypre_parasails_sym",ssStream.str()); - PetscOptionsSetValue("-pc_hypre_parasails_nlevels", - Utilities::int_to_string( - additional_data.n_levels - ).c_str()); + set_option_value ("-pc_hypre_parasails_nlevels", + Utilities::to_string(additional_data.n_levels)); ssStream.str(""); // empty the stringstream ssStream << additional_data.threshold; - PetscOptionsSetValue("-pc_hypre_parasails_thresh", ssStream.str().c_str()); + set_option_value("-pc_hypre_parasails_thresh", ssStream.str()); ssStream.str(""); // empty the stringstream ssStream << additional_data.filter; - PetscOptionsSetValue("-pc_hypre_parasails_filter", ssStream.str().c_str()); + set_option_value("-pc_hypre_parasails_filter", ssStream.str()); ierr = PCSetFromOptions (pc); AssertThrow (ierr == 0, ExcPETScError(ierr)); diff --git a/tests/slepc/solve_01.cc b/tests/slepc/solve_01.cc index f1da212bb3..74a089c6e8 100644 --- a/tests/slepc/solve_01.cc +++ b/tests/slepc/solve_01.cc @@ -27,6 +27,7 @@ #include #include #include +#include #include #include #include @@ -110,9 +111,9 @@ int main(int argc, char **argv) PETScWrappers::Vector(dim)); std::vector v(n_eigenvalues); - PetscOptionsSetValue("-st_ksp_type","cg"); - PetscOptionsSetValue("-st_pc_type", "jacobi"); - PetscOptionsSetValue("-st_ksp_tol", "1e-15"); + PETScWrappers::set_option_value("-st_ksp_type","cg"); + PETScWrappers::set_option_value("-st_pc_type", "jacobi"); + PETScWrappers::set_option_value("-st_ksp_tol", "1e-15"); { SLEPcWrappers::SolverKrylovSchur solver(control); @@ -129,7 +130,7 @@ int main(int argc, char **argv) check_solve (solver, control, A,B,u,v); } - PetscOptionsSetValue("-st_ksp_type","preonly"); + PETScWrappers::set_option_value("-st_ksp_type","preonly"); { SLEPcWrappers::SolverGeneralizedDavidson solver(control); check_solve (solver, control, A,B,u,v); @@ -141,8 +142,8 @@ int main(int argc, char **argv) check_solve (solver, control, A,B,u,v); } - PetscOptionsSetValue("-st_ksp_type","cg"); - PetscOptionsSetValue("-st_ksp_max_it", "10"); + PETScWrappers::set_option_value("-st_ksp_type","cg"); + PETScWrappers::set_option_value("-st_ksp_max_it", "10"); { SLEPcWrappers::SolverJacobiDavidson solver(control); check_solve (solver, control, A,B,u,v); diff --git a/tests/slepc/solve_04.cc b/tests/slepc/solve_04.cc index 6fe702ba43..0cf353a712 100644 --- a/tests/slepc/solve_04.cc +++ b/tests/slepc/solve_04.cc @@ -27,6 +27,7 @@ #include #include #include +#include #include #include #include @@ -133,9 +134,9 @@ int main(int argc, char **argv) // set extra settings for JD; Otherwise, at least on OSX, // the number of eigensolver iterations is different between debug and release modes! - PetscOptionsSetValue("-st_ksp_type","cg"); - PetscOptionsSetValue("-st_pc_type", "jacobi"); - PetscOptionsSetValue("-st_ksp_max_it", "10"); + PETScWrappers::set_option_value("-st_ksp_type","cg"); + PETScWrappers::set_option_value("-st_pc_type", "jacobi"); + PETScWrappers::set_option_value("-st_ksp_max_it", "10"); { SLEPcWrappers::SolverJacobiDavidson solver(control); check_solve (solver, control, A,u,v);