SD BatchOptimizer: Implement cached evaluation

diff --git a/doc/news/changes/minor/20210313Jean-PaulPelteret-4 b/doc/news/changes/minor/20210313Jean-PaulPelteret-4
new file mode 100644 (file)
index 0000000..94b846c
--- /dev/null
+++ b/doc/news/changes/minor/20210313Jean-PaulPelteret-4
@@ -0,0 +1,8 @@
+New: The new method Differentiation::SD::BatchOptimizer::extract() allows one to
+extract results from previous evaluation of the symbolic expressions. 
+This means that one can safely use a single instance of a batch optimizer to
+pre-compute a set of results, cache them, and then later extract some results
+from the cached data without re-evaluating any symbolic expressions.
+<br>
+(Jean-Paul Pelteret, 2021/03/13)
+

diff --git a/include/deal.II/differentiation/sd/symengine_optimizer.h b/include/deal.II/differentiation/sd/symengine_optimizer.h
index 440f15761b10f7338500b717e91362bf4be59359..55dff430ad314e07533063e3c414dd034fbd7da9 100644 (file)
--- a/include/deal.II/differentiation/sd/symengine_optimizer.h
+++ b/include/deal.II/differentiation/sd/symengine_optimizer.h
@@ -1124,6 +1124,11 @@ namespace Differentiation
        * @tparam TensorType The type of tensor to be evaluated and returned
        *         (i.e. Tensor or SymmetricTensor).
        * @param[in] symbol_tensor The symbolic tensor that is to be evaluated.
+       * @param[in] cached_evaluation A vector that stores the numerical
+       *            values of all dependent variables referenced by the
+       *            @p optimizer. This vector is, most typically, first attained
+       *            by a call to the BatchOptimizer::evaluate() variant that
+       *            takes no arguments.
        * @param[in] optimizer The optimizer that can evaluate the input
        *       @p symbol_tensor.
        * @return TensorType<rank, dim, NumberType> The numeric result that the
@@ -1136,6 +1141,7 @@ namespace Differentiation
       TensorType<rank, dim, NumberType>
       tensor_evaluate_optimized(
         const TensorType<rank, dim, Expression> &symbol_tensor,
+        const std::vector<NumberType> &          cached_evaluation,
         const BatchOptimizer<NumberType> &       optimizer)
       {
         TensorType<rank, dim, NumberType> out;
@@ -1143,7 +1149,8 @@ namespace Differentiation
           {
             const TableIndices<rank> indices(
               out.unrolled_to_component_indices(i));
-            out[indices] = optimizer.evaluate(symbol_tensor[indices]);
+            out[indices] =
+              optimizer.extract(symbol_tensor[indices], cached_evaluation);
           }
         return out;
       }
@@ -1161,6 +1168,11 @@ namespace Differentiation
        * @tparam TensorType The type of tensor to be evaluated and returned
        *         (i.e. Tensor or SymmetricTensor).
        * @param[in] symbol_tensor The symbolic tensor that is to be evaluated.
+       * @param[in] cached_evaluation A vector that stores the numerical
+       *            values of all dependent variables referenced by the
+       *            @p optimizer. This vector is, most typically, first attained
+       *            by a call to the BatchOptimizer::evaluate() variant that
+       *            takes no arguments.
        * @param[in] optimizer The optimizer that can evaluate the input
        *       @p symbol_tensor.
        * @return TensorType<rank, dim, NumberType> The numeric result that the
@@ -1170,6 +1182,7 @@ namespace Differentiation
       SymmetricTensor<4, dim, NumberType>
       tensor_evaluate_optimized(
         const SymmetricTensor<4, dim, Expression> &symbol_tensor,
+        const std::vector<NumberType> &            cached_evaluation,
         const BatchOptimizer<NumberType> &         optimizer)
       {
         SymmetricTensor<4, dim, NumberType> out;
@@ -1182,7 +1195,8 @@ namespace Differentiation
             {
               const TableIndices<4> indices =
                 make_rank_4_tensor_indices<dim>(i, j);
-              out[indices] = optimizer.evaluate(symbol_tensor[indices]);
+              out[indices] =
+                optimizer.extract(symbol_tensor[indices], cached_evaluation);
             }
         return out;
       }
