]> https://gitweb.dealii.org/ - dealii.git/commitdiff
AD Helpers: Introduce base class for QP-level helper classes 6979/head
authorJean-Paul Pelteret <jppelteret@gmail.com>
Mon, 11 Feb 2019 17:47:04 +0000 (18:47 +0100)
committerJean-Paul Pelteret <jppelteret@gmail.com>
Mon, 11 Feb 2019 18:07:22 +0000 (19:07 +0100)
include/deal.II/differentiation/ad/ad_helpers.h
source/differentiation/ad/ad_helpers.cc
source/differentiation/ad/ad_helpers.inst1.in
source/differentiation/ad/ad_helpers.inst2.in

index 5460362a72bf36625cb16c2cfac6dbd042fc2010..4f36c1dbfe1104d6a9898a868ccc572c3f91c100 100644 (file)
@@ -2617,6 +2617,323 @@ namespace Differentiation
       }
 
     } // namespace internal
+
+
+
+    /**
+     * A base helper class that facilitates the evaluation of point-wise defined
+     * functions. This is the point-wise counterpart of the
+     * ADHelperCellLevelBase class, and was conceived for computations at a
+     * continuum point, or quadrature point, rather than for finite-element
+     * level calculations. That being said, the interface to this and the
+     * derived classes are sufficiently generic that the dependent function(s)
+     * and their argument(s), that are the independent variables, can be
+     * interpreted in any manner that the user may choose.
+     *
+     * As it offers a field-based interface, this class would
+     * typically be used to compute the derivatives of a constitutive law
+     * defined at a quadrature point; however, it may also be used in other
+     * contexts, such as to compute the linearization of a set of local
+     * nonlinear equations.
+     *
+     * @warning ADOL-C does not support the standard threading models used by
+     * deal.II, so this class should @b not be embedded within a multithreaded
+     * function when using ADOL-C number types. It is, however, suitable for use
+     * in both serial and MPI routines.
+     *
+     * @author Jean-Paul Pelteret, 2016, 2017, 2018
+     */
+    template <int                  dim,
+              enum AD::NumberTypes ADNumberTypeCode,
+              typename ScalarType = double>
+    class ADHelperPointLevelFunctionsBase
+      : public ADHelperBase<ADNumberTypeCode, ScalarType>
+    {
+    public:
+      /**
+       * Type definition for the dimension of the associated input and output
+       * tensor types.
+       */
+      static const unsigned int dimension = dim;
+
+      /**
+       * Type definition for the floating point number type that is used in,
+       * and results from, all computations.
+       */
+      using scalar_type =
+        typename ADHelperBase<ADNumberTypeCode, ScalarType>::scalar_type;
+
+      /**
+       * Type definition for the auto-differentiation number type that is used
+       * in all computations.
+       */
+      using ad_type =
+        typename ADHelperBase<ADNumberTypeCode, ScalarType>::ad_type;
+
+      /**
+       * @name Constructor / destructor
+       */
+      //@{
+
+      /**
+       * The constructor for the class.
+       *
+       * @param[in] n_independent_variables The number of independent variables
+       * that will be used in the definition of the functions that it is
+       * desired to compute the sensitivities of. In the computation of
+       * $\mathbf{f}(\mathbf{X})$, this will be the number of inputs
+       * $\mathbf{X}$, i.e. the dimension of the domain space.
+       * @param[in] n_dependent_variables The number of scalar functions to be
+       * defined that will have a sensitivity to the given independent
+       * variables. In the computation of $\mathbf{f}(\mathbf{X})$, this will
+       * be the number of outputs $\mathbf{f}$, i.e. the dimension of the
+       * image space.
+       */
+      ADHelperPointLevelFunctionsBase(
+        const unsigned int n_independent_variables,
+        const unsigned int n_dependent_variables);
+
+      /**
+       * Destructor
+       */
+      virtual ~ADHelperPointLevelFunctionsBase() = default;
+
+      //@}
+
+      /**
+       * @name Independent variables
+       */
+      //@{
+
+      /**
+       * @copydoc ADHelperBase::reset()
+       */
+      virtual void
+      reset(const unsigned int n_independent_variables =
+              dealii::numbers::invalid_unsigned_int,
+            const unsigned int n_dependent_variables =
+              dealii::numbers::invalid_unsigned_int,
+            const bool clear_registered_tapes = true) override;
+
+      /**
+       * Register the complete set of independent variables $\mathbf{X}$.
+       *
+       * @param[in] values A field that defines the values of all independent
+       * variables. When considering taped AD numbers with branching functions,
+       * to avoid potential issues with branch switching it may be a good idea
+       * to choose these values close or equal to those that will be later
+       * evaluated and differentiated around.
+       *
+       * @note The input value type must correspond to this class's @p scalar_type.
+       * Depending on the selected @p ADNumberTypeCode, this may or may not
+       * correspond with the @p ScalarType prescribed as a template argument.
+       *
+       * @note For taped AD numbers, this operation is only valid in recording mode.
+       */
+      void
+      register_independent_variables(const std::vector<scalar_type> &values);
+
+      /**
+       * Register the subset of independent variables
+       * $\mathbf{A} \subset \mathbf{X}$.
+       *
+       * @param[in] value A field that defines a number of independent
+       * variables. When considering taped AD numbers with branching functions,
+       * to avoid potential issues with branch switching it may be a good idea
+       * to choose these values close or equal to those that will be later
+       * evaluated and differentiated around.
+       * @param[in] extractor An extractor associated with the input field
+       * variables. This effectively defines which components of the global set
+       * of independent variables this field is associated with.
+       *
+       * @note The input value type must correspond to this class's @p scalar_type.
+       * Depending on the selected @p ADNumberTypeCode, this may or may not
+       * correspond with the @p ScalarType prescribed as a template argument.
+       *
+       * @note The input extractor must correspond to the input @p ValueType.
+       * So, for example, if a value is a rank-1 tensor
+       * (i.e. of type Tensor<1,dim,scalar_type>), then the extractor must
+       * be an FEValuesExtractors::Vector or FEValuesExtractors::Tensor<1>.
+       *
+       * @note This function may be repeatedly used until a call to
+       * finalize_sensitive_independent_variables() or
+       * get_sensitive_variables() is made.
+       *
+       * @note For taped AD numbers, this operation is only valid in recording mode.
+       */
+      template <typename ValueType, typename ExtractorType>
+      void
+      register_independent_variable(const ValueType &    value,
+                                    const ExtractorType &extractor);
+
+      /**
+       * Return the complete set of independent variables as represented by
+       * auto-differentiable numbers. These are the independent
+       * variables $\mathbf{X}$ at which the dependent values are evaluated
+       * and differentiated.
+       *
+       * This function indicates to the AD library that implements the
+       * auto-differentiable number type that operations performed on these
+       * numbers are to be tracked so they are considered "sensitive"
+       * variables. This is, therefore, the set of variables with which one
+       * would then perform computations, and based on which one can then
+       * extract both the value of the function and its derivatives with the
+       * member functions below. The values of the components of the returned
+       * object are initialized to the values set with
+       * register_independent_variable().
+       *
+       * @return An array of auto-differentiable type numbers.
+       *
+       * @note For taped AD numbers, this operation is only valid in recording mode.
+       */
+      const std::vector<ad_type> &
+      get_sensitive_variables() const;
+
+      /*
+       * Extract a subset of the independent variables as represented by
+       * auto-differentiable numbers. These are the independent
+       * variables $\mathbf{A} \subset \mathbf{X}$ at which the dependent values
+       * are evaluated and differentiated.
+       *
+       * This function indicates to the AD library that implements the
+       * auto-differentiable number type that operations performed on these
+       * numbers are to be tracked so they are considered "sensitive"
+       * variables. This is, therefore, the set of variables with which one
+       * would then perform computations, and based on which one can then
+       * extract both the value of the function and its derivatives with the
+       * member functions below. The values of the components of the returned
+       * object are initialized to the values set with
+       * register_independent_variable().
+       *
+       * @param[in] extractor An extractor associated with the input field
+       * variables. This effectively defines which components of the global set
+       * of independent variables this field is associated with.
+       * @return An object of auto-differentiable type numbers. The return type is
+       * based on the input extractor, and will be either a scalar,
+       * Tensor<1,dim>, Tensor<2,dim>, or SymmetricTensor<2,dim>.
+       *
+       * @note For taped AD numbers, this operation is only valid in recording mode.
+       */
+      template <typename ExtractorType>
+      typename internal::Extractor<dim,
+                                   ExtractorType>::template tensor_type<ad_type>
+      get_sensitive_variables(const ExtractorType &extractor) const;
+
+      //@}
+
+      /**
+       * @name Operations specific to taped mode: Reusing tapes
+       */
+      //@{
+
+      /**
+       * Set the values for the independent variables $\mathbf{X}$.
+       *
+       * @param[in] values A vector that defines the values of all
+       * independent variables.
+       *
+       * @note The input value type must correspond to this class's @p scalar_type.
+       * Depending on the selected @p ADNumberTypeCode, this may or may not
+       * correspond with the @p ScalarType prescribed as a template argument.
+       *
+       * @note If the @p keep_independent_values flag has been set when
+       * ADHelperBase::start_recording_operations() is called then the tape is
+       * immediately usable after creation, and the values of the independent
+       * variables set by register_independent_variables() are those at which
+       * the function is to be evaluated. In this case, a separate call to this
+       * function is not strictly necessary.
+       */
+      void
+      set_independent_variables(const std::vector<scalar_type> &values);
+
+      /**
+       * Set the values for a subset of independent variables
+       * $\mathbf{A} \subset \mathbf{X}$.
+       *
+       * @param[in] value A field that defines the values of a number of
+       * independent variables.
+       * @param[in] extractor An extractor associated with the input field
+       * variables. This effectively defines which components of the global set
+       * of independent variables this field is associated with.
+       *
+       * @note The input value type must correspond to this class's @p scalar_type.
+       * Depending on the selected @p ADNumberTypeCode, this may or may not
+       * correspond with the @p ScalarType prescribed as a template argument.
+       *
+       * @note The input extractor must correspond to the input @p ValueType.
+       * So, for example, if a value is a rank-1 tensor
+       * (i.e. of type Tensor<1,dim,scalar_type>), then the extractor must
+       * be an FEValuesExtractors::Vector or FEValuesExtractors::Tensor<1>.
+       *
+       * @note If the @p keep_independent_values flag has been set when
+       * ADHelperBase::start_recording_operations() is called then the tape is
+       * immediately usable after creation, and the values of the independent
+       * variables set by register_independent_variable() are those at which the
+       * function is to be evaluated. In this case, a separate call to this
+       * function is not strictly necessary.
+       */
+      template <typename ValueType, typename ExtractorType>
+      void
+      set_independent_variable(const ValueType &    value,
+                               const ExtractorType &extractor);
+
+      //@}
+
+    protected:
+      /**
+       * @name Independent variables
+       */
+      //@{
+
+      /**
+       * Set the actual value of the independent variable $X_{i}$.
+       *
+       * @param[in] index The index in the vector of independent variables.
+       * @param[in] symmetric_component Mark whether this index relates to a
+       * component of a field that has a symmetric counterpart
+       * (e.g. if @p index represents an off-diagonal entry in a symmetric
+       * tensor).
+       * @param[in] value The value to set the index'd independent variable to.
+       */
+      void
+      set_sensitivity_value(const unsigned int index,
+                            const bool         symmetric_component,
+                            const scalar_type &value);
+
+      /**
+       * Return whether the @p index'th independent variables is one for which
+       * we must take into account symmetry when extracting their gradient or
+       * Hessian values.
+       */
+      bool
+      is_symmetric_independent_variable(const unsigned int index) const;
+
+      /**
+       * Return the number of independent variables that have been marked as
+       * being components of a symmetric field.
+       */
+      unsigned int
+      n_symmetric_independent_variables() const;
+
+      //@}
+
+    private:
+      /**
+       * @name Independent variables
+       */
+      //@{
+
+      /**
+       * The independent variables for which we must take into account symmetry
+       * when extracting their gradient or Hessian values.
+       */
+      std::vector<bool> symmetric_independent_variables;
+
+      //@}
+
+    }; // class ADHelperPointLevelFunctionsBase
+
+
   } // namespace AD
 } // namespace Differentiation
 
