]> https://gitweb.dealii.org/ - dealii.git/commitdiff
PETSc: Fix backward compatibility 15368/head
authorStefano Zampini <stefano.zampini@gmail.com>
Thu, 15 Jun 2023 22:26:28 +0000 (01:26 +0300)
committerStefano Zampini <stefano.zampini@gmail.com>
Thu, 15 Jun 2023 22:51:28 +0000 (01:51 +0300)
include/deal.II/lac/petsc_compatibility.h
include/deal.II/lac/petsc_snes.templates.h
include/deal.II/lac/petsc_ts.templates.h
source/lac/petsc_compatibility.cc
source/lac/petsc_precondition.cc
tests/petsc/vector_wrap_01.cc
tests/petsc_complex/17.cc
tests/petsc_complex/vector_wrap_01.cc
tests/tests.h

index 1b94af490336944d7a24dfdf7b9ad7d010462f31..bd061ea8fb9078b53df0536b9ebd0e05ed88f9d5 100644 (file)
 /*
  * Rather than using ifdefs everywhere, try to wrap older versions of PETSc
  * functions in one place.
+ *
+ * Functions that are not inlined are:
+ * - Functions returning PetscErrorCode that are supposed to be called within
+ *   PETSc callbacks.
+ * - Functions that need access to internal PETSc headers that we don't want
+ *   to expose to deal.II users
  */
 #ifndef dealii_petsc_compatibility_h
 #define dealii_petsc_compatibility_h
@@ -98,24 +104,47 @@ namespace PETScWrappers
     set_matrix_option(matrix, MAT_KEEP_NONZERO_PATTERN, PETSC_TRUE);
   }
 
+
+
   /**
    * Tell PETSc that the status of the vector has changed.
    */
   void
   petsc_increment_state_counter(Vec v);
 
+
+
   /**
    * Tell PETSc that the status of the matrix has changed.
    */
   void
   petsc_increment_state_counter(Mat A);
 
+
+  /**
+   * Set the failed reason for the preconditioner.
+   */
+  PetscErrorCode
+  pc_set_failed_reason(PC pc, PCFailedReason reason);
+
+
+
   /**
    * Resets internal domain error flags in the SNES object.
    */
   void
   snes_reset_domain_flags(SNES snes);
 
+
+
+  /**
+   * Resets internal domain error flags in the SNES object.
+   */
+  void
+  snes_set_jacobian_domain_error(SNES snes);
+
+
+
   /**
    * Tell PETSc nonlinear solver to use matrix free finite differencing (MFFD).
    *
@@ -128,6 +157,8 @@ namespace PETScWrappers
   void
   set_use_matrix_free(SNES snes, const bool mf_operator, const bool mf);
 
+
+
   /**
    * Tell PETSc ODE solver to use matrix free finite differencing (MFFD).
    *
@@ -140,24 +171,32 @@ namespace PETScWrappers
   void
   set_use_matrix_free(TS ts, const bool mf_operator, const bool mf);
 
+
+
   /**
    * Set final time for ODE integration.
    */
   void
   ts_set_max_time(TS ts, const PetscReal maxtime);
 
+
+
   /**
    * Set maximum number of steps for ODE integration.
    */
   void
   ts_set_max_steps(TS ts, const PetscInt maxsteps);
 
+
+
   /**
    * Return current step number.
    */
   unsigned int
   ts_get_step_number(TS ts);
 
+
+
   /**
    * Return true if the TS has a SNES object.
    */
index cf5af01c4d551f0a89e66b5f4a50462eba2a4285..f9dbfb674e641b7eaab1e17c4607a68830c0c945 100644 (file)
@@ -226,7 +226,9 @@ namespace PETScWrappers
     AssertPETSc(SNESReset(snes));
     // By default PETSc does not check for Jacobian errors in optimized
     // mode. Here we do it unconditionally.
+#  if DEAL_II_PETSC_VERSION_GTE(3, 11, 0)
     AssertPETSc(SNESSetCheckJacobianDomainError(snes, PETSC_TRUE));
