]> https://gitweb.dealii.org/ - dealii.git/commitdiff
ARKode: remove/deprecate unused functionality 11661/head
authorSebastian Proell <proell@lnm.mw.tum.de>
Mon, 1 Feb 2021 14:52:58 +0000 (15:52 +0100)
committerSebastian Proell <proell@lnm.mw.tum.de>
Tue, 2 Feb 2021 11:24:15 +0000 (12:24 +0100)
remove internal solver reference in operators
deprecate unnecessary reinit_vector

doc/news/changes/incompatibilities/20210201Proell [new file with mode: 0644]
include/deal.II/sundials/arkode.h
source/sundials/arkode.cc
tests/sundials/arkode_01.cc
tests/sundials/arkode_02.cc
tests/sundials/arkode_03.cc
tests/sundials/arkode_04.cc
tests/sundials/arkode_05.cc
tests/sundials/arkode_06.cc
tests/sundials/arkode_07.cc
tests/sundials/arkode_08.cc

diff --git a/doc/news/changes/incompatibilities/20210201Proell b/doc/news/changes/incompatibilities/20210201Proell
new file mode 100644 (file)
index 0000000..503fdaa
--- /dev/null
@@ -0,0 +1,3 @@
+Deprecated: ARKode does no longer need a `reinit_vector` function.
+<br>
+(Sebastian Proell, 2021/02/01)
index bd5bda0dcf80850fd6d83a24b7a4f5ca55a68236..0196c33a81ba17875d92b6fc8e0099e87284d212 100644 (file)
@@ -280,10 +280,8 @@ namespace SUNDIALS
    * $z_j$ where $j < i$). Additional information on the specific predictor
    * algorithms implemented in ARKode is provided in ARKode documentation.
    *
-   * The user has to provide the implementation of the following
-   *`std::function`s:
-   *  - reinit_vector()
-   * and either one or both of
+   * The user has to provide the implementation of at least one (or both) of the
+   * following `std::function`s:
    *  - implicit_function()
    *  - explicit_function()
    *
@@ -363,12 +361,7 @@ namespace SUNDIALS
    *
    * SUNDIALS::ARKode<VectorType> ode;
    *
-   * ode.reinit_vector = [] (VectorType&v)
-   * {
-   *   v.reinit(2);
-   * };
-   *
-   * double kappa = 1.0;
+   * const double kappa = 1.0;
    *
    * ode.explicit_function = [kappa] (double,
    *                                  const VectorType &y,
@@ -634,9 +627,14 @@ namespace SUNDIALS
     get_arkode_memory() const;
 
     /**
-     * A function object that users need to supply and that is intended to
-     * reinit the given vector.
+     * A function object that was used to `reinit` the given vector. Setting
+     * this field does no longer have any effect and all auxiliary vectors are
+     * reinit-ed automatically based on the user-supplied vector in solve_ode().
+     *
+     * @deprecated This function is no longer used and can be safely removed in
+     *   user code.
      */
+    DEAL_II_DEPRECATED
     std::function<void(VectorType &)> reinit_vector;
 
     /**
@@ -1373,19 +1371,12 @@ namespace SUNDIALS
     /**
      * Constructor.
      *
-     * @param solver The ARKode solver that uses this operator
      * @param A_data Data required by @p a_times_fn
      * @param a_times_fn A function pointer to the function that computes A*v
      */
-    SundialsOperator(ARKode<VectorType> &solver,
-                     void *              A_data,
-                     ATimesFn            a_times_fn);
+    SundialsOperator(void *A_data, ATimesFn a_times_fn);
 
   private:
-    /**
-     * Reference to the ARKode object that uses this SundialsOperator.
-     */
-    ARKode<VectorType> &solver;
     /**
      * Data necessary to evaluate a_times_fn.
      */
@@ -1421,23 +1412,14 @@ namespace SUNDIALS
     /**
      * Constructor.
      *
-     * @param solver The ARKode solver that uses this operator
      * @param P_data Data required by @p p_solve_fn
      * @param p_solve_fn A function pointer to the function that computes A*v
      * @param tol Tolerance, that an iterative solver should use to judge
      *   convergence
      */
-    SundialsPreconditioner(ARKode<VectorType> &solver,
-                           void *              P_data,
-                           PSolveFn            p_solve_fn,
-                           double              tol);
+    SundialsPreconditioner(void *P_data, PSolveFn p_solve_fn, double tol);
 
   private:
-    /**
-     * Reference to the ARKode object that uses this SundialsPreconditioner.
-     */
-    ARKode<VectorType> &solver;
-
     /**
      * Data necessary to calls p_solve_fn
      */
index 7cf6d94a4d39c8cc09d7e52e3d168b14f5a7ca8a..41e88ffa2f15401bb74c1d26d91a4e7ea82acd6d 100644 (file)
@@ -384,11 +384,8 @@ namespace SUNDIALS
 
       LinearSolveFunction<VectorType> lsolve;
 
-      //! solver reference to forward all solver related calls to user-specified
-      //! functions
-      ARKode<VectorType> *solver;
-      void *              P_data;
-      void *              A_data;
+      void *P_data;
+      void *A_data;
     };
 
 