@@ -2675,6 +2992,113 @@ namespace Differentiation
     }
 
 
+
+    /* ----------------- ADHelperPointLevelFunctionsBase ----------------- */
+
+
+
+    template <int                  dim,
+              enum AD::NumberTypes ADNumberTypeCode,
+              typename ScalarType>
+    template <typename ValueType, typename ExtractorType>
+    void
+    ADHelperPointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::
+      register_independent_variable(const ValueType &    value,
+                                    const ExtractorType &extractor)
+    {
+      // This is actually the same thing as the set_independent_variable
+      // function, in the sense that we simply populate our array of independent
+      // values with a meaningful number. However, in this case we need to
+      // double check that we're not registering these variables twice
+#    ifdef DEBUG
+      const std::vector<unsigned int> index_set(
+        internal::extract_field_component_indices<dim>(extractor));
+      for (unsigned int i = 0; i < index_set.size(); ++i)
+        {
+          Assert(
+            this->registered_independent_variable_values[index_set[i]] == false,
+            ExcMessage(
+              "Overlapping indices for independent variables. "
+              "One or more indices associated with the field that "
+              "is being registered as an independent variable have "
+              "already been associated with another field. This suggests "
+              "that the component offsets used to construct their counterpart "
+              "extractors are incompatible with one another. Make sure that "
+              "the first component for each extractor properly takes into "
+              "account the dimensionality of the preceeding fields."));
+        }
+#    endif
+      set_independent_variable(value, extractor);
+    }
+
+
+
+    template <int                  dim,
+              enum AD::NumberTypes ADNumberTypeCode,
+              typename ScalarType>
+    template <typename ValueType, typename ExtractorType>
+    void
+    ADHelperPointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::
+      set_independent_variable(const ValueType &    value,
+                               const ExtractorType &extractor)
+    {
+      const std::vector<unsigned int> index_set(
+        internal::extract_field_component_indices<dim>(extractor));
+      for (unsigned int i = 0; i < index_set.size(); ++i)
+        {
+          set_sensitivity_value(
+            index_set[i],
+            internal::Extractor<dim, ExtractorType>::symmetric_component(i),
+            internal::get_tensor_entry(value, i));
+        }
+    }
+
+
+
+    template <int                  dim,
+              enum AD::NumberTypes ADNumberTypeCode,
+              typename ScalarType>
+    template <typename ExtractorType>
+    typename internal::Extractor<dim, ExtractorType>::template tensor_type<
+      typename ADHelperBase<ADNumberTypeCode, ScalarType>::ad_type>
+    ADHelperPointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::
+      get_sensitive_variables(const ExtractorType &extractor) const
+    {
+      if (ADNumberTraits<ad_type>::is_taped == true)
+        {
+          Assert(this->active_tape_index() !=
+                   Numbers<ad_type>::invalid_tape_index,
+                 ExcMessage("Invalid tape index"));
+        }
+
+      // If necessary, finalize the internally stored vector of
+      // AD numbers that represents the independent variables
+      this->finalize_sensitive_independent_variables();
+      Assert(this->independent_variables.size() ==
+               this->n_independent_variables(),
+             ExcDimensionMismatch(this->independent_variables.size(),
+                                  this->n_independent_variables()));
+
+      const std::vector<unsigned int> index_set(
+        internal::extract_field_component_indices<dim>(extractor));
+      typename internal::Extractor<dim,
+                                   ExtractorType>::template tensor_type<ad_type>
+        out;
+
+      for (unsigned int i = 0; i < index_set.size(); ++i)
+        {
+          const unsigned int index = index_set[i];
+          Assert(index < this->n_independent_variables(), ExcInternalError());
+          Assert(this->registered_independent_variable_values[index] == true,
+                 ExcInternalError());
+          internal::get_tensor_entry(out, i) =
+            this->independent_variables[index];
+        }
+
+      return out;
+    }
+
+
   } // namespace AD
 } // namespace Differentiation
 
