]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Introduce sundials_types.h 16271/head
authorDaniel Arndt <arndtd@ornl.gov>
Sun, 19 Nov 2023 15:34:38 +0000 (10:34 -0500)
committerDaniel Arndt <arndtd@ornl.gov>
Sun, 19 Nov 2023 19:04:29 +0000 (14:04 -0500)
include/deal.II/sundials/arkode.h
include/deal.II/sundials/ida.h
include/deal.II/sundials/kinsol.h
include/deal.II/sundials/n_vector.templates.h
include/deal.II/sundials/sundials_types.h [new file with mode: 0644]
source/sundials/arkode.cc
source/sundials/ida.cc
source/sundials/kinsol.cc
source/sundials/sunlinsol_wrapper.cc
tests/sundials/arkode_06.cc
tests/sundials/arkode_07.cc

index 7e254952cbbbab7f434614f346a59eb0e153f0f1..f969eb80e700c39828b7299ea321bc1bff79fe57 100644 (file)
 #  include <deal.II/base/discrete_time.h>
 
 #  include <deal.II/sundials/n_vector.h>
+#  include <deal.II/sundials/sundials_types.h>
 #  include <deal.II/sundials/sunlinsol_wrapper.h>
 
 #  include <boost/signals2.hpp>
 
 #  include <sundials/sundials_linearsolver.h>
 #  include <sundials/sundials_math.h>
-#  include <sundials/sundials_types.h>
 
 #  include <exception>
 #  include <memory>
index 45dcc0d141ef45b197b9bb9d023d94b3399be2ab..9e6acf22b2d13791ce71d226ca3e82cba46a3a3d 100644 (file)
@@ -39,6 +39,7 @@
 #    include <ida/ida.h>
 #  endif
 
+#  include <deal.II/sundials/sundials_types.h>
 #  include <deal.II/sundials/sunlinsol_wrapper.h>
 
 #  include <boost/signals2.hpp>
@@ -46,7 +47,6 @@
 #  include <nvector/nvector_serial.h>
 #  include <sundials/sundials_config.h>
 #  include <sundials/sundials_math.h>
-#  include <sundials/sundials_types.h>
 
 #  include <memory>
 
index 1c980951dad14c711718f00988c42485ce3f28bc..3f098c7d775c2228a4504e411def1d52eae25baa 100644 (file)
 #  include <deal.II/lac/vector.h>
 #  include <deal.II/lac/vector_memory.h>
 
+#  include <deal.II/sundials/sundials_types.h>
+
 #  include <boost/signals2.hpp>
 
 #  include <kinsol/kinsol.h>
 #  include <nvector/nvector_serial.h>
 #  include <sundials/sundials_math.h>
-#  include <sundials/sundials_types.h>
 
 #  include <exception>
 #  include <memory>
index 1fe1c0cf9a702bc7dd8790f36b70adc77e4bdaa3..12455a9c2222f102f24c36226dd177344dc79509 100644 (file)
@@ -20,6 +20,7 @@
 #include <deal.II/base/config.h>
 
 #include <deal.II/sundials/n_vector.h>
+#include <deal.II/sundials/sundials_types.h>
 
 #ifdef DEAL_II_WITH_SUNDIALS
 
