From: Jean-Paul Pelteret <jppelteret@gmail.com>
Date: Fri, 3 Nov 2017 18:56:01 +0000 (+0100)
Subject: Add support and implement ProductType for Trilinos Sacado::Rad numbers
X-Git-Tag: v9.0.0-rc1~803^2~2
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=67968eaab8e24f9ea4249577030ef43f9c92da72;p=dealii.git

Add support and implement ProductType for Trilinos Sacado::Rad numbers
---

diff --git a/cmake/configure/configure_2_trilinos.cmake b/cmake/configure/configure_2_trilinos.cmake
index f409ba96d2..5c252d69a5 100644
--- a/cmake/configure/configure_2_trilinos.cmake
+++ b/cmake/configure/configure_2_trilinos.cmake
@@ -198,7 +198,11 @@ MACRO(FEATURE_TRILINOS_CONFIGURE_EXTERNAL)
       "Sacado::Fad::DFad<double>"
       "Sacado::Fad::DFad<float>"
       "Sacado::Fad::DFad<Sacado::Fad::DFad<double>>"
-      "Sacado::Fad::DFad<Sacado::Fad::DFad<float>>")
+      "Sacado::Fad::DFad<Sacado::Fad::DFad<float>>"
+      "Sacado::Rad::ADvar<double>"
+      "Sacado::Rad::ADvar<float>"
+      "Sacado::Rad::ADvar<Sacado::Fad::DFad<double>>"
+      "Sacado::Rad::ADvar<Sacado::Fad::DFad<float>>")
 
   IF (TRILINOS_WITH_MPI)
     SET(DEAL_II_EXPAND_EPETRA_VECTOR "LinearAlgebra::EpetraWrappers::Vector")
diff --git a/include/deal.II/differentiation/ad/sacado_product_types.h b/include/deal.II/differentiation/ad/sacado_product_types.h
index 0ff6ce32fc..f5b0cf8d5f 100644
--- a/include/deal.II/differentiation/ad/sacado_product_types.h
+++ b/include/deal.II/differentiation/ad/sacado_product_types.h
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2015 by the deal.II authors
+// Copyright (C) 2015, 2017 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -25,10 +25,21 @@
 
 DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
 #  include <Sacado.hpp>
+// It appears that some versions of Trilinos do not directly or indirectly
+// include all the headers for all forward and reverse Sacado AD types.
+// So we directly include these both here as a precaution.
+// Standard forward AD classes (templated)
+#  include <Sacado_Fad_DFad.hpp>
+// Reverse AD classes (templated)
+#  include <Sacado_trad.hpp>
 DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
 
 DEAL_II_NAMESPACE_OPEN
 
+
+/* -------------- Sacado::Fad::DFad (Differentiation::AD::NumberTypes::[sacado_fad/sacado_fad_fad]) -------------- */
+
+
 namespace internal
 {
 
@@ -71,11 +82,12 @@ namespace internal
   template <typename T, typename U>
   struct ProductTypeImpl<Sacado::Fad::DFad<T>, Sacado::Fad::DFad<U> >
   {
-    typedef Sacado::Fad::DFad<typename ProductType<T,U>::type > type;
+    typedef Sacado::Fad::DFad<typename ProductType<T,U>::type> type;
   };
 
 }
 
+
 template <typename T>
 struct EnableIfScalar<Sacado::Fad::DFad<T> >
 {
@@ -83,6 +95,129 @@ struct EnableIfScalar<Sacado::Fad::DFad<T> >
 };
 
 
+/* -------------- Sacado::Rad::ADvar (Differentiation::AD::NumberTypes::[sacado_rad/sacado_rad_fad]) -------------- */
+
+
+namespace internal
+{
+
+  template <typename T>
+  struct ProductTypeImpl<Sacado::Rad::ADvar<T>, float>
+  {
+    typedef Sacado::Rad::ADvar<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<float, Sacado::Rad::ADvar<T> >
+  {
+    typedef Sacado::Rad::ADvar<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<Sacado::Rad::ADvar<T>, double>
+  {
+    typedef Sacado::Rad::ADvar<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<double, Sacado::Rad::ADvar<T> >
+  {
+    typedef Sacado::Rad::ADvar<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<Sacado::Rad::ADvar<T>, int>
+  {
+    typedef Sacado::Rad::ADvar<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<int, Sacado::Rad::ADvar<T> >
+  {
+    typedef Sacado::Rad::ADvar<T> type;
+  };
+
+  template <typename T, typename U>
+  struct ProductTypeImpl<Sacado::Rad::ADvar<T>, Sacado::Rad::ADvar<U> >
+  {
+    typedef Sacado::Rad::ADvar<typename ProductType<T,U>::type> type;
+  };
+
+  /* --- Sacado::Rad::ADvar: Temporary type --- */
+
+  template <typename T>
+  struct ProductTypeImpl<Sacado::Rad::ADvari<T>, float>
+  {
+    typedef Sacado::Rad::ADvari<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<float, Sacado::Rad::ADvari<T> >
+  {
+    typedef Sacado::Rad::ADvari<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<Sacado::Rad::ADvari<T>, double>
+  {
+    typedef Sacado::Rad::ADvari<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<double, Sacado::Rad::ADvari<T> >
+  {
+    typedef Sacado::Rad::ADvari<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<Sacado::Rad::ADvari<T>, int>
+  {
+    typedef Sacado::Rad::ADvari<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<int, Sacado::Rad::ADvari<T> >
+  {
+    typedef Sacado::Rad::ADvari<T> type;
+  };
+
+  template <typename T, typename U>
+  struct ProductTypeImpl<Sacado::Rad::ADvari<T>, Sacado::Rad::ADvari<U> >
+  {
+    typedef Sacado::Rad::ADvari<typename ProductType<T,U>::type> type;
+  };
+
+  /* --- Sacado::Rad::ADvar: Main and temporary type --- */
+
+  template <typename T, typename U>
+  struct ProductTypeImpl<Sacado::Rad::ADvar<T>, Sacado::Rad::ADvari<U> >
+  {
+    typedef Sacado::Rad::ADvar<typename ProductType<T,U>::type> type;
+  };
+
+  template <typename T, typename U>
+  struct ProductTypeImpl<Sacado::Rad::ADvari<T>, Sacado::Rad::ADvar<U> >
+  {
+    typedef Sacado::Rad::ADvar<typename ProductType<T,U>::type> type;
+  };
+
+}
+
+
+template <typename T>
+struct EnableIfScalar<Sacado::Rad::ADvar<T> >
+{
+  typedef Sacado::Rad::ADvar<T> type;
+};
+
+
+template <typename T>
+struct EnableIfScalar<Sacado::Rad::ADvari<T> >
+{
+  typedef Sacado::Rad::ADvari<T> type;
+};
+
+
 
 /**
  * Compute the scalar product $a:b=\sum_{i,j} a_{ij}b_{ij}$ between two