From 55417b5493eb932acfa10ef756d10f9ed4454681 Mon Sep 17 00:00:00 2001 From: Jean-Paul Pelteret Date: Fri, 3 Nov 2017 14:18:11 +0100 Subject: [PATCH] Implement type traits and helper classes for Sacado numbers --- include/deal.II/differentiation/ad.h | 3 + .../differentiation/ad/sacado_number_types.h | 804 ++++++++++++++++++ .../differentiation/ad/sacado_product_types.h | 6 +- source/differentiation/ad/CMakeLists.txt | 3 + .../differentiation/ad/sacado_number_types.cc | 34 + .../ad/sacado_number_types.inst1.in | 53 ++ .../ad/sacado_number_types.inst2.in | 27 + 7 files changed, 927 insertions(+), 3 deletions(-) create mode 100644 include/deal.II/differentiation/ad/sacado_number_types.h create mode 100644 source/differentiation/ad/sacado_number_types.cc create mode 100644 source/differentiation/ad/sacado_number_types.inst1.in create mode 100644 source/differentiation/ad/sacado_number_types.inst2.in diff --git a/include/deal.II/differentiation/ad.h b/include/deal.II/differentiation/ad.h index 3e7b5f9c62..aa8a7b8df7 100644 --- a/include/deal.II/differentiation/ad.h +++ b/include/deal.II/differentiation/ad.h @@ -26,6 +26,9 @@ #include #include +#include +#include + DEAL_II_NAMESPACE_OPEN /** diff --git a/include/deal.II/differentiation/ad/sacado_number_types.h b/include/deal.II/differentiation/ad/sacado_number_types.h new file mode 100644 index 0000000000..b0cb3ed175 --- /dev/null +++ b/include/deal.II/differentiation/ad/sacado_number_types.h @@ -0,0 +1,804 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2016 - 2017 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#ifndef dealii_differentiation_ad_sacado_number_types_h +#define dealii_differentiation_ad_sacado_number_types_h + +#include + +#ifdef DEAL_II_WITH_TRILINOS + +DEAL_II_DISABLE_EXTRA_DIAGNOSTICS +#include +// 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 +// Reverse AD classes (templated) +#include +DEAL_II_ENABLE_EXTRA_DIAGNOSTICS + +#include +#include + +#include +#include + +#include +#include + +DEAL_II_NAMESPACE_OPEN + + +namespace Differentiation +{ + namespace AD + { + + /** + * A struct to indicate whether a given @p NumberType is a + * Sacado number or not. By default, numbers are not considered to + * have the necessary characteristics to fulfill this condition. + * + * @author Jean-Paul Pelteret, 2017 + */ + template + struct is_sacado_number; + + + /** + * A struct to indicate whether a given @p NumberType is a supported Sacado::Fad + * number or not. By default, numbers are not considered to have the necessary + * characteristics to fulfill this condition. + * + * @author Jean-Paul Pelteret, 2017 + */ + template + struct is_sacado_dfad_number; + + + /** + * A struct to indicate whether a given @p NumberType is a supported Sacado::Rad + * number or not. By default, numbers are not considered to have the necessary + * characteristics to fulfill this condition. + * + * @author Jean-Paul Pelteret, 2017 + */ + template + struct is_sacado_rad_number; + + + namespace internal + { + + /** + * A struct that provides a uniform interface to critical + * implementation details of a @p SacadoNumber. It defines + * various number types, and records how many levels of + * differentiation the number is able to support. + * + * @author Jean-Paul Pelteret, 2017 + */ + template + struct SacadoNumberInfo; + + } + + + + } // namespace AD +} // namespace Differentiation + + +/* --------------------------- inline and template functions and specializations ------------------------- */ + + +#ifndef DOXYGEN + +namespace Differentiation +{ + namespace AD + { + + namespace internal + { + + // The documentation on Sacado numbers is pretty sparse and/or hard to + // navigate. As a point of reference, see + // https://trilinos.org/docs/dev/packages/sacado/doc/html/classSacado_1_1Fad_1_1SimpleFad.html + // for semi-applicable documentation for the Sacado::Fad::Dfad class. + // and the examples in + // https://github.com/trilinos/Trilinos/tree/master/packages/sacado/example + // + // If one dares to venture there, the relevant files for the classes supported here are: + // + // Forward-mode auto-differentiable types: + // https://github.com/trilinos/Trilinos/blob/master/packages/sacado/src/sacado_dfad_DFad.hpp + // https://github.com/trilinos/Trilinos/blob/master/packages/sacado/src/sacado_dfad_GeneralFad.hpp + // + // Reverse-mode auto-differentiable types: + // https://github.com/trilinos/Trilinos/blob/master/packages/sacado/src/Sacado_trad.hpp + + + /** + * Specialization for Sacado::Fad numbers + */ + template + struct SacadoNumberInfo >::value + >::type> + { + typedef SacadoNumber ad_type; + typedef typename ad_type::scalar_type scalar_type; + typedef typename ad_type::value_type value_type; + typedef typename ad_type::value_type derivative_type; + + static const unsigned int n_supported_derivative_levels + = 1 + SacadoNumberInfo::n_supported_derivative_levels; + }; + + + /** + * Specialization for Sacado::Rad numbers + */ + template + struct SacadoNumberInfo >::value + >::type> + { + typedef SacadoNumber ad_type; + typedef typename ad_type::ADVari::scalar_type scalar_type; + typedef typename ad_type::ADVari::value_type value_type; + typedef typename ad_type::ADVari::value_type derivative_type; + + static const unsigned int n_supported_derivative_levels + = 1 + SacadoNumberInfo::n_supported_derivative_levels; + }; + + + /** + * Specialization for floating point numbers. + * + * This is required as a termination point for the recursive + * templates used in the above specializations. + */ + template + struct SacadoNumberInfo::type>::value + >::type> + { + static const unsigned int n_supported_derivative_levels = 0; + }; + + + /** + * A specialization for the information struct for Sacado dynamic forward + * auto-differentiable numbers. + */ + template + struct ADNumberInfoFromEnum::value>::type + > + { + static const bool is_taped = false; + typedef Sacado::Fad::DFad real_type; + typedef typename SacadoNumberInfo::derivative_type derivative_type; + static const unsigned int n_supported_derivative_levels + = SacadoNumberInfo::n_supported_derivative_levels; + }; + + + /** + * A specialization for the information struct for nested Sacado dynamic + * forward auto-differentiable numbers. + */ + template + struct ADNumberInfoFromEnum::value>::type + > + { + static const bool is_taped = false; + typedef Sacado::Fad::DFad< Sacado::Fad::DFad > real_type; + typedef typename SacadoNumberInfo::derivative_type derivative_type; + static const unsigned int n_supported_derivative_levels + = SacadoNumberInfo::n_supported_derivative_levels; + }; + + + /** + * A specialization for the information struct for Sacado dynamic reverse + * auto-differentiable numbers. + */ + template + struct ADNumberInfoFromEnum::value>::type + > + { + static const bool is_taped = false; + typedef Sacado::Rad::ADvar real_type; + typedef typename SacadoNumberInfo::derivative_type derivative_type; + static const unsigned int n_supported_derivative_levels + = SacadoNumberInfo::n_supported_derivative_levels; + }; + + + /** + * A specialization for the information struct for Sacado dynamic nested + * reverse-forward auto-differentiable numbers. + */ + template + struct ADNumberInfoFromEnum::value>::type + > + { + static const bool is_taped = false; + typedef Sacado::Rad::ADvar< Sacado::Fad::DFad > real_type; + typedef typename SacadoNumberInfo::derivative_type derivative_type; + static const unsigned int n_supported_derivative_levels + = SacadoNumberInfo::n_supported_derivative_levels; + }; + + + /** + * Specialization of the marking strategy for Sacado::Fad::DFad + * auto-differentiable numbers + */ + template + struct Marking< Sacado::Fad::DFad > + { + typedef typename SacadoNumberInfo< Sacado::Fad::DFad >::ad_type ad_type; + typedef typename SacadoNumberInfo< Sacado::Fad::DFad >::derivative_type derivative_type; + typedef typename SacadoNumberInfo< Sacado::Fad::DFad >::scalar_type scalar_type; + + /* + * Initialize the state of an independent variable. + */ + static void + independent_variable(const scalar_type &in, + const unsigned int index, + const unsigned int n_independent_variables, + ad_type &out) + { + // It is required that we first initialise the outer number before + // any of the nested ones. + out = ad_type(n_independent_variables, index, in); + + // Initialize potential nested directional derivatives + Marking::independent_variable( + in, index, n_independent_variables, out.val()); + } + + /* + * Initialize the state of a dependent variable. + */ + static void + dependent_variable(ad_type &out, + const ad_type &func) + { + out = func; + } + }; + + + /** + * Specialization of the marking strategy for Sacado::Rad::ADvar + * auto-differentiable numbers. + */ + template + struct Marking< Sacado::Rad::ADvar > + { + typedef typename SacadoNumberInfo< Sacado::Rad::ADvar >::ad_type ad_type; + typedef typename SacadoNumberInfo< Sacado::Rad::ADvar >::derivative_type derivative_type; + typedef typename SacadoNumberInfo< Sacado::Rad::ADvar >::scalar_type scalar_type; + + /* + * Initialize the state of an independent variable. + */ + static void + independent_variable(const scalar_type &in, + const unsigned int index, + const unsigned int n_independent_variables, + ad_type &out) + { + // For Sacado::Rad::ADvar numbers, we have to initialize the + // ADNumber with an already fully-configured value. This means + // that if this nests another ADNumber then the nested number + // must already be setup and ready for use. + + // Initialize potential nested directional derivatives + derivative_type derivative_initializer; + Marking::independent_variable( + in, index, n_independent_variables, + derivative_initializer); + + // Initialize the outer ad_type + out = derivative_initializer; + } + + /* + * Initialize the state of a dependent variable. + */ + static void + dependent_variable(ad_type &out, + const ad_type &func) + { + out = func; + } + }; + + + /** + * A struct to help extract certain information associated with + * Sacado dynamic reverse auto-differentiable numbers. The @p NumberType + * can be either a floating point number or another Sacado type. + * + * @author Jean-Paul Pelteret, 2017 + */ + template + struct ExtractData< Sacado::Fad::DFad > + { + typedef typename SacadoNumberInfo< Sacado::Fad::DFad >::derivative_type derivative_type; + typedef typename SacadoNumberInfo< Sacado::Fad::DFad >::scalar_type scalar_type; + typedef typename SacadoNumberInfo< Sacado::Fad::DFad >::value_type value_type; + + /** + * Extract the real scalar value. + */ + static scalar_type + value (const Sacado::Fad::DFad &x) + { + return ExtractData::value(x.val()); + } + + + /** + * Extract the number of directional derivatives. + */ + static unsigned int + n_directional_derivatives (const Sacado::Fad::DFad &x) + { + return x.size(); + } + + + /** + * Extract the directional derivative in the specified @p direction. + */ + static derivative_type + directional_derivative (const Sacado::Fad::DFad &x, + const unsigned int direction) + { + if (x.hasFastAccess()) + return x.fastAccessDx(direction); + else + return x.dx(direction); + } + }; + + + /** + * A struct to help extract certain information associated with + * Sacado dynamic reverse auto-differentiable numbers. The @p NumberType + * can be either a floating point number or another Sacado type. + * + * @author Jean-Paul Pelteret, 2017 + */ + template + struct ExtractData< Sacado::Rad::ADvar > + { + typedef typename SacadoNumberInfo< Sacado::Rad::ADvar >::derivative_type derivative_type; + typedef typename SacadoNumberInfo< Sacado::Rad::ADvar >::scalar_type scalar_type; + typedef typename SacadoNumberInfo< Sacado::Rad::ADvar >::value_type value_type; + + /** + * Extract the real scalar value. + */ + static scalar_type + value (const Sacado::Rad::ADvar &x) + { + return ExtractData::value(x.val()); + } + + + /** + * Extract the number of directional derivatives. + */ + static unsigned int + n_directional_derivatives (const Sacado::Rad::ADvar &) + { + // There are as many directional derivatives as there are + // independent variables, but each independent variable can + // only return one directional derivative. + return 1; + } + + + /** + * Extract the directional derivative in the specified @p direction. + * + * @note For reverse-mode AD, @p x should represent an independent + * variable with respect to which one wishes to take the derivative + * @p df/dx of a dependent function @p f(x). + */ + static derivative_type + directional_derivative (const Sacado::Rad::ADvar &x, + const unsigned int ) + { + return x.adj(); + } + }; + + } // namespace internal + + + /* -------------- NumberTypes::sacado_dfad -------------- */ + + + /** + * Specialization of the general ADNumberTraits class that + * provides relevant information for auto-differentiable numbers. + * This specialization is for the case where @p ADNumberType is an + * Sacado::Fad::DFad number templated on a floating point type. + */ + template + struct ADNumberTraits >::value + >::type> + : NumberTraits + {}; + + + /** + * Specialization of the general ADNumberTraits class that + * provides relevant information for auto-differentiable numbers. + * This specialization is for the case where @p ADNumberType is an + * complex Sacado::Fad::DFad number templated on a floating point type. + */ + template + struct ADNumberTraits > >::value + >::type> + : NumberTraits, NumberTypes::sacado_dfad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) Sacado number type. + */ + template<> + struct NumberTraits,NumberTypes::sacado_dfad> + : NumberTraits >::scalar_type,NumberTypes::sacado_dfad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) complex Sacado number type. + */ + template<> + struct NumberTraits >,NumberTypes::sacado_dfad> + : NumberTraits > >::scalar_type,NumberTypes::sacado_dfad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) Sacado number type. + */ + template<> + struct NumberTraits,NumberTypes::sacado_dfad> + : NumberTraits >::scalar_type,NumberTypes::sacado_dfad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) complex Sacado number type. + */ + template<> + struct NumberTraits >,NumberTypes::sacado_dfad> + : NumberTraits > >::scalar_type,NumberTypes::sacado_dfad> + {}; + + + /* -------------- NumberTypes::sacado_rad -------------- */ + + + /** + * Specialization of the general ADNumberTraits class that + * provides relevant information for auto-differentiable numbers. + * This specialization is for the case where @p ADNumberType is an + * Sacado::Rad::ADvar number templated on a floating point type. + */ + template + struct ADNumberTraits >::value + >::type> + : NumberTraits + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) Sacado number type. + */ + template<> + struct NumberTraits,NumberTypes::sacado_rad> + : NumberTraits >::scalar_type,NumberTypes::sacado_rad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) Sacado number type. + */ + template<> + struct NumberTraits,NumberTypes::sacado_rad> + : NumberTraits >::scalar_type,NumberTypes::sacado_rad> + {}; + + +#ifdef DEAL_II_TRILINOS_CXX_SUPPORTS_SACADO_COMPLEX_RAD + + + /** + * Specialization of the general ADNumberTraits class that + * provides relevant information for auto-differentiable numbers. + * This specialization is for the case where @p ADNumberType is an + * complex Sacado::Rad::ADvar number templated on a float type. + */ + template + struct ADNumberTraits > >::value + >::type> + : NumberTraits, NumberTypes::sacado_rad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) complex Sacado number type. + */ + template<> + struct NumberTraits >,NumberTypes::sacado_rad> + : NumberTraits > >::scalar_type,NumberTypes::sacado_rad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) complex Sacado number type. + */ + template<> + struct NumberTraits >,NumberTypes::sacado_rad> + : NumberTraits > >::scalar_type,NumberTypes::sacado_rad> + {}; + + +#endif + + + /* -------------- NumberTypes::sacado_dfad_dfad -------------- */ + + /** + * Specialization of the general ADNumberTraits class that + * provides relevant information for auto-differentiable numbers. + * This specialization is for the case where @p ADNumberType is an + * once nested Sacado::Fad::DFad number, with the inner number templated on + * a floating point type. + */ + template + struct ADNumberTraits > >::value + >::type> + : NumberTraits + {}; + + + /** + * Specialization of the general ADNumberTraits class that + * provides relevant information for auto-differentiable numbers. + * This specialization is for the case where @p ADNumberType is an + * complex nested Sacado::Fad::DFad number templated on a floating + * point type. + */ + template + struct ADNumberTraits > > >::value + >::type> + : NumberTraits, NumberTypes::sacado_dfad_dfad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) nested Sacado number type. + */ + template<> + struct NumberTraits >,NumberTypes::sacado_dfad_dfad> + : NumberTraits > >::scalar_type,NumberTypes::sacado_dfad_dfad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) nested complex Sacado number type. + */ + template<> + struct NumberTraits > >,NumberTypes::sacado_dfad_dfad> + : NumberTraits > > >::scalar_type,NumberTypes::sacado_dfad_dfad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) nested Sacado number type. + */ + template<> + struct NumberTraits >,NumberTypes::sacado_dfad_dfad> + : NumberTraits > >::scalar_type,NumberTypes::sacado_dfad_dfad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) nested complex Sacado number type. + */ + template<> + struct NumberTraits > >,NumberTypes::sacado_dfad_dfad> + : NumberTraits > > >::scalar_type,NumberTypes::sacado_dfad_dfad> + {}; + + + /* -------------- NumberTypes::sacado_rad_dfad -------------- */ + + /** + * Specialization of the general ADNumberTraits class that + * provides relevant information for auto-differentiable numbers. + * This specialization is for the case where @p ADNumberType is an + * once nested Sacado::Fad::DFad number, with the inner number templated on + * a floating point type. + */ + template + struct ADNumberTraits > >::value + >::type> + : NumberTraits + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) nested Sacado number type. + */ + template<> + struct NumberTraits >,NumberTypes::sacado_rad_dfad> + : NumberTraits > >::scalar_type,NumberTypes::sacado_rad_dfad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) nested Sacado number type. + */ + template<> + struct NumberTraits >,NumberTypes::sacado_rad_dfad> + : NumberTraits > >::scalar_type,NumberTypes::sacado_rad_dfad> + {}; + + +#ifdef DEAL_II_TRILINOS_CXX_SUPPORTS_SACADO_COMPLEX_RAD + + + /** + * Specialization of the general ADNumberTraits class that + * provides relevant information for auto-differentiable numbers. + * This specialization is for the case where @p ADNumberType is an + * complex nested Sacado::Fad::DFad number templated on a floating + * point type. + */ + template + struct ADNumberTraits > > >::value + >::type> + : NumberTraits, NumberTypes::sacado_rad_dfad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) nested complex Sacado number type. + */ + template<> + struct NumberTraits > >,NumberTypes::sacado_rad_dfad> + : NumberTraits > > >::scalar_type,NumberTypes::sacado_rad_dfad> + {}; + + + /** + * Specialization of the NumberTraits struct for + * the (otherwise disabled) nested complex Sacado number type. + */ + template<> + struct NumberTraits > >,NumberTypes::sacado_rad_dfad> + : NumberTraits > > >::scalar_type,NumberTypes::sacado_rad_dfad> + {}; + + +#endif + + + /* -------------- Additional type traits -------------- */ + + + template + struct is_sacado_dfad_number + : std::false_type + {}; + + + template + struct is_sacado_dfad_number::type>::type_code == NumberTypes::sacado_dfad || + ADNumberTraits::type>::type_code == NumberTypes::sacado_dfad_dfad + >::type> + : std::true_type + {}; + + + template + struct is_sacado_rad_number + : std::false_type + {}; + + + template + struct is_sacado_rad_number::type>::type_code == NumberTypes::sacado_rad || + ADNumberTraits::type>::type_code == NumberTypes::sacado_rad_dfad + >::type> + : std::true_type + {}; + + + template + struct is_sacado_number + : std::false_type + {}; + + + template + struct is_sacado_number::value || + is_sacado_rad_number::value + >::type> + : std::true_type + {}; + + } +} + +#endif // DOXYGEN + + + +DEAL_II_NAMESPACE_CLOSE + + +#endif // DEAL_II_WITH_TRILINOS + +#endif diff --git a/include/deal.II/differentiation/ad/sacado_product_types.h b/include/deal.II/differentiation/ad/sacado_product_types.h index f5b0cf8d5f..756c53d571 100644 --- a/include/deal.II/differentiation/ad/sacado_product_types.h +++ b/include/deal.II/differentiation/ad/sacado_product_types.h @@ -24,14 +24,14 @@ #ifdef DEAL_II_WITH_TRILINOS DEAL_II_DISABLE_EXTRA_DIAGNOSTICS -# include +#include // 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 +#include // Reverse AD classes (templated) -# include +#include DEAL_II_ENABLE_EXTRA_DIAGNOSTICS DEAL_II_NAMESPACE_OPEN diff --git a/source/differentiation/ad/CMakeLists.txt b/source/differentiation/ad/CMakeLists.txt index d445d49661..babc220ee5 100644 --- a/source/differentiation/ad/CMakeLists.txt +++ b/source/differentiation/ad/CMakeLists.txt @@ -17,10 +17,13 @@ INCLUDE_DIRECTORIES(BEFORE ${CMAKE_CURRENT_BINARY_DIR}) SET(_src adolc_number_types.cc + sacado_number_types.cc ) SET(_inst adolc_number_types.inst.in + sacado_number_types.inst1.in + sacado_number_types.inst2.in ) FILE(GLOB _header diff --git a/source/differentiation/ad/sacado_number_types.cc b/source/differentiation/ad/sacado_number_types.cc new file mode 100644 index 0000000000..b5ddd93a27 --- /dev/null +++ b/source/differentiation/ad/sacado_number_types.cc @@ -0,0 +1,34 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +#include + +#ifdef DEAL_II_WITH_TRILINOS + +#include + +DEAL_II_NAMESPACE_OPEN + +/*---------------------- Explicit Instantiations ----------------------*/ + +#include "sacado_number_types.inst1" +#ifdef DEAL_II_TRILINOS_CXX_SUPPORTS_SACADO_COMPLEX_RAD +#include "sacado_number_types.inst2" +#endif + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/source/differentiation/ad/sacado_number_types.inst1.in b/source/differentiation/ad/sacado_number_types.inst1.in new file mode 100644 index 0000000000..d297c16dfd --- /dev/null +++ b/source/differentiation/ad/sacado_number_types.inst1.in @@ -0,0 +1,53 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +for (Number : DIFFERENTIABLE_TRILINOS_SACADO_REAL_SCALARS) +{ + namespace Differentiation + \{ + namespace AD + \{ + template struct ADNumberTraits; + \} + \} +} + + +for (Number : REAL_SCALARS) +{ + namespace Differentiation + \{ + namespace AD + \{ + template struct NumberTraits; + template struct NumberTraits; + template struct NumberTraits; + template struct NumberTraits; + \} + \} +} + + +for (Number : COMPLEX_SCALARS) +{ + namespace Differentiation + \{ + namespace AD + \{ + template struct NumberTraits; + template struct NumberTraits; + \} + \} +} diff --git a/source/differentiation/ad/sacado_number_types.inst2.in b/source/differentiation/ad/sacado_number_types.inst2.in new file mode 100644 index 0000000000..580a215633 --- /dev/null +++ b/source/differentiation/ad/sacado_number_types.inst2.in @@ -0,0 +1,27 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +for (Number : COMPLEX_SCALARS) +{ + namespace Differentiation + \{ + namespace AD + \{ + template struct NumberTraits; + template struct NumberTraits; + \} + \} +} -- 2.39.5