]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix compilation with PETSc 3.7.0. 2754/head
authorDavid Wells <wellsd2@rpi.edu>
Thu, 7 Jul 2016 09:50:59 +0000 (05:50 -0400)
committerDavid Wells <wellsd2@rpi.edu>
Thu, 7 Jul 2016 12:59:56 +0000 (08:59 -0400)
doc/news/changes.h
include/deal.II/lac/petsc_compatibility.h [new file with mode: 0644]
source/lac/petsc_precondition.cc
tests/slepc/solve_01.cc
tests/slepc/solve_04.cc

index 6e393a1905122e50d3680b40adbc9515b29026b5..4a537c072058682c59e93447160840ccb1823b13 100644 (file)
@@ -150,6 +150,12 @@ inconvenience this causes.
 <h3>General</h3>
 
 <ol>
+ <li> New: The library is now compatible with PETSc 3.7.0. Part of this change
+ included adding a new header, <tt>petsc_compatibility.h</tt>, which provides
+ some version-independent functions for using common PETSc functions.
+ <br>
+ (David Wells, 2016/07/07)
+ </li>
 
  <li> 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 (file)
index 0000000..705ffa0
--- /dev/null
@@ -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 <deal.II/base/config.h>
+
+#ifdef DEAL_II_WITH_PETSC
+
+#include <petscconf.h>
+#include <petscpc.h>
+
+#include <string>
+
+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
index 7dd465f16ead90ab3d65093ab43377f0c674e385..96bd524b527e92fb1552b0452466a75d781cba07 100644 (file)
@@ -18,6 +18,7 @@
 #ifdef DEAL_II_WITH_PETSC
 
 #  include <deal.II/base/utilities.h>
+#  include <deal.II/lac/petsc_compatibility.h>
 #  include <deal.II/lac/petsc_matrix_base.h>
 #  include <deal.II/lac/petsc_vector_base.h>
 #  include <deal.II/lac/petsc_solver.h>
@@ -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));
index f1da212bb3ee1a959ee7c60d6b46e3436822fb45..74a089c6e87b4d7e145a60164114a128e7a407aa 100644 (file)
@@ -27,6 +27,7 @@
 #include <iostream>
 #include <iomanip>
 #include <deal.II/base/logstream.h>
+#include <deal.II/lac/petsc_compatibility.h>
 #include <deal.II/lac/petsc_sparse_matrix.h>
 #include <deal.II/lac/petsc_vector.h>
 #include <deal.II/lac/petsc_solver.h>
@@ -110,9 +111,9 @@ int main(int argc, char **argv)
                                           PETScWrappers::Vector(dim));
     std::vector<PetscScalar> 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);
index 6fe702ba4350e80a136257a6dd7aa5718e951f67..0cf353a7127028aff5bb5e5d093aec0c6e01129d 100644 (file)
@@ -27,6 +27,7 @@
 #include <iostream>
 #include <iomanip>
 #include <deal.II/base/logstream.h>
+#include <deal.II/lac/petsc_compatibility.h>
 #include <deal.II/lac/petsc_sparse_matrix.h>
 #include <deal.II/lac/petsc_vector.h>
 #include <deal.II/lac/petsc_solver.h>
@@ -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);

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.