index 2c3a9a5a47bd970894f0738802a48b9fbaee5285..7789a475430b8c88434bdb958fdb5779c33ded6e 100644 (file)
@@ -1143,6 +1143,175 @@ namespace Differentiation
     }
 
 
+
+    /* ----------------- ADHelperPointLevelFunctionsBase  ----------------- */
+
+
+
+    template <int                  dim,
+              enum AD::NumberTypes ADNumberTypeCode,
+              typename ScalarType>
+    ADHelperPointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::
+      ADHelperPointLevelFunctionsBase(
+        const unsigned int n_independent_variables,
+        const unsigned int n_dependent_variables)
+      : ADHelperBase<ADNumberTypeCode, ScalarType>(n_independent_variables,
+                                                   n_dependent_variables)
+      , symmetric_independent_variables(n_independent_variables, false)
+    {}
+
+
+
+    template <int                  dim,
+              enum AD::NumberTypes ADNumberTypeCode,
+              typename ScalarType>
+    void
+    ADHelperPointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::reset(
+      const unsigned int n_independent_variables,
+      const unsigned int n_dependent_variables,
+      const bool         clear_registered_tapes)
+    {
+      ADHelperBase<ADNumberTypeCode, ScalarType>::reset(n_independent_variables,
+                                                        n_dependent_variables,
+                                                        clear_registered_tapes);
+
+      const unsigned int new_n_independent_variables =
+        (n_independent_variables != dealii::numbers::invalid_unsigned_int ?
+           n_independent_variables :
+           this->n_independent_variables());
+      symmetric_independent_variables =
+        std::vector<bool>(new_n_independent_variables, false);
+    }
+
+
+
+    template <int                  dim,
+              enum AD::NumberTypes ADNumberTypeCode,
+              typename ScalarType>
+    bool
+    ADHelperPointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::
+      is_symmetric_independent_variable(const unsigned int index) const
+    {
+      Assert(index < symmetric_independent_variables.size(),
+             ExcInternalError());
+      return symmetric_independent_variables[index];
+    }
+
+
+
+    template <int                  dim,
+              enum AD::NumberTypes ADNumberTypeCode,
+              typename ScalarType>
+    unsigned int
+    ADHelperPointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::
+      n_symmetric_independent_variables() const
+    {
+      return std::count(symmetric_independent_variables.begin(),
+                        symmetric_independent_variables.end(),
+                        true);
+    }
+
+
+
+    template <int                  dim,
+              enum AD::NumberTypes ADNumberTypeCode,
+              typename ScalarType>
+    void
+    ADHelperPointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::
+      register_independent_variables(const std::vector<scalar_type> &values)
+    {
+      // This is actually the same thing the set_independent_variable function,
+      // in the sense that we simply populate our array of independent values
+      // with a meaningful number. However, in this case we need to double check
+      // that we're not registering these variables twice
+      Assert(values.size() == this->n_independent_variables(),
+             ExcMessage(
+               "Vector size does not match number of independent variables"));
+      for (unsigned int i = 0; i < this->n_independent_variables(); ++i)
+        {
+          Assert(this->registered_independent_variable_values[i] == false,
+                 ExcMessage("Independent variable value already registered."));
+        }
+      set_independent_variables(values);
+    }
+
+
+
+    template <int                  dim,
+              enum AD::NumberTypes ADNumberTypeCode,
+              typename ScalarType>
+    const std::vector<
+      typename ADHelperPointLevelFunctionsBase<dim,
+                                               ADNumberTypeCode,
+                                               ScalarType>::ad_type> &
+    ADHelperPointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::
+      get_sensitive_variables() const
+    {
+      if (ADNumberTraits<ad_type>::is_taped == true)
+        {
+          Assert(this->active_tape_index() !=
+                   Numbers<ad_type>::invalid_tape_index,
+                 ExcMessage("Invalid tape index"));
+        }
+
+      // Just in case the user has not done so, we repeat the call to
+      // initialize the internally stored vector of AD numbers that
+      // represents the independent variables.
+      this->finalize_sensitive_independent_variables();
+      Assert(this->independent_variables.size() ==
+               this->n_independent_variables(),
+             ExcDimensionMismatch(this->independent_variables.size(),
+                                  this->n_independent_variables()));
+
+      return this->independent_variables;
+    }
+
+
+
+    template <int                  dim,
+              enum AD::NumberTypes ADNumberTypeCode,
+              typename ScalarType>
+    void
+    ADHelperPointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::
+      set_sensitivity_value(const unsigned int index,
+                            const bool         symmetric_component,
+                            const scalar_type &value)
+    {
+      ADHelperBase<ADNumberTypeCode, ScalarType>::set_sensitivity_value(index,
+                                                                        value);
+      Assert(
+        index < this->n_independent_variables(),
+        ExcMessage(
+          "Trying to set the symmetry flag of a non-existent independent variable."));
+      Assert(index < symmetric_independent_variables.size(),
+             ExcInternalError());
+      symmetric_independent_variables[index] = symmetric_component;
+    }
+
+
+
+    template <int                  dim,
+              enum AD::NumberTypes ADNumberTypeCode,
+              typename ScalarType>
+    void
+    ADHelperPointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::
+      set_independent_variables(const std::vector<scalar_type> &values)
+    {
+      if (ADNumberTraits<ad_type>::is_taped == true)
+        {
+          Assert(this->active_tape_index() !=
+                   Numbers<ad_type>::invalid_tape_index,
+                 ExcMessage("Invalid tape index"));
+        }
+      Assert(values.size() == this->n_independent_variables(),
+             ExcMessage(
+               "Vector size does not match number of independent variables"));
+      for (unsigned int i = 0; i < this->n_independent_variables(); ++i)
+        ADHelperBase<ADNumberTypeCode, ScalarType>::set_sensitivity_value(
+          i, values[i]);
+    }
+
+
   } // namespace AD
 } // namespace Differentiation
 
