]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use correct types for sundials 7.
authorMarc Fehling <mafehling.git@gmail.com>
Fri, 12 Apr 2024 09:02:39 +0000 (11:02 +0200)
committerMarc Fehling <mafehling.git@gmail.com>
Fri, 12 Apr 2024 10:06:10 +0000 (12:06 +0200)
include/deal.II/sundials/sunlinsol_wrapper.h
source/sundials/arkode.cc
source/sundials/ida.cc
source/sundials/kinsol.cc
source/sundials/sunlinsol_wrapper.cc

index eeb98898ccb8c0cea4d68a25dc01118f1ad02812..107b768f875ae3a3a5dc967e90e861d01f823676 100644 (file)
@@ -59,7 +59,9 @@ namespace SUNDIALS
      * @param a_times_fn A function pointer to the function that computes A*v
      * @param linsol_ctx The context object used to set up the linear solver and all vectors
      */
-    SundialsOperator(void *A_data, ATimesFn a_times_fn, SUNContext linsol_ctx);
+    SundialsOperator(void       *A_data,
+                     SUNATimesFn a_times_fn,
+                     SUNContext  linsol_ctx);
 #  else
     /**
      * Constructor.
@@ -76,17 +78,23 @@ namespace SUNDIALS
      */
     void *A_data;
 
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
     /**
      * %Function pointer declared by SUNDIALS to evaluate the matrix vector
      * product.
      */
-    ATimesFn a_times_fn;
+    SUNATimesFn a_times_fn;
 
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
     /**
      * Context object used for SUNDIALS logging.
      */
     SUNContext linsol_ctx;
+#  else
+    /**
+     * %Function pointer declared by SUNDIALS to evaluate the matrix vector
+     * product.
+     */
+    ATimesFn a_times_fn;
 #  endif
   };
 
@@ -124,10 +132,10 @@ namespace SUNDIALS
      * 6.0.0. If you are using an earlier version of SUNDIALS then you need to
      * use the other constructor.
      */
-    SundialsPreconditioner(void      *P_data,
-                           PSolveFn   p_solve_fn,
-                           SUNContext linsol_ctx,
-                           double     tol);
+    SundialsPreconditioner(void       *P_data,
+                           SUNPSolveFn p_solve_fn,
+                           SUNContext  linsol_ctx,
+                           double      tol);
 #  else
     /**
      * Constructor.
@@ -150,17 +158,23 @@ namespace SUNDIALS
      */
     void *P_data;
 
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
     /**
      * %Function pointer to a function that computes the preconditioner
      * application.
      */
-    PSolveFn p_solve_fn;
+    SUNPSolveFn p_solve_fn;
 
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
     /**
      * Context object used for SUNDIALS logging.
      */
     SUNContext linsol_ctx;
+#  else
+    /**
+     * %Function pointer to a function that computes the preconditioner
+     * application.
+     */
+    PSolveFn p_solve_fn;
 #  endif
 
     /**
index 38d8ca06f3b4a000622f5d143cff3b456938eb20..a1a909290147ea98182580867ee4e110121f4663 100644 (file)
@@ -69,12 +69,19 @@ namespace SUNDIALS
   {
     set_functions_to_trigger_an_assert();
 
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
     // SUNDIALS will always duplicate communicators if we provide them. This
     // can cause problems if SUNDIALS is configured with MPI and we pass along
     // MPI_COMM_SELF in a serial application as MPI won't be
     // initialized. Hence, work around that by just not providing a
     // communicator in that case.
+#  if DEAL_II_SUNDIALS_VERSION_GTE(7, 0, 0)
+    const int status =
+      SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? SUN_COMM_NULL :
+                                                            mpi_communicator,
+                        &arkode_ctx);
+    (void)status;
+    AssertARKode(status);
+#  elif DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
     const int status =
       SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
                                                             &mpi_communicator,
@@ -239,9 +246,20 @@ namespace SUNDIALS
     int status;
     (void)status;
 
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+#  if DEAL_II_SUNDIALS_VERSION_GTE(7, 0, 0)
     status = SUNContext_Free(&arkode_ctx);
     AssertARKode(status);
+
+    // Same comment applies as in class constructor:
+    status =
+      SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? SUN_COMM_NULL :
+                                                            mpi_communicator,
+                        &arkode_ctx);
+    AssertARKode(status);
+#  elif DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+    status = SUNContext_Free(&arkode_ctx);
+    AssertARKode(status);
+
     // Same comment applies as in class constructor:
     status =
       SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
@@ -409,7 +427,7 @@ namespace SUNDIALS
 #  else
             sun_linear_solver =
               SUNLinSol_SPGMR(y_template,
-                              PREC_NONE,
+                              SUN_PREC_NONE,
                               0 /*krylov subvectors, 0 uses default*/,
                               arkode_ctx);
 #  endif
