From 56913d7e791b248ccb3887cd11efe94344cb6e93 Mon Sep 17 00:00:00 2001 From: Jean-Paul Pelteret Date: Fri, 10 May 2019 22:29:34 +0200 Subject: [PATCH] Add dummy AD::NumberTypes::none for floating point number types These number types are scalar arithmetic types, that hold no derivatives. They are implemented primarily to facilitate the use of template meta-programming techniques to switch between different AD types. This covers the case when the user does not want an AD type at all, but rather a primitive type. --- .../differentiation/ad/ad_number_traits.h | 357 ++++++++++++++++-- .../differentiation/ad/ad_number_types.h | 16 +- 2 files changed, 346 insertions(+), 27 deletions(-) diff --git a/include/deal.II/differentiation/ad/ad_number_traits.h b/include/deal.II/differentiation/ad/ad_number_traits.h index 4675c493bf..8cf3c9866b 100644 --- a/include/deal.II/differentiation/ad/ad_number_traits.h +++ b/include/deal.II/differentiation/ad/ad_number_traits.h @@ -17,6 +17,7 @@ #define dealii_differentiation_ad_ad_number_traits_h #include +#include #include @@ -334,6 +335,24 @@ namespace Differentiation {}; + /** + * A specialization for the information struct for arithmetic / floating + * point numbers. + */ + template + struct ADNumberInfoFromEnum< + ScalarType, + Differentiation::AD::NumberTypes::none, + typename std::enable_if< + std::is_floating_point::value>::type> + { + static const bool is_taped = false; + using real_type = ScalarType; + using derivative_type = ScalarType; + static const unsigned int n_supported_derivative_levels = 0; + }; + + /** * A dummy specialization for floating point numbers. This is helpful * for nested auto-differentiable numbers, where a recursive marking @@ -1068,50 +1087,338 @@ namespace Differentiation # endif + /** * A dummy specialization for floating point numbers. This is necessary to - * deal with the tricky case of higher-order derivative extraction from - * AD numbers. + * deal with the special case of the ADNumberTypeCode that represents a + * floating point number. * - * Sacado nests the directional derivatives within Sacado numbers, while - * ADOL-C does not. So, starting off from a scalar function f(x), the number - * type resulting from the computation of a directional derivative f'(x) of - * that function is a floating point number for ADOL-C number types and - * Sacado::FAD, but that of a Sacado::FAD> is a - * Sacado::FAD. - * - * For this reason, when trying to extract higher-order derivatives of a - * number that does not support them the input value to this function may be - * a scalar type. + * This specific version specializes the generic case where the template + * @p ScalarType represents a floating point type, or complex number that + * is templated on a floating point type. */ template - struct ADNumberTraits< + struct NumberTraits< ScalarType, - typename std::enable_if::value>::type> + NumberTypes::none, + typename std::enable_if< + std::is_floating_point::value || + (boost::is_complex::value && + std::is_floating_point::type>::value)>::type> { + /** + * The internal number type code. + */ + static constexpr enum NumberTypes type_code = NumberTypes::none; + + // The clang compiler does not seem to like these + // variables being defined as constant expressions + // (the tests /ad_number_traits_02 will + // fail with linking errors). However, GCC complains + // about the use of non-constant expressions in + // std::conditional. +# ifdef __clang__ + + /** + * A flag to indicate whether the number is of + * the taped variety or not. + */ + static const bool is_taped; + + + /** + * A flag to indicate whether the number is of + * the tapeless variety or not. + */ + static const bool is_tapeless; + + + /** + * A flag to indicate whether the number represents + * a real value. + */ + static const bool is_real_valued; + + + /** + * A flag to indicate whether the number represents + * a complex value. + */ + static const bool is_complex_valued; + + + /** + * The number of directional derivatives that can be + * taken with this (floating point) number. + */ + static const unsigned int n_supported_derivative_levels; + +# else + + /** + * A flag to indicate whether the number is of + * the taped variety or not. + */ + static constexpr bool is_taped = false; + + + /** + * A flag to indicate whether the number is of + * the tapeless variety or not. + */ + static constexpr bool is_tapeless = false; + + + /** + * A flag to indicate whether the number represents + * a real value. + */ + static constexpr bool is_real_valued = + (!boost::is_complex::value); + + + /** + * A flag to indicate whether the number represents + * a complex value. + */ + static constexpr bool is_complex_valued = !is_real_valued; + + + /** + * The number of directional derivatives that can be + * taken with this (floating point) number. + */ + static constexpr unsigned int n_supported_derivative_levels = 0; + +# endif + + /** * Underlying floating point value type. * This could real-valued or complex-valued. */ using scalar_type = ScalarType; - static ScalarType - get_directional_derivative(const ScalarType & /*x*/, - const unsigned int /*direction*/) + + /** + * Type for real numbers. + */ + using real_type = + typename dealii::numbers::NumberTraits::real_type; + + + /** + * Type for complex numbers. + */ + using complex_type = std::complex; + + + /** + * The actual auto-differentiable number type. + */ + using ad_type = ScalarType; + + /** + * The actual auto-differentiable number directional derivative type. + */ + using derivative_type = ScalarType; + + + /** + * Extract the value of an auto-differentiable number. + */ + static scalar_type + get_scalar_value(const ad_type &x) + { + return x; + } + + + /** + * Extract the derivative value of an auto-differentiable number. + */ + static derivative_type + get_directional_derivative(const ad_type &, const unsigned int) { - // If the AD drivers are correctly implemented then we should not get - // here. This is essentially a dummy for when the ADNumberTypeCode for - // the original AD number (from which one is getting a derivative >= 2) - // is one that specified ADOL-C taped and tapeless numbers, or a - // non-nested Sacado number. - AssertThrow( + Assert( false, ExcMessage( - "Floating point numbers have no directional derivatives.")); - return 0.0; + "Floating point/arithmetic numbers have no directional derivatives.")); + return derivative_type(); + } + + + /** + * Extract the number of directional derivatives value tracked by + * an auto-differentiable number. + */ + static unsigned int + n_directional_derivatives(const ad_type &) + { + Assert( + false, + ExcMessage( + "Floating point/arithmetic numbers have no directional derivatives.")); + return 0; } }; +# ifdef __clang__ + + template + const bool NumberTraits< + ScalarType, + NumberTypes::none, + typename std::enable_if< + std::is_floating_point::value || + (boost::is_complex::value && + std::is_floating_point::type>::value)>::type>::is_taped = false; + + + template + const bool NumberTraits< + ScalarType, + NumberTypes::none, + typename std::enable_if< + std::is_floating_point::value || + (boost::is_complex::value && + std::is_floating_point::type>::value)>::type>::is_tapeless = false; + + + template + const bool NumberTraits< + ScalarType, + NumberTypes::none, + typename std::enable_if< + std::is_floating_point::value || + (boost::is_complex::value && + std::is_floating_point::type>::value)>::type>::is_real_valued = + (!boost::is_complex::value); + + + template + const bool NumberTraits< + ScalarType, + NumberTypes::none, + typename std::enable_if< + std::is_floating_point::value || + (boost::is_complex::value && + std::is_floating_point::type>::value)>::type>::is_complex_valued = + !(NumberTraits::is_real_valued); + + + template + const unsigned int NumberTraits< + ScalarType, + NumberTypes::none, + typename std::enable_if< + std::is_floating_point::value || + (boost::is_complex::value && + std::is_floating_point::type>::value)>::type>::n_supported_derivative_levels = + 0; + +# else + + template + constexpr bool NumberTraits< + ScalarType, + NumberTypes::none, + typename std::enable_if< + std::is_floating_point::value || + (boost::is_complex::value && + std::is_floating_point::type>::value)>::type>::is_taped; + + + template + constexpr bool NumberTraits< + ScalarType, + NumberTypes::none, + typename std::enable_if< + std::is_floating_point::value || + (boost::is_complex::value && + std::is_floating_point::type>::value)>::type>::is_tapeless; + + + template + constexpr bool NumberTraits< + ScalarType, + NumberTypes::none, + typename std::enable_if< + std::is_floating_point::value || + (boost::is_complex::value && + std::is_floating_point::type>::value)>::type>::is_real_valued; + + + template + constexpr bool NumberTraits< + ScalarType, + NumberTypes::none, + typename std::enable_if< + std::is_floating_point::value || + (boost::is_complex::value && + std::is_floating_point::type>::value)>::type>::is_complex_valued; + + + template + constexpr unsigned int NumberTraits< + ScalarType, + NumberTypes::none, + typename std::enable_if< + std::is_floating_point::value || + (boost::is_complex::value && + std::is_floating_point::type>::value)>::type>::n_supported_derivative_levels; + +# endif + + + /** + * A dummy specialization for floating point numbers. + * + * This is necessary to deal with two situations: + * 1. when the dummy "none" number type is to be used, and + * 2. the tricky case of higher-order derivative extraction from AD numbers. + * + * @note Sacado nests the directional derivatives within Sacado numbers, while + * ADOL-C does not. So, starting off from a scalar function f(x), the number + * type resulting from the computation of a directional derivative f'(x) of + * that function is a floating point number for ADOL-C number types and + * Sacado::FAD, but that of a Sacado::FAD> is a + * Sacado::FAD. + * + * For this reason, when trying to extract higher-order derivatives of a + * number that does not support them the input value to this function may be + * a scalar type. + */ + template + struct ADNumberTraits< + ScalarType, + typename std::enable_if::value>::type> + : NumberTraits + {}; + + /** + * A dummy specialization for complex floating point numbers. + */ + template + struct ADNumberTraits< + ComplexScalarType, + typename std::enable_if< + boost::is_complex::value && + std::is_floating_point::value>:: + type> : NumberTraits + {}; + } // namespace AD } // namespace Differentiation diff --git a/include/deal.II/differentiation/ad/ad_number_types.h b/include/deal.II/differentiation/ad/ad_number_types.h index e199ec25bf..f20d1b4130 100644 --- a/include/deal.II/differentiation/ad/ad_number_types.h +++ b/include/deal.II/differentiation/ad/ad_number_types.h @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2016 - 2017 by the deal.II authors +// Copyright (C) 2016 - 2017, 2019 by the deal.II authors // // This file is part of the deal.II library. // @@ -31,10 +31,22 @@ namespace Differentiation * is to be used for computations. If a type that is selected for use * is not available in the library, a run-time error will be thrown. * - * @author Jean-Paul Pelteret, 2017 + * @author Jean-Paul Pelteret, 2017, 2019 */ enum class NumberTypes { + /** + * A dummy type for floating point numbers (i.e., non-differentiable + * scalar types). + * + * This option exists to facilitate the use of template meta-programming + * techniques to select, based on this enumeration, which + * auto-differentiable number type will be used to perform calculations. + * It will not permit any computations because the underlying number types + * resulting from its selection are floating point types, rather than + * auto-differentiable numbers. + */ + none, /** * Taped forward and reverse-mode ADOL-C number type (n-differentiable). * -- 2.39.5