@@ -422,18 +419,15 @@ namespace SUNDIALS
                         N_Vector b,
                         realtype tol)
     {
-      auto                content = access_content<VectorType>(LS);
-      ARKode<VectorType> &solver  = *(content->solver);
+      auto content = access_content<VectorType>(LS);
 
       auto *src_b = internal::unwrap_nvector_const<VectorType>(b);
       auto *dst_x = internal::unwrap_nvector<VectorType>(x);
 
-      SundialsOperator<VectorType> op(solver,
-                                      content->A_data,
-                                      content->a_times_fn);
+      SundialsOperator<VectorType> op(content->A_data, content->a_times_fn);
 
       SundialsPreconditioner<VectorType> preconditioner(
-        solver, content->P_data, content->preconditioner_solve, tol);
+        content->P_data, content->preconditioner_solve, tol);
 
       return content->lsolve(op, preconditioner, *dst_x, *src_b, tol);
     }
@@ -503,8 +497,7 @@ namespace SUNDIALS
   class SundialsLinearSolverWrapper
   {
   public:
-    SundialsLinearSolverWrapper(ARKode<VectorType> &            solver,
-                                LinearSolveFunction<VectorType> lsolve)
+    SundialsLinearSolverWrapper(LinearSolveFunction<VectorType> lsolve)
     {
       sun_linear_solver                  = SUNLinSolNewEmpty();
       sun_linear_solver->ops->gettype    = arkode_linsol_get_type;
@@ -515,7 +508,6 @@ namespace SUNDIALS
       sun_linear_solver->ops->setpreconditioner =
         arkode_linsol_set_preconditioner<VectorType>;
 
-      content.solver             = &solver;
       content.lsolve             = lsolve;
       sun_linear_solver->content = &content;
     }
@@ -811,7 +803,7 @@ namespace SUNDIALS
           {
             linear_solver =
               std::make_unique<SundialsLinearSolverWrapper<VectorType>>(
-                *this, solve_linearized_system);
+                solve_linearized_system);
             sun_linear_solver = linear_solver->get_wrapped_solver();
           }
         else
@@ -882,7 +874,7 @@ namespace SUNDIALS
           {
             mass_solver =
               std::make_unique<SundialsLinearSolverWrapper<VectorType>>(
-                *this, solve_mass);
+                solve_mass);
             sun_mass_linear_solver = mass_solver->get_wrapped_solver();
           }
         else
@@ -949,11 +941,9 @@ namespace SUNDIALS
 #  if DEAL_II_SUNDIALS_VERSION_GTE(4, 0, 0)
 
   template <typename VectorType>
-  SundialsOperator<VectorType>::SundialsOperator(ARKode<VectorType> &solver,
-                                                 void *              A_data,
-                                                 ATimesFn            a_times_fn)
-    : solver(solver)
-    , A_data(A_data)
+  SundialsOperator<VectorType>::SundialsOperator(void *   A_data,
+                                                 ATimesFn a_times_fn)
+    : A_data(A_data)
     , a_times_fn(a_times_fn)
 
   {
@@ -978,12 +968,10 @@ namespace SUNDIALS
 
   template <typename VectorType>
   SundialsPreconditioner<VectorType>::SundialsPreconditioner(
-    ARKode<VectorType> &solver,
-    void *              P_data,
-    PSolveFn            p_solve_fn,
-    double              tol)
-    : solver(solver)
-    , P_data(P_data)
+    void *   P_data,
+    PSolveFn p_solve_fn,
+    double   tol)
+    : P_data(P_data)
     , p_solve_fn(p_solve_fn)
     , tol(tol)
   {}