@@ -1855,7 +1869,9 @@ namespace Differentiation
       /**
        * Returns the result of a value substitution into the optimized
        * counterpart of all dependent functions. This function fetches all of
-       * those cached values.
+       * those cached values. If these results are stored elsewhere, and you
+       * want to use this optimizer to extract components of the cached results
+       * then you can make use of the extract() functions listed below.
        *
        * These values were computed by substituting a @p substitution_values map
        * during substitute() call.
@@ -1921,6 +1937,56 @@ namespace Differentiation
       SymmetricTensor<rank, dim, ReturnType>
       evaluate(const SymmetricTensor<rank, dim, Expression> &funcs) const;
 
+
+      /**
+       * Returns the result of a value substitution into the optimized
+       * counterpart of @p func. This function fetches that one value
+       * stored in the @p cached_evaluation vector. This vector is, most
+       * typically, first attained by a call to the evaluate() variant that
+       * takes no arguments.
+       */
+      ReturnType
+      extract(const Expression &             func,
+              const std::vector<ReturnType> &cached_evaluation) const;
+
+
+      /**
+       * Returns the result of a value substitution into the optimized
+       * counterpart of @p funcs. This function fetches that subset of
+       * values stored in the @p cached_evaluation vector. This vector is, most
+       * typically, first attained by a call to the evaluate() variant that
+       * takes no arguments.
+       */
+      std::vector<ReturnType>
+      extract(const std::vector<Expression> &funcs,
+              const std::vector<ReturnType> &cached_evaluation) const;
+
+
+      /**
+       * Returns the result of a tensor value substitution into the optimized
+       * counterpart of @p funcs. This function fetches those tensor components
+       * stored in the @p cached_evaluation vector. This vector is, most
+       * typically, first attained by a call to the evaluate() variant that
+       * takes no arguments.
+       */
+      template <int rank, int dim>
+      Tensor<rank, dim, ReturnType>
+      extract(const Tensor<rank, dim, Expression> &funcs,
+              const std::vector<ReturnType> &      cached_evaluation) const;
+
+
+      /**
+       * Returns the result of a tensor value substitution into the optimized
+       * counterpart of @p funcs. This function fetches those symmetric
+       * tensor components stored in the @p cached_evaluation vector. This vector
+       * is, most typically, first attained by a call to the evaluate() variant
+       * that takes no arguments.
+       */
+      template <int rank, int dim>
+      SymmetricTensor<rank, dim, ReturnType>
+      extract(const SymmetricTensor<rank, dim, Expression> &funcs,
+              const std::vector<ReturnType> &cached_evaluation) const;
+
       //@}
 
     private:
@@ -2445,8 +2511,9 @@ namespace Differentiation
     template <typename ReturnType>
     template <int rank, int dim>
     Tensor<rank, dim, ReturnType>
-    BatchOptimizer<ReturnType>::evaluate(
-      const Tensor<rank, dim, Expression> &funcs) const
+    BatchOptimizer<ReturnType>::extract(
+      const Tensor<rank, dim, Expression> &funcs,
+      const std::vector<ReturnType> &      cached_evaluation) const
     {
       Assert(
         values_substituted() == true,
@@ -2454,16 +2521,30 @@ namespace Differentiation
           "The optimizer is not configured to perform evaluation. "
           "This action can only performed after substitute() has been called."));
 
-      return internal::tensor_evaluate_optimized(funcs, *this);
+      return internal::tensor_evaluate_optimized(funcs,
+                                                 cached_evaluation,
+                                                 *this);
     }
 
 
 
     template <typename ReturnType>
     template <int rank, int dim>