@@ -176,12 +177,6 @@ namespace SUNDIALS
      */
     namespace NVectorOperations
     {
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-      using realtype = ::sunrealtype;
-#  else
-      using realtype = ::realtype;
-#  endif
-
       N_Vector_ID
       get_vector_id(N_Vector v);
 
@@ -202,18 +197,22 @@ namespace SUNDIALS
 
       template <typename VectorType>
       void
-      linear_sum(realtype a, N_Vector x, realtype b, N_Vector y, N_Vector z);
+      linear_sum(SUNDIALS::realtype a,
+                 N_Vector           x,
+                 SUNDIALS::realtype b,
+                 N_Vector           y,
+                 N_Vector           z);
 
       template <typename VectorType>
-      realtype
+      SUNDIALS::realtype
       dot_product(N_Vector x, N_Vector y);
 
       template <typename VectorType>
-      realtype
+      SUNDIALS::realtype
       weighted_l2_norm(N_Vector x, N_Vector y);
 
       template <typename VectorType>
-      realtype
+      SUNDIALS::realtype
       l1_norm(N_Vector x);
 
       template <typename VectorType>
@@ -251,47 +250,47 @@ namespace SUNDIALS
       elementwise_abs(N_Vector x, N_Vector z);
 
       template <typename VectorType>
-      realtype
+      SUNDIALS::realtype
       weighted_rms_norm(N_Vector x, N_Vector w);
 
       template <typename VectorType>
-      realtype
+      SUNDIALS::realtype
       weighted_rms_norm_mask(N_Vector x, N_Vector w, N_Vector mask);
 
       template <typename VectorType>
-      realtype
+      SUNDIALS::realtype
       max_norm(N_Vector x);
 
       template <typename VectorType,
                 std::enable_if_t<is_serial_vector<VectorType>::value, int> = 0>
-      realtype
+      SUNDIALS::realtype
       min_element(N_Vector x);
 
       template <typename VectorType,
                 std::enable_if_t<!is_serial_vector<VectorType>::value &&
                                    !IsBlockVector<VectorType>::value,
                                  int> = 0>
-      realtype
+      SUNDIALS::realtype
       min_element(N_Vector x);
 
       template <typename VectorType,
                 std::enable_if_t<!is_serial_vector<VectorType>::value &&
                                    IsBlockVector<VectorType>::value,
                                  int> = 0>
-      realtype
+      SUNDIALS::realtype
       min_element(N_Vector x);
 
       template <typename VectorType>
       void
-      scale(realtype c, N_Vector x, N_Vector z);
+      scale(SUNDIALS::realtype c, N_Vector x, N_Vector z);
 
       template <typename VectorType>
       void
-      set_constant(realtype c, N_Vector v);
+      set_constant(SUNDIALS::realtype c, N_Vector v);
 
       template <typename VectorType>
       void
-      add_constant(N_Vector x, realtype b, N_Vector z);
+      add_constant(N_Vector x, SUNDIALS::realtype b, N_Vector z);
 
       template <typename VectorType>
       const MPI_Comm &
@@ -652,7 +651,11 @@ namespace SUNDIALS
 
       template <typename VectorType>
       void
-      linear_sum(realtype a, N_Vector x, realtype b, N_Vector y, N_Vector z)
+      linear_sum(SUNDIALS::realtype a,
+                 N_Vector           x,
+                 SUNDIALS::realtype b,
+                 N_Vector           y,
+                 N_Vector           z)
       {
         auto *x_dealii = unwrap_nvector_const<VectorType>(x);
         auto *y_dealii = unwrap_nvector_const<VectorType>(y);
@@ -672,7 +675,7 @@ namespace SUNDIALS
 
 
       template <typename VectorType>
-      realtype
+      SUNDIALS::realtype
       dot_product(N_Vector x, N_Vector y)
       {
         return *unwrap_nvector_const<VectorType>(x) *
@@ -682,7 +685,7 @@ namespace SUNDIALS
 
 
       template <typename VectorType>
-      realtype
+      SUNDIALS::realtype
       weighted_l2_norm(N_Vector x, N_Vector w)
       {
         // TODO copy can be avoided by a custom kernel
@@ -695,7 +698,7 @@ namespace SUNDIALS
 
 
       template <typename VectorType>
-      realtype
+      SUNDIALS::realtype
       l1_norm(N_Vector x)
       {
         return unwrap_nvector_const<VectorType>(x)->l1_norm();
@@ -705,7 +708,7 @@ namespace SUNDIALS
 
       template <typename VectorType>
       void
-      set_constant(realtype c, N_Vector v)
+      set_constant(SUNDIALS::realtype c, N_Vector v)
       {
         auto *v_dealii = unwrap_nvector<VectorType>(v);
         *v_dealii      = c;
@@ -715,7 +718,7 @@ namespace SUNDIALS
 
       template <typename VectorType>
       void
-      add_constant(N_Vector x, realtype b, N_Vector z)
+      add_constant(N_Vector x, SUNDIALS::realtype b, N_Vector z)
       {
         auto *x_dealii = unwrap_nvector_const<VectorType>(x);
         auto *z_dealii = unwrap_nvector<VectorType>(z);
@@ -752,7 +755,7 @@ namespace SUNDIALS
 
 
       template <typename VectorType>
-      realtype
+      SUNDIALS::realtype
       weighted_rms_norm(N_Vector x, N_Vector w)
       {
         // TODO copy can be avoided by a custom kernel
@@ -766,7 +769,7 @@ namespace SUNDIALS
 
 
       template <typename VectorType>
-      realtype
+      SUNDIALS::realtype
       weighted_rms_norm_mask(N_Vector x, N_Vector w, N_Vector mask)
       {
         // TODO copy can be avoided by a custom kernel
@@ -782,7 +785,7 @@ namespace SUNDIALS
 
 
       template <typename VectorType>
-      realtype
+      SUNDIALS::realtype
       max_norm(N_Vector x)
       {
         return unwrap_nvector_const<VectorType>(x)->linfty_norm();
@@ -792,7 +795,7 @@ namespace SUNDIALS
 
       template <typename VectorType,
                 std::enable_if_t<is_serial_vector<VectorType>::value, int>>
-      realtype
+      SUNDIALS::realtype
       min_element(N_Vector x)
       {
         auto *vector = unwrap_nvector_const<VectorType>(x);
@@ -805,7 +808,7 @@ namespace SUNDIALS
                 std::enable_if_t<!is_serial_vector<VectorType>::value &&
                                    !IsBlockVector<VectorType>::value,
                                  int>>
-      realtype
+      SUNDIALS::realtype
       min_element(N_Vector x)
       {
         auto *vector = unwrap_nvector_const<VectorType>(x);
@@ -831,7 +834,7 @@ namespace SUNDIALS
                 std::enable_if_t<!is_serial_vector<VectorType>::value &&
                                    IsBlockVector<VectorType>::value,
                                  int>>
-      realtype
+      SUNDIALS::realtype
       min_element(N_Vector x)
       {
         auto *vector = unwrap_nvector_const<VectorType>(x);
@@ -869,7 +872,7 @@ namespace SUNDIALS
 
       template <typename VectorType>
       void
-      scale(realtype c, N_Vector x, N_Vector z)
+      scale(SUNDIALS::realtype c, N_Vector x, N_Vector z)
       {
         auto *x_dealii = unwrap_nvector_const<VectorType>(x);
         auto *z_dealii = unwrap_nvector<VectorType>(z);
diff --git a/include/deal.II/sundials/sundials_types.h b/include/deal.II/sundials/sundials_types.h
new file mode 100644 (file)
index 0000000..66db338
--- /dev/null
@@ -0,0 +1,41 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2023 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_sundials_types_h
+#define dealii_sundials_types_h
+
+#include <deal.II/base/config.h>
+
+#ifdef DEAL_II_WITH_SUNDIALS
+#  include <sundials/sundials_types.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace SUNDIALS
+{
+/**
+ * Alias for the real type used by SUNDIALS.
+ */
+#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
+  using realtype = ::sunrealtype;
+#  else
+  using realtype = ::realtype;
+#  endif
+} // namespace SUNDIALS
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif // DEAL_II_WITH_SUNDIALS
+#endif // dealii_sundials_types_h
index a27a52a06e6a72b508a4136c54756cb60b508c05..abb320af767086a043a85b9cc7881ab64bd8250b 100644 (file)
@@ -269,15 +269,10 @@ namespace SUNDIALS
     Assert(explicit_function || implicit_function,
            ExcFunctionNotProvided("explicit_function || implicit_function"));
 
-    auto explicit_function_callback = [](
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                                        sunrealtype tt,
-#  else
-                                        realtype tt,
-#  endif
-                                        N_Vector yy,
-                                        N_Vector yp,
-                                        void    *user_data) -> int {
+    auto explicit_function_callback = [](SUNDIALS::realtype tt,
+                                         N_Vector           yy,
+                                         N_Vector           yp,
+                                         void              *user_data) -> int {
       Assert(user_data != nullptr, ExcInternalError());
       ARKode<VectorType> &solver =
         *static_cast<ARKode<VectorType> *>(user_data);
@@ -294,15 +289,10 @@ namespace SUNDIALS
     };
 
 
-    auto implicit_function_callback = [](
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                                        sunrealtype tt,
-#  else
-                                        realtype tt,
-#  endif
-                                        N_Vector yy,
-                                        N_Vector yp,
-                                        void    *user_data) -> int {
+    auto implicit_function_callback = [](SUNDIALS::realtype tt,
+                                         N_Vector           yy,
+                                         N_Vector           yp,
+                                         void              *user_data) -> int {
       Assert(user_data != nullptr, ExcInternalError());
       ARKode<VectorType> &solver =
         *static_cast<ARKode<VectorType> *>(user_data);
@@ -428,16 +418,12 @@ namespace SUNDIALS
         status = ARKStepSetLinearSolver(arkode_mem, sun_linear_solver, nullptr);
         AssertARKode(status);
 
-        auto jacobian_times_vector_callback = [](N_Vector v,
-                                                 N_Vector Jv,
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                                                 sunrealtype t,
-#  else
-                                                 realtype      t,
-#  endif
-                                                 N_Vector y,
-                                                 N_Vector fy,
-                                                 void    *user_data,
+        auto jacobian_times_vector_callback = [](N_Vector           v,
+                                                 N_Vector           Jv,
+                                                 SUNDIALS::realtype t,
+                                                 N_Vector           y,
+                                                 N_Vector           fy,
+                                                 void              *user_data,
                                                  N_Vector) -> int {
           Assert(user_data != nullptr, ExcInternalError());
           ARKode<VectorType> &solver =
@@ -459,15 +445,10 @@ namespace SUNDIALS
             *src_fy);
         };
 
-        auto jacobian_times_vector_setup_callback = [](
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                                                      sunrealtype t,
-#  else
-                                                      realtype t,
-#  endif
-                                                      N_Vector y,
-                                                      N_Vector fy,
-                                                      void *user_data) -> int {
+        auto jacobian_times_vector_setup_callback = [](SUNDIALS::realtype t,
+                                                       N_Vector           y,
+                                                       N_Vector           fy,
+                                                       void *user_data) -> int {
           Assert(user_data != nullptr, ExcInternalError());
           ARKode<VectorType> &solver =
             *static_cast<ARKode<VectorType> *>(user_data);
@@ -490,25 +471,15 @@ namespace SUNDIALS
         AssertARKode(status);
         if (jacobian_preconditioner_solve)
           {
-            auto solve_with_jacobian_callback = [](
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                                                  sunrealtype t,
-#  else
-                                                  realtype   t,
-#  endif
-                                                  N_Vector y,
-                                                  N_Vector fy,
-                                                  N_Vector r,
-                                                  N_Vector z,
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                                                  sunrealtype gamma,
-                                                  sunrealtype delta,
-#  else
-                                                  realtype   gamma,
-                                                  realtype   delta,
-#  endif
-                                                  int   lr,
-                                                  void *user_data) -> int {
+            auto solve_with_jacobian_callback = [](SUNDIALS::realtype t,
+                                                   N_Vector           y,
+                                                   N_Vector           fy,
+                                                   N_Vector           r,
+                                                   N_Vector           z,
+                                                   SUNDIALS::realtype gamma,
+                                                   SUNDIALS::realtype delta,
+                                                   int                lr,
+                                                   void *user_data) -> int {
               Assert(user_data != nullptr, ExcInternalError());
               ARKode<VectorType> &solver =
                 *static_cast<ARKode<VectorType> *>(user_data);
@@ -532,22 +503,13 @@ namespace SUNDIALS
                 lr);
             };
 
-            auto jacobian_solver_setup_callback = [](
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                                                    sunrealtype t,
-#  else
-                                                    realtype t,
-#  endif
-                                                    N_Vector     y,
-                                                    N_Vector     fy,
-                                                    booleantype  jok,
-                                                    booleantype *jcurPtr,
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                                                    sunrealtype gamma,
-#  else
-                                                    realtype gamma,
-#  endif
-                                                    void *user_data) -> int {
+            auto jacobian_solver_setup_callback = [](SUNDIALS::realtype t,
+                                                     N_Vector           y,
+                                                     N_Vector           fy,
+                                                     booleantype        jok,
+                                                     booleantype       *jcurPtr,
+                                                     SUNDIALS::realtype gamma,
+                                                     void *user_data) -> int {
               Assert(user_data != nullptr, ExcInternalError());
               ARKode<VectorType> &solver =
                 *static_cast<ARKode<VectorType> *>(user_data);
@@ -665,13 +627,7 @@ namespace SUNDIALS
         AssertARKode(status);
 
         auto mass_matrix_times_vector_setup_callback =
-          [](
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-            sunrealtype t,
-#  else
-            realtype                                         t,
-#  endif
-            void *mtimes_data) -> int {
+          [](SUNDIALS::realtype t, void *mtimes_data) -> int {
           Assert(mtimes_data != nullptr, ExcInternalError());
           ARKode<VectorType> &solver =
             *static_cast<ARKode<VectorType> *>(mtimes_data);
@@ -680,13 +636,9 @@ namespace SUNDIALS
             solver.mass_times_setup, solver.pending_exception, t);
         };
 
-        auto mass_matrix_times_vector_callback = [](N_Vector v,
-                                                    N_Vector Mv,
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                                                    sunrealtype t,
-#  else
-                                                    realtype t,
-#  endif
+        auto mass_matrix_times_vector_callback = [](N_Vector           v,
+                                                    N_Vector           Mv,
+                                                    SUNDIALS::realtype t,
                                                     void *mtimes_data) -> int {
           Assert(mtimes_data != nullptr, ExcInternalError());
           ARKode<VectorType> &solver =
@@ -713,13 +665,8 @@ namespace SUNDIALS
 
         if (mass_preconditioner_solve)
           {
-            auto mass_matrix_solver_setup_callback = [](
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                                                       sunrealtype t,
-#  else
-                                                       realtype t,
-#  endif
-                                                       void *user_data) -> int {
+            auto mass_matrix_solver_setup_callback =
+              [](SUNDIALS::realtype t, void *user_data) -> int {
               Assert(user_data != nullptr, ExcInternalError());
               ARKode<VectorType> &solver =
                 *static_cast<ARKode<VectorType> *>(user_data);
@@ -728,21 +675,12 @@ namespace SUNDIALS
                 solver.mass_preconditioner_setup, solver.pending_exception, t);
             };
 
-            auto solve_with_mass_matrix_callback = [](
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                                                     sunrealtype t,
-#  else
-                                                     realtype   t,
-#  endif
-                                                     N_Vector r,
-                                                     N_Vector z,
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                                                     sunrealtype delta,
-#  else
-                                                     realtype   delta,
-#  endif
-                                                     int   lr,
-                                                     void *user_data) -> int {
+            auto solve_with_mass_matrix_callback = [](SUNDIALS::realtype t,
+                                                      N_Vector           r,
+                                                      N_Vector           z,
+                                                      SUNDIALS::realtype delta,
+                                                      int                lr,
+                                                      void *user_data) -> int {
               Assert(user_data != nullptr, ExcInternalError());
               ARKode<VectorType> &solver =
                 *static_cast<ARKode<VectorType> *>(user_data);
index ee26b25fa6977bed6018beaa2eef1888d1b98543..6d57319cb7a4c3550fdffa159fbabac907e3ad06 100644 (file)
@@ -227,16 +227,11 @@ namespace SUNDIALS
 
     status = IDAInit(
       ida_mem,
-      [](
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-        sunrealtype tt,
-#  else
-        realtype tt,
-#  endif
-        N_Vector yy,
-        N_Vector yp,
-        N_Vector rr,
-        void    *user_data) -> int {
+      [](SUNDIALS::realtype tt,
+         N_Vector           yy,
+         N_Vector           yp,
+         N_Vector           rr,
+         void              *user_data) -> int {
         IDA<VectorType> &solver = *static_cast<IDA<VectorType> *>(user_data);
 
         auto *src_yy   = internal::unwrap_nvector_const<VectorType>(yy);
@@ -321,7 +316,7 @@ namespace SUNDIALS
 #  if DEAL_II_SUNDIALS_VERSION_LT(6, 0, 0)
     LS = SUNLinSolNewEmpty();
 #  else
-    LS = SUNLinSolNewEmpty(ida_ctx);
+    LS      = SUNLinSolNewEmpty(ida_ctx);
 #  endif
 
     LS->content = this;
@@ -349,14 +344,9 @@ namespace SUNDIALS
                 ExcFunctionNotProvided("solve_with_jacobian"));
     LS->ops->solve = [](SUNLinearSolver LS,
                         SUNMatrix /*ignored*/,
-                        N_Vector x,
-                        N_Vector b,
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                        sunrealtype tol
-#  else
-                        realtype tol
-#  endif
-                        ) -> int {
+                        N_Vector           x,
+                        N_Vector           b,
+                        SUNDIALS::realtype tol) -> int {
       IDA<VectorType> &solver = *static_cast<IDA<VectorType> *>(LS->content);
 
       auto *src_b = internal::unwrap_nvector_const<VectorType>(b);
@@ -391,7 +381,7 @@ namespace SUNDIALS
 #  if DEAL_II_SUNDIALS_VERSION_LT(6, 0, 0)
     J = SUNMatNewEmpty();
 #  else
-    J = SUNMatNewEmpty(ida_ctx);
+    J       = SUNMatNewEmpty(ida_ctx);
 #  endif
     J->content = this;
 
@@ -424,22 +414,16 @@ namespace SUNDIALS
     // calling IDASetLinearSolver
     status = IDASetJacFn(
       ida_mem,
-      [](
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-        sunrealtype tt,
-        sunrealtype cj,
-#  else
-        realtype tt,
-        realtype cj,
-#  endif
-        N_Vector yy,
-        N_Vector yp,
-        N_Vector /* residual */,
-        SUNMatrix /* ignored */,
-        void *user_data,
-        N_Vector /* tmp1 */,
-        N_Vector /* tmp2 */,
-        N_Vector /* tmp3 */) -> int {
+      [](SUNDIALS::realtype tt,
+         SUNDIALS::realtype cj,
+         N_Vector           yy,
+         N_Vector           yp,
+         N_Vector /* residual */,
+         SUNMatrix /* ignored */,
+         void *user_data,
+         N_Vector /* tmp1 */,
+         N_Vector /* tmp2 */,
+         N_Vector /* tmp3 */) -> int {
         Assert(user_data != nullptr, ExcInternalError());
         IDA<VectorType> &solver = *static_cast<IDA<VectorType> *>(user_data);
 
index f80fe1d73c973f4498f8235ee40c64813cd133b4..483844536c0c185148e4ea1da95923f4e1200b0c 100644 (file)
@@ -362,14 +362,9 @@ namespace SUNDIALS
 
         LS->ops->solve = [](SUNLinearSolver LS,
                             SUNMatrix /*ignored*/,
-                            N_Vector x,
-                            N_Vector b,
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                            sunrealtype tol
-#  else
-                            realtype tol
-#  endif
-                            ) -> int {
+                            N_Vector           x,
+                            N_Vector           b,
+                            SUNDIALS::realtype tol) -> int {
           // Receive the object that describes the linear solver and
           // unpack the pointer to the KINSOL object from which we can then
           // get the 'reinit' and 'solve' functions.
@@ -399,7 +394,7 @@ namespace SUNDIALS
 #  if DEAL_II_SUNDIALS_VERSION_LT(6, 0, 0)
         J = SUNMatNewEmpty();
 #  else
-        J = SUNMatNewEmpty(kinsol_ctx);
+        J  = SUNMatNewEmpty(kinsol_ctx);
 #  endif
         J->content = this;
 
index 45101dcf5b76ba781267c87fd092a232f89090e5..d3877f20cac06046784c8326f1acaf91f9c0562e 100644 (file)
@@ -36,6 +36,7 @@
 #  endif
 
 #  include <deal.II/sundials/n_vector.h>
+#  include <deal.II/sundials/sundials_types.h>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -195,14 +196,9 @@ namespace SUNDIALS
     int
     arkode_linsol_solve(SUNLinearSolver LS,
                         SUNMatrix /*ignored*/,
-                        N_Vector x,
-                        N_Vector b,
-#  if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                        sunrealtype tol
-#  else
-                        realtype tol
-#  endif
-    )
+                        N_Vector           x,
+                        N_Vector           b,
+                        SUNDIALS::realtype tol)
     {
       LinearSolverContent<VectorType> *content = access_content<VectorType>(LS);
 
index 701cb3829b5c791ddeedf72798f633f3069c928a..ab54be76c32ebebd759094d3c9e4eb8a7d167fb6 100644 (file)
@@ -95,17 +95,11 @@ main()
   };
 
 
-  ode.jacobian_times_setup = [&](
-#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                               sunrealtype t,
-#else
-                               realtype t,
-#endif
-                               const VectorType &y,
-                               const VectorType &fy) {
-    J       = 0;
-    J(2, 2) = -1.0 / eps;
-  };
+  ode.jacobian_times_setup =
+    [&](SUNDIALS::realtype t, const VectorType &y, const VectorType &fy) {
+      J       = 0;
+      J(2, 2) = -1.0 / eps;
+    };
 
   ode.jacobian_times_vector = [&](const VectorType &v,
                                   VectorType       &Jv,
index 1e715c40f89d6600b91dcd986987fd89e5d35df5..8d9745bd413fc20485e9186bf7c8d8b166d219b1 100644 (file)
@@ -95,17 +95,11 @@ main()
   };
 
 
-  ode.jacobian_times_setup = [&](
-#if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0)
-                               sunrealtype t,
-#else
-                               realtype t,
-#endif
-                               const VectorType &y,
-                               const VectorType &fy) {
-    J       = 0;
-    J(2, 2) = -1.0 / eps;
-  };
+  ode.jacobian_times_setup =
+    [&](SUNDIALS::realtype t, const VectorType &y, const VectorType &fy) {
+      J       = 0;
+      J(2, 2) = -1.0 / eps;
+    };
 
   ode.jacobian_times_vector = [&](const VectorType &v,
                                   VectorType       &Jv,

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.