]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Support both SUNDIALS 2.7.0 and 3.0.0
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Fri, 3 Nov 2017 16:27:49 +0000 (17:27 +0100)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Fri, 3 Nov 2017 16:51:29 +0000 (17:51 +0100)
include/deal.II/base/config.h.in
include/deal.II/sundials/ida.h
source/sundials/arkode.cc
source/sundials/ida.cc
source/sundials/kinsol.cc

index fcf1065a62c94e89b8a551a20c511ee6fd3e0ab2..1074270ed6c52bf9c71f99e2d51c139aaf1d38d9 100644 (file)
     (major)*1000000 + (minor)*10000 + (subminor)*100 + (patch))
 #endif
 
+ /*
+  * SUNDIALS:
+  *
+  *  Note: The following definitions will be set in sundials_config.h,
+  *        so we don't repeat them here.
+  *
+  *  SUNDIALS_VERSION_MAJOR
+  *  SUNDIALS_VERSION_MINOR
+  *  SUNDIALS_VERSION_PATCH
+  */
+
+ #define DEAL_II_SUNDIALS_VERSION_GTE(major,minor,subminor) \
+   ((SUNDIALS_VERSION_MAJOR * 10000 + \
+     SUNDIALS_VERSION_MINOR * 100 + \
+     SUNDIALS_VERSION_SUBMINOR) \
+     >=  \
+     (major)*10000 + (minor)*100 + (subminor))
+
+ #define DEAL_II_SUNDIALS_VERSION_LT(major,minor,subminor) \
+   ((SUNDIALS_VERSION_MAJOR * 10000 + \
+     SUNDIALS_VERSION_MINOR * 100 + \
+     SUNDIALS_VERSION_SUBMINOR) \
+     <  \
+     (major)*10000 + (minor)*100 + (subminor))
+
 /*
  * PETSc:
  *
index b270ebe0dffbb6d2a5d7263336038d2b60b25b2a..5441c4675f773d178e519cccab1324973a8fedb2 100644 (file)
 
 #include <ida/ida.h>
 #include <ida/ida_spils.h>
-#include <ida/ida_spgmr.h>
-#include <ida/ida_spbcgs.h>
-#include <ida/ida_sptfqmr.h>
+#include <sundials/sundials_config.h>
+#if DEAL_II_SUNDIALS_VERSION_LT(3,0,0)
+#  include <ida/ida_spgmr.h>
+#  include <ida/ida_spbcgs.h>
+#  include <ida/ida_sptfqmr.h>
+#endif
 #include <nvector/nvector_serial.h>
 #include <sundials/sundials_math.h>
 #include <sundials/sundials_types.h>
index ff2afbad10a95a5a1e5f8ff97847092471fdec35..33d62dbf70fe9f09da90518193719fdd36cee8ac 100644 (file)
@@ -32,6 +32,8 @@
 #include <deal.II/base/utilities.h>
 #include <deal.II/sundials/copy.h>
 
+#include <sundials/sundials_config.h>
+
 #include <iostream>
 #include <iomanip>
 #include <arkode/arkode_impl.h>
@@ -137,7 +139,9 @@ namespace SUNDIALS
     template<typename VectorType>
     int t_arkode_solve_jacobian(ARKodeMem arkode_mem,
                                 N_Vector b,
+#if DEAL_II_SUNDIALS_VERSION_LT(3,0,0)
                                 N_Vector,
+#endif
                                 N_Vector ycur,
                                 N_Vector fcur)
     {
@@ -186,8 +190,13 @@ namespace SUNDIALS
 
     template<typename VectorType>
     int t_arkode_solve_mass(ARKodeMem arkode_mem,
+#if DEAL_II_SUNDIALS_VERSION_LT(3,0,0)
                             N_Vector b,
-                            N_Vector)
+                            N_Vector
+#else
+                            N_Vector b
+#endif
+                           )
     {
       ARKode<VectorType> &solver = *static_cast<ARKode<VectorType> *>(arkode_mem->ark_user_data);
       GrowingVectorMemory<VectorType> mem;
@@ -424,7 +433,9 @@ namespace SUNDIALS
         if (setup_jacobian)
           {
             ARKode_mem->ark_lsetup = t_arkode_setup_jacobian<VectorType>;
+#if DEAL_II_SUNDIALS_VERSION_LT(3,0,0)
             ARKode_mem->ark_setupNonNull = true;
+#endif
           }
       }
     else
@@ -441,7 +452,9 @@ namespace SUNDIALS
         if (setup_mass)
           {
             ARKode_mem->ark_msetup = t_arkode_setup_mass<VectorType>;
+#if DEAL_II_SUNDIALS_VERSION_LT(3,0,0)
             ARKode_mem->ark_MassSetupNonNull = true;
+#endif
           }
       }
 
index 188a44cc97fba8a6314b3154d4aa834a98422b44..3bbc1eada68b57b6e2d2520d1a8449ede09fb859 100644 (file)
@@ -387,7 +387,9 @@ namespace SUNDIALS
 
     IDA_mem->ida_lsetup = t_dae_lsetup<VectorType>;
     IDA_mem->ida_lsolve = t_dae_solve<VectorType>;
+#if DEAL_II_SUNDIALS_VERSION_LT(3,0,0)
     IDA_mem->ida_setupNonNull = true;
+#endif
 
     status = IDASetMaxOrd(ida_mem, data.maximum_order);
     AssertIDA(status);
index 178940a478283e0c6a10c482608ec3b18f18037d..d807fbfda068d5312f80102c4c5d8615a175152b 100644 (file)
 #include <deal.II/base/utilities.h>
 #include <deal.II/lac/block_vector.h>
 #ifdef DEAL_II_WITH_TRILINOS
-#include <deal.II/lac/trilinos_parallel_block_vector.h>
-#include <deal.II/lac/trilinos_vector.h>
+#  include <deal.II/lac/trilinos_parallel_block_vector.h>
+#  include <deal.II/lac/trilinos_vector.h>
 #endif
 #ifdef DEAL_II_WITH_PETSC
-#include <deal.II/lac/petsc_parallel_block_vector.h>
-#include <deal.II/lac/petsc_parallel_vector.h>
+#  include <deal.II/lac/petsc_parallel_block_vector.h>
+#  include <deal.II/lac/petsc_parallel_vector.h>
 #endif
 #include <deal.II/base/utilities.h>
 #include <deal.II/sundials/copy.h>
 
-#include <kinsol/kinsol_dense.h>
+#include <sundials/sundials_config.h>
+#if DEAL_II_SUNDIALS_VERSION_GTE(3,0,0)
+#  include <sunmatrix/sunmatrix_dense.h>
+#  include <sunlinsol/sunlinsol_dense.h>
+#  include <kinsol/kinsol_direct.h>
+#else
+#  include <kinsol/kinsol_dense.h>
+#endif
 
 #include <iostream>
 #include <iomanip>
@@ -254,12 +261,20 @@ namespace SUNDIALS
         if (setup_jacobian)
           {
             KIN_mem->kin_lsetup = t_kinsol_setup_jacobian<VectorType>;
+#if DEAL_II_SUNDIALS_VERSION_LT(3,0,0)
             KIN_mem->kin_setupNonNull = true;
+#endif
           }
       }
     else
       {
+#if DEAL_II_SUNDIALS_VERSION_GTE(3,0,0)
+        const auto J = SUNDenseMatrix(system_size, system_size);
+        const auto LS = SUNDenseLinearSolver(u_scale, J);
+        status = KINDlsSetLinearSolver(kinsol_mem, LS, J);
+#else
         status = KINDense(kinsol_mem, system_size);
+#endif
         AssertKINSOL(status);
       }
 

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.