index 28c49b0378ae161caf748c52e7d6530369a71383..279f2ee6ad1960d6e28cee42ca538f0b4b6e42eb 100644 (file)
@@ -76,8 +76,6 @@ main(int argc, char **argv)
 
   SUNDIALS::ARKode<VectorType> ode(data);
 
-  ode.reinit_vector = [&](VectorType &v) { v.reinit(2); };
-
   double kappa = 1.0;
 
   ode.explicit_function =
index afae0d877b7d4a90784e9aa811e6b320436c3d74..8e4d395981c279d2694d9fb843101c2039f88092 100644 (file)
@@ -69,8 +69,6 @@ main(int argc, char **argv)
 
   SUNDIALS::ARKode<VectorType> ode(data);
 
-  ode.reinit_vector = [&](VectorType &v) { v.reinit(2); };
-
   double kappa = 1.0;
 
   ode.implicit_function =
index db0a27ce609d05646caa5f190641fbd91782f849..0714ac9751fe0e48040c424ed9393baf138eeef8 100644 (file)
@@ -72,11 +72,6 @@ main(int argc, char **argv)
 
   SUNDIALS::ARKode<VectorType> ode(data);
 
-  ode.reinit_vector = [&](VectorType &v) {
-    // Three independent variables
-    v.reinit(3);
-  };
-
   // Parameters
   double u0 = 3.9, v0 = 1.1, w0 = 2.8, a = 1.2, b = 2.5, eps = 1e-5;
 
index 3198c8a57e37ee718ff127f8adfdcac37aaeeb80..04116ad06ca1945e94b73b34e7879aa904591fec 100644 (file)
@@ -72,11 +72,6 @@ main(int argc, char **argv)
 
   SUNDIALS::ARKode<VectorType> ode(data);
 
-  ode.reinit_vector = [&](VectorType &v) {
-    // Three independent variables
-    v.reinit(3);
-  };
-
   // Parameters
   double u0 = 3.9, v0 = 1.1, w0 = 2.8, a = 1.2, b = 2.5, eps = 1e-5;
   // Explicit jacobian.
index 40339066bd92a327887693a747e17abf2724c05e..66d37344f3c8a29e7e663ddca639495b3eca310e 100644 (file)
@@ -72,11 +72,6 @@ main(int argc, char **argv)
 
   SUNDIALS::ARKode<VectorType> ode(data);
 
-  ode.reinit_vector = [&](VectorType &v) {
-    // Three independent variables
-    v.reinit(3);
-  };
-
   // Parameters
   double u0 = 3.9, v0 = 1.1, w0 = 2.8, a = 1.2, b = 2.5, eps = 1e-5;
   // Explicit jacobian.
index 8bd82b01ba4a005095e9f7fc2ebbe986bb42996e..d0f5a69810c539aded40203761b89eb31a710984 100644 (file)
@@ -79,11 +79,6 @@ main(int argc, char **argv)
 
   SUNDIALS::ARKode<VectorType> ode(data);
 
-  ode.reinit_vector = [&](VectorType &v) {
-    // Three independent variables
-    v.reinit(3);
-  };
-
   // Parameters
   double u0 = 3.9, v0 = 1.1, w0 = 2.8, a = 1.2, b = 2.5, eps = 1e-5;
   // Explicit jacobian.
index cbbb08378a1c054fd4bcdc74156095c1ffc43820..264099ce4f8d8d57652f7c1b91cf9e95dc820da3 100644 (file)
@@ -79,11 +79,6 @@ main(int argc, char **argv)
 
   SUNDIALS::ARKode<VectorType> ode(data);
 
-  ode.reinit_vector = [&](VectorType &v) {
-    // Three independent variables
-    v.reinit(3);
-  };
-
   // Parameters
   double u0 = 3.9, v0 = 1.1, w0 = 2.8, a = 1.2, b = 2.5, eps = 1e-5;
   // Explicit jacobian.
index 97f8a863b2a5a0fe197e42bac214de16dba8a37d..7abd387302aa9b21ab37a54f01d3d486fea891eb 100644 (file)
@@ -69,8 +69,6 @@ main(int argc, char **argv)
 
   SUNDIALS::ARKode<VectorType> ode(data);
 
-  ode.reinit_vector = [&](VectorType &v) { v.reinit(2); };
-
   // Explicit jacobian = stiffness matrix
   FullMatrix<double> K(2, 2);
   K(0, 0) = K(1, 1) = 0.5;

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.