@@ -505,8 +523,13 @@ namespace SUNDIALS
             auto jacobian_solver_setup_callback = [](SUNDIALS::realtype t,
                                                      N_Vector           y,
                                                      N_Vector           fy,
-                                                     booleantype        jok,
-                                                     booleantype       *jcurPtr,
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                                     sunbooleantype  jok,
+                                                     sunbooleantype *jcurPtr,
+#  else
+                                                     booleantype  jok,
+                                                     booleantype *jcurPtr,
+#  endif
                                                      SUNDIALS::realtype gamma,
                                                      void *user_data) -> int {
               Assert(user_data != nullptr, ExcInternalError());
@@ -612,13 +635,20 @@ namespace SUNDIALS
 #  else
             sun_mass_linear_solver =
               SUNLinSol_SPGMR(y_template,
-                              PREC_NONE,
+                              SUN_PREC_NONE,
                               0 /*krylov subvectors, 0 uses default*/,
                               arkode_ctx);
 #  endif
           }
+
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+        sunbooleantype mass_time_dependent =
+          data.mass_is_time_independent ? SUNFALSE : SUNTRUE;
+#  else
         booleantype mass_time_dependent =
           data.mass_is_time_independent ? SUNFALSE : SUNTRUE;
+#  endif
+
         status = ARKStepSetMassLinearSolver(arkode_mem,
                                             sun_mass_linear_solver,
                                             nullptr,
index 65a7093ab7845402cef6a24985af13025f94523c..e9d137a6294f34f16685db32345f642371c22e2c 100644 (file)
@@ -61,12 +61,21 @@ namespace SUNDIALS
     , mpi_communicator(mpi_comm)
     , pending_exception(nullptr)
   {
+    set_functions_to_trigger_an_assert();
+
     // SUNDIALS will always duplicate communicators if we provide them. This
     // can cause problems if SUNDIALS is configured with MPI and we pass along
     // MPI_COMM_SELF in a serial application as MPI won't be
     // initialized. Hence, work around that by just not providing a
     // communicator in that case.
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+#  if DEAL_II_SUNDIALS_VERSION_GTE(7, 0, 0)
+    const int status =
+      SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? SUN_COMM_NULL :
+                                                            mpi_communicator,
+                        &ida_ctx);
+    (void)status;
+    AssertIDA(status);
+#  elif DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
     const int status =
       SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
                                                             &mpi_communicator,
@@ -74,7 +83,6 @@ namespace SUNDIALS
     (void)status;
     AssertIDA(status);
 #  endif
-    set_functions_to_trigger_an_assert();
   }
 
 
@@ -187,7 +195,17 @@ namespace SUNDIALS
     int status;
     (void)status;
 
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+#  if DEAL_II_SUNDIALS_VERSION_GTE(7, 0, 0)
+    status = SUNContext_Free(&ida_ctx);
+    AssertIDA(status);
+
+    // Same comment applies as in class constructor:
+    status =
+      SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? SUN_COMM_NULL :
+                                                            mpi_communicator,
+                        &ida_ctx);
+    AssertIDA(status);
+#  elif DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
     status = SUNContext_Free(&ida_ctx);
     AssertIDA(status);
 
index 34fbc8cc80e59a99d79c9ec1b72ba78454cb32bb..7ad8d70afda1b003570aee77e26201fcac395931 100644 (file)
 // Make sure we #include the SUNDIALS config file...
 #  include <sundials/sundials_config.h>
 // ...before the rest of the SUNDIALS files:
-#  include <kinsol/kinsol_direct.h>
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+#    include <kinsol/kinsol_ls.h>
+#  else
+#    include <kinsol/kinsol_direct.h>
+#  endif
 #  include <sunlinsol/sunlinsol_dense.h>
 #  include <sunmatrix/sunmatrix_dense.h>
 
@@ -175,12 +179,19 @@ namespace SUNDIALS
   {
     set_functions_to_trigger_an_assert();
 
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
     // SUNDIALS will always duplicate communicators if we provide them. This
     // can cause problems if SUNDIALS is configured with MPI and we pass along
     // MPI_COMM_SELF in a serial application as MPI won't be
     // initialized. Hence, work around that by just not providing a
     // communicator in that case.
+#  if DEAL_II_SUNDIALS_VERSION_GTE(7, 0, 0)
+    const int status =
+      SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? SUN_COMM_NULL :
+                                                            mpi_communicator,
+                        &kinsol_ctx);
+    (void)status;
+    AssertKINSOL(status);
+#  elif DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
     const int status =
       SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
                                                             &mpi_communicator,