index bf91f88855b7a56105830a0b6919153331ede771..ef44c9398496eb0a73bcfc8516f34d6175f16b45 100644 (file)
@@ -101,3 +101,42 @@ for ()
     \}
     \}
 }
+
+
+for (deal_II_dimension : DIMENSIONS ; number : REAL_SCALARS)
+{
+  namespace Differentiation
+  \{
+  namespace AD
+  \{
+    // -------------------------- ADHelperPointLevelFunctionsBase ----------------------
+
+    template
+    class ADHelperPointLevelFunctionsBase<deal_II_dimension,NumberTypes::adolc_taped,number>;
+
+    template
+    class ADHelperPointLevelFunctionsBase<deal_II_dimension,NumberTypes::adolc_tapeless,number>;
+
+    \}
+    \}
+}
+
+// Instantiations for ADHelpers for which the underlying number type is fixed
+for (deal_II_dimension : DIMENSIONS)
+{
+    namespace Differentiation
+    \{
+    namespace AD
+    \{
+
+    // -------------------------- ADHelperPointLevelFunctionsBase ----------------------
+
+    template
+    class ADHelperPointLevelFunctionsBase<deal_II_dimension,NumberTypes::adolc_taped,typename NumberTraits<double,NumberTypes::adolc_taped>::ad_type>;
+
+    template
+    class ADHelperPointLevelFunctionsBase<deal_II_dimension,NumberTypes::adolc_tapeless,typename NumberTraits<double,NumberTypes::adolc_tapeless>::ad_type>;
+
+    \}
+    \}
+}
index 30c30a054ce39d354e41219747075a6d9d30a154..5fa4b2511e377467f632bc6f5a6d861451211acd 100644 (file)
@@ -149,3 +149,55 @@ for ()
     \}
     \}
 }