+#  endif
   }
 
 
@@ -378,7 +380,7 @@ namespace PETScWrappers
       const int err    = call_and_possibly_capture_snes_exception(
         user->jacobian,
         user->pending_exception,
-        [snes]() -> void { AssertPETSc(SNESSetJacobianDomainError(snes)); },
+        [snes]() -> void { snes_set_jacobian_domain_error(snes); },
         xdealii,
         Adealii,
         Pdealii);
@@ -424,7 +426,7 @@ namespace PETScWrappers
       const int err    = call_and_possibly_capture_snes_exception(
         user->setup_jacobian,
         user->pending_exception,
-        [snes]() -> void { AssertPETSc(SNESSetJacobianDomainError(snes)); },
+        [snes]() -> void { snes_set_jacobian_domain_error(snes); },
         xdealii);
       if (err)
         return PetscError(
index 3061a1776453b5023035ad664f61b60b2fa7172a..f071222ad9fe433ac962a90fc21bb8272b0ea666 100644 (file)
@@ -492,7 +492,7 @@ namespace PETScWrappers
         [ts]() -> void {
           SNES snes;
           AssertPETSc(TSGetSNES(ts, &snes));
-          AssertPETSc(SNESSetJacobianDomainError(snes));
+          snes_set_jacobian_domain_error(snes);
         },
         t,
         xdealii,
@@ -553,7 +553,7 @@ namespace PETScWrappers
         [ts]() -> void {
           SNES snes;
           AssertPETSc(TSGetSNES(ts, &snes));
-          AssertPETSc(SNESSetJacobianDomainError(snes));
+          snes_set_jacobian_domain_error(snes);
         },
         t,
         xdealii,
@@ -653,7 +653,7 @@ namespace PETScWrappers
         [ts]() -> void {
           SNES snes;
           AssertPETSc(TSGetSNES(ts, &snes));
-          AssertPETSc(SNESSetJacobianDomainError(snes));
+          snes_set_jacobian_domain_error(snes);
         },
         t,
         xdealii,
@@ -869,13 +869,14 @@ namespace PETScWrappers
 
     // By default PETSc does not check for Jacobian errors in optimized
     // mode. Here we do it unconditionally.
+#  if DEAL_II_PETSC_VERSION_GTE(3, 11, 0)
     if (ts_has_snes(ts))
       {
         SNES snes;
         AssertPETSc(TSGetSNES(ts, &snes));
         AssertPETSc(SNESSetCheckJacobianDomainError(snes, PETSC_TRUE));
       }
-
+#  endif
     // Having set everything up, now do the actual work
     // and let PETSc do the time stepping. If there is
     // a pending exception, then one of the user callbacks
index ad702fce50cb8c847f872972beff5828bf987b0c..b787052001c1fe6fc5571de05bcc54cd556a4f87 100644 (file)
@@ -22,6 +22,7 @@
 
 #ifdef DEAL_II_WITH_PETSC
 
+#  include <petsc/private/pcimpl.h>
 #  include <petsc/private/petscimpl.h>
 #  include <petsc/private/snesimpl.h>
 #  include <petsc/private/tsimpl.h>
@@ -52,11 +53,36 @@ namespace PETScWrappers
     AssertPETSc(PetscObjectStateIncrease(reinterpret_cast<PetscObject>(A)));
   }
 
+  PetscErrorCode
+  pc_set_failed_reason(PC pc, PCFailedReason reason)
+  {
+#  if DEAL_II_PETSC_VERSION_GTE(3, 14, 0)
+    return PCSetFailedReason(pc, reason);
+#  else
+    pc->failedreason = reason;
+    return 0;
+#  endif
+  }
+
   void
   snes_reset_domain_flags(SNES snes)
   {
+#  if DEAL_II_PETSC_VERSION_GTE(3, 11, 0)
     snes->jacobiandomainerror = PETSC_FALSE;
-    snes->domainerror         = PETSC_FALSE;
+#  endif
+    snes->domainerror = PETSC_FALSE;
+  }
+
+  void
+  snes_set_jacobian_domain_error(SNES snes)
+  {
+#  if DEAL_II_PETSC_VERSION_GTE(3, 11, 0)
+    snes->jacobiandomainerror = PETSC_TRUE;
+#  else
+    // There is no equivalent, and since this used to stop
+    // computations, we opt to set the converged reason
+    snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN;
+#  endif
   }
 
   void