@@ -234,9 +245,17 @@ namespace SUNDIALS
     AssertKINSOL(status);
 #  endif
 
-#  if DEAL_II_SUNDIALS_VERSION_LT(6, 0, 0)
-    kinsol_mem = KINCreate();
-#  else
+
+#  if DEAL_II_SUNDIALS_VERSION_GTE(7, 0, 0)
+    // Same comment applies as in class constructor:
+    status =
+      SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? SUN_COMM_NULL :
+                                                            mpi_communicator,
+                        &kinsol_ctx);
+    AssertKINSOL(status);
+
+    kinsol_mem = KINCreate(kinsol_ctx);
+#  elif DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
     // Same comment applies as in class constructor:
     status =
       SUNContext_Create(mpi_communicator == MPI_COMM_SELF ? nullptr :
@@ -245,6 +264,8 @@ namespace SUNDIALS
     AssertKINSOL(status);
 
     kinsol_mem = KINCreate(kinsol_ctx);
+#  else
+    kinsol_mem = KINCreate();
 #  endif
 
     status = KINSetUserData(kinsol_mem, static_cast<void *>(this));
index e20124eedfb67aa238f3697fb0e4c7864127e554..4d728b7e78a4806dc901ce6c9d98e30e72cff16b 100644 (file)
@@ -144,12 +144,16 @@ namespace SUNDIALS
         , pending_exception(pending_exception)
       {}
 
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+      SUNATimesFn a_times_fn;
+      SUNPSetupFn preconditioner_setup;
+      SUNPSolveFn preconditioner_solve;
+
+      SUNContext linsol_ctx;
+#  else
       ATimesFn a_times_fn;
       PSetupFn preconditioner_setup;
       PSolveFn preconditioner_solve;
-
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-      SUNContext linsol_ctx;
 #  endif
 
       LinearSolveFunction<VectorType> lsolve;
@@ -254,10 +258,15 @@ namespace SUNDIALS
     }
 
 
-
     template <typename VectorType>
     int
-    arkode_linsol_set_a_times(SUNLinearSolver LS, void *A_data, ATimesFn ATimes)
+    arkode_linsol_set_a_times(SUNLinearSolver LS,
+                              void           *A_data,
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                              SUNATimesFn ATimes)
+#  else
+                              ATimesFn        ATimes)
+#  endif
     {
       LinearSolverContent<VectorType> *content = access_content<VectorType>(LS);
 
@@ -272,8 +281,13 @@ namespace SUNDIALS
     int
     arkode_linsol_set_preconditioner(SUNLinearSolver LS,
                                      void           *P_data,
-                                     PSetupFn        p_setup,
-                                     PSolveFn        p_solve)
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+                                     SUNPSetupFn p_setup,
+                                     SUNPSolveFn p_solve)
+#  else
+                                     PSetupFn p_setup,
+                                     PSolveFn p_solve)
+#  endif
     {
       LinearSolverContent<VectorType> *content = access_content<VectorType>(LS);
 
@@ -342,9 +356,9 @@ namespace SUNDIALS
 
 #  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
   template <typename VectorType>
-  SundialsOperator<VectorType>::SundialsOperator(void      *A_data,
-                                                 ATimesFn   a_times_fn,
-                                                 SUNContext linsol_ctx)
+  SundialsOperator<VectorType>::SundialsOperator(void       *A_data,
+                                                 SUNATimesFn a_times_fn,
+                                                 SUNContext  linsol_ctx)
     : A_data(A_data)
     , a_times_fn(a_times_fn)
     , linsol_ctx(linsol_ctx)
@@ -394,10 +408,10 @@ namespace SUNDIALS
 #  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
   template <typename VectorType>
   SundialsPreconditioner<VectorType>::SundialsPreconditioner(
-    void      *P_data,
-    PSolveFn   p_solve_fn,
-    SUNContext linsol_ctx,
-    double     tol)
+    void       *P_data,
+    SUNPSolveFn p_solve_fn,
+    SUNContext  linsol_ctx,
+    double      tol)
     : P_data(P_data)
     , p_solve_fn(p_solve_fn)
     , linsol_ctx(linsol_ctx)

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.