+
+
+for (deal_II_dimension : DIMENSIONS ; number : REAL_SCALARS)
+{
+  namespace Differentiation
+  \{
+  namespace AD
+  \{
+
+    // -------------------------- ADHelperPointLevelFunctionsBase ----------------------
+
+    template
+    class ADHelperPointLevelFunctionsBase<deal_II_dimension,NumberTypes::sacado_dfad_dfad,number>;
+
+    template
+    class ADHelperPointLevelFunctionsBase<deal_II_dimension,NumberTypes::sacado_dfad,number>;
+
+    template
+    class ADHelperPointLevelFunctionsBase<deal_II_dimension,NumberTypes::sacado_rad,number>;
+
+    template
+    class ADHelperPointLevelFunctionsBase<deal_II_dimension,NumberTypes::sacado_rad_dfad,number>;
+
+    \}
+    \}
+}
+
+// Instantiations for ADHelpers for which the underlying number type is fixed
+for (deal_II_dimension : DIMENSIONS)
+{
+    namespace Differentiation
+    \{
+    namespace AD
+    \{
+
+    // -------------------------- ADHelperPointLevelFunctionsBase ----------------------
+
+    template
+    class ADHelperPointLevelFunctionsBase<deal_II_dimension,NumberTypes::sacado_dfad_dfad,typename NumberTraits<double,NumberTypes::sacado_dfad_dfad>::ad_type>;
+
+    template
+    class ADHelperPointLevelFunctionsBase<deal_II_dimension,NumberTypes::sacado_dfad,typename NumberTraits<double,NumberTypes::sacado_dfad>::ad_type>;
+
+    template
+    class ADHelperPointLevelFunctionsBase<deal_II_dimension,NumberTypes::sacado_rad,typename NumberTraits<double,NumberTypes::sacado_rad>::ad_type>;
+
+    template
+    class ADHelperPointLevelFunctionsBase<deal_II_dimension,NumberTypes::sacado_rad_dfad,typename NumberTraits<double,NumberTypes::sacado_rad_dfad>::ad_type>;
+
+    \}
+    \}
+}

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.