index 94f1122572f83be715d6a6c33afdc4ea1aad6085..f018a7bdc0481a402ca73f1aed4126e2df6bd169 100644 (file)
@@ -1113,7 +1113,7 @@ namespace PETScWrappers
     // Failed reason is not reset uniformly within the
     // interface code of PCSetUp in PETSc.
     // We handle it here.
-    PetscCall(PCSetFailedReason(ppc, PC_NOERROR));
+    PetscCall(pc_set_failed_reason(ppc, PC_NOERROR));
     PetscFunctionReturn(0);
   }
 
@@ -1141,7 +1141,7 @@ namespace PETScWrappers
       }
     catch (const RecoverableUserCallbackError &)
       {
-        PetscCall(PCSetFailedReason(ppc, PC_SUBPC_ERROR));
+        PetscCall(pc_set_failed_reason(ppc, PC_SUBPC_ERROR));
       }
     catch (...)
       {
@@ -1182,7 +1182,7 @@ namespace PETScWrappers
       }
     catch (const RecoverableUserCallbackError &)
       {
-        PetscCall(PCSetFailedReason(ppc, PC_SUBPC_ERROR));
+        PetscCall(pc_set_failed_reason(ppc, PC_SUBPC_ERROR));
       }
     catch (...)
       {
index dff5c14c73ae690920856c72d9b16328c906314b..2dc83293b4c2d7f72f0b92da1a434a31c2b6b955 100644 (file)
@@ -61,11 +61,7 @@ main(int argc, char **argv)
         test(v, w);
       }
 
-#if DEAL_II_PETSC_VERSION_LT(3, 2, 0)
-      ierr = VecDestroy(vpetsc);
-#else
       ierr = VecDestroy(&vpetsc);
-#endif
 
       AssertThrow(ierr == 0, ExcPETScError(ierr));
     }
index fb1544922969884d4f8a7d7c5f833709abbd9a0e..17653fbaa013d47442a9db3493b08c4dad6f447e 100644 (file)
@@ -34,11 +34,7 @@ test(PETScWrappers::MPI::Vector &v)
     {
       const PetscScalar s(1. * k, 2. * k);
       v(k) = s;
-#if DEAL_II_PETSC_VERSION_LT(3, 7, 0)
-      norm += std::fabs(1. * k) + std::fabs(2. * k);
-#else
       norm += std::abs(s);
-#endif
     }
   v.compress(VectorOperation::insert);
 
index 7acb96b363da18fa352074f6d3e233bc48f804a3..044cde359b06a22c74545fa6f1185efc8db1c182 100644 (file)
@@ -61,11 +61,7 @@ main(int argc, char **argv)
         test(v, w);
       }
 
-#if DEAL_II_PETSC_VERSION_LT(3, 2, 0)
-      ierr = VecDestroy(vpetsc);
-#else
       ierr = VecDestroy(&vpetsc);
-#endif
 
       AssertThrow(ierr == 0, ExcPETScError(ierr));
     }
index 16185387d93096a08af533e4f20117ec858a2982..d5c6250d5f54cf77a47aef57bde86b86b54ea6db 100644 (file)
@@ -461,7 +461,6 @@ namespace
   void
   check_petsc_allocations()
   {
-#  if DEAL_II_PETSC_VERSION_GTE(3, 2, 0)
     PetscStageLog stageLog;
     PetscLogGetStageLog(&stageLog);
 
@@ -492,7 +491,6 @@ namespace
 
     if (errors)
       throw dealii::ExcMessage("PETSc memory leak");
-#  endif
   }
 } // namespace
 #endif

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.