-    SymmetricTensor<rank, dim, ReturnType>
+    Tensor<rank, dim, ReturnType>
     BatchOptimizer<ReturnType>::evaluate(
-      const SymmetricTensor<rank, dim, Expression> &funcs) const
+      const Tensor<rank, dim, Expression> &funcs) const
+    {
+      return extract(funcs, dependent_variables_output);
+    }
+
+
+
+    template <typename ReturnType>
+    template <int rank, int dim>
+    SymmetricTensor<rank, dim, ReturnType>
+    BatchOptimizer<ReturnType>::extract(
+      const SymmetricTensor<rank, dim, Expression> &funcs,
+      const std::vector<ReturnType> &               cached_evaluation) const
     {
       Assert(
         values_substituted() == true,
@@ -2471,7 +2552,20 @@ namespace Differentiation
           "The optimizer is not configured to perform evaluation. "
           "This action can only performed after substitute() has been called."));
 
-      return internal::tensor_evaluate_optimized(funcs, *this);
+      return internal::tensor_evaluate_optimized(funcs,
+                                                 cached_evaluation,
+                                                 *this);
+    }
+
+
+
+    template <typename ReturnType>
+    template <int rank, int dim>
+    SymmetricTensor<rank, dim, ReturnType>
+    BatchOptimizer<ReturnType>::evaluate(
+      const SymmetricTensor<rank, dim, Expression> &funcs) const
+    {
+      return extract(funcs, dependent_variables_output);
     }
 
 #  endif // DOXYGEN

diff --git a/source/differentiation/sd/symengine_optimizer.cc b/source/differentiation/sd/symengine_optimizer.cc
index e90516e92873a87c5369545f633aca553e677ea9..09b2c2b74ca5006e1f0f0108d80003d0d9216ce5 100644 (file)
--- a/source/differentiation/sd/symengine_optimizer.cc
+++ b/source/differentiation/sd/symengine_optimizer.cc
@@ -523,7 +523,9 @@ namespace Differentiation
 
     template <typename ReturnType>
     ReturnType
-    BatchOptimizer<ReturnType>::evaluate(const Expression &func) const
+    BatchOptimizer<ReturnType>::extract(
+      const Expression &             func,
+      const std::vector<ReturnType> &cached_evaluation) const
     {
       Assert(
         values_substituted() == true,
@@ -613,27 +615,47 @@ namespace Differentiation
              ExcMessage("Function has not been registered."));
       Assert(it->second < n_dependent_variables(), ExcInternalError());
 
-      return dependent_variables_output[it->second];
+      return cached_evaluation[it->second];
+    }
+
+
+
+    template <typename ReturnType>
+    ReturnType
+    BatchOptimizer<ReturnType>::evaluate(const Expression &func) const
+    {
+      return extract(func, dependent_variables_output);
     }
 
 
 
     template <typename ReturnType>
     std::vector<ReturnType>
-    BatchOptimizer<ReturnType>::evaluate(
-      const std::vector<Expression> &funcs) const
+    BatchOptimizer<ReturnType>::extract(
+      const std::vector<Expression> &funcs,
+      const std::vector<ReturnType> &cached_evaluation) const
     {
       std::vector<ReturnType> out;
       out.reserve(funcs.size());
 
       for (const auto &func : funcs)
-        out.emplace_back(evaluate(func));
+        out.emplace_back(extract(func, cached_evaluation));
 
       return out;
     }
 
 
 
+    template <typename ReturnType>
+    std::vector<ReturnType>
+    BatchOptimizer<ReturnType>::evaluate(
+      const std::vector<Expression> &funcs) const
+    {
+      return extract(funcs, dependent_variables_output);
+    }
+
+
+
     template <typename ReturnType>
     bool
     BatchOptimizer<ReturnType>::is_valid_nonunique_dependent_variable(