Reference documentation for deal.II version 9.2.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
sacado_number_types.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2016 - 2020 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_differentiation_ad_sacado_number_types_h
17 #define dealii_differentiation_ad_sacado_number_types_h
18 
19 #include <deal.II/base/config.h>
20 
21 #include <type_traits>
22 
23 
25 
26 
27 namespace Differentiation
28 {
29  namespace AD
30  {
38  template <typename NumberType, typename = void>
39  struct is_sacado_number : std::false_type
40  {};
41 
42 
50  template <typename NumberType, typename = void>
51  struct is_sacado_dfad_number : std::false_type
52  {};
53 
54 
62  template <typename NumberType, typename = void>
63  struct is_sacado_rad_number : std::false_type
64  {};
65 
66  } // namespace AD
67 } // namespace Differentiation
68 
69 
71 
72 
73 
74 #ifdef DEAL_II_TRILINOS_WITH_SACADO
75 
77 # include <Sacado.hpp>
78 // It appears that some versions of Trilinos do not directly or indirectly
79 // include all the headers for all forward and reverse Sacado AD types.
80 // So we directly include these both here as a precaution.
81 // Standard forward AD classes (templated)
82 # include <Sacado_Fad_DFad.hpp>
83 // Reverse AD classes (templated)
84 # include <Sacado_trad.hpp>
86 
87 # include <deal.II/base/exceptions.h>
88 # include <deal.II/base/numbers.h>
89 
92 
93 # include <complex>
94 
96 
97 
98 namespace Differentiation
99 {
100  namespace AD
101  {
102  namespace internal
103  {
112  template <typename SacadoNumber, typename = void>
114 
115  } // namespace internal
116 
117 
118 
119  } // namespace AD
120 } // namespace Differentiation
121 
122 
123 /* ----------- inline and template functions and specializations ----------- */
124 
125 
126 # ifndef DOXYGEN
127 
128 namespace Differentiation
129 {
130  namespace AD
131  {
132  namespace internal
133  {
134  // The documentation on Sacado numbers is pretty sparse and/or hard to
135  // navigate. As a point of reference, see
136  // https://trilinos.org/docs/dev/packages/sacado/doc/html/classSacado_1_1Fad_1_1SimpleFad.html
137  // for semi-applicable documentation for the Sacado::Fad::Dfad class.
138  // and the examples in
139  // https://github.com/trilinos/Trilinos/tree/master/packages/sacado/example
140  //
141  // If one dares to venture there, the relevant files for the classes
142  // supported here are:
143  //
144  // Forward-mode auto-differentiable types:
145  // https://github.com/trilinos/Trilinos/blob/master/packages/sacado/src/sacado_dfad_DFad.hpp
146  // https://github.com/trilinos/Trilinos/blob/master/packages/sacado/src/sacado_dfad_GeneralFad.hpp
147  //
148  // Reverse-mode auto-differentiable types:
149  // https://github.com/trilinos/Trilinos/blob/master/packages/sacado/src/Sacado_trad.hpp
150 
151 
155  template <typename SacadoNumber>
156  struct SacadoNumberInfo<
157  SacadoNumber,
158  typename std::enable_if<std::is_same<
159  SacadoNumber,
160  Sacado::Fad::DFad<typename SacadoNumber::value_type>>::value>::type>
161  {
162  using ad_type = SacadoNumber;
163  using scalar_type = typename ad_type::scalar_type;
164  using value_type = typename ad_type::value_type;
165  using derivative_type = typename ad_type::value_type;
166 
167  static const unsigned int n_supported_derivative_levels =
169  };
170 
171 
175  template <typename SacadoNumber>
176  struct SacadoNumberInfo<
177  SacadoNumber,
178  typename std::enable_if<std::is_same<
179  SacadoNumber,
180  Sacado::Rad::ADvar<typename SacadoNumber::value_type>>::value>::type>
181  {
182  using ad_type = SacadoNumber;
183  using scalar_type = typename ad_type::ADVari::scalar_type;
184  using value_type = typename ad_type::ADVari::value_type;
185  using derivative_type = typename ad_type::ADVari::value_type;
186 
187  static const unsigned int n_supported_derivative_levels =
189  };
190 
191 
198  template <typename Number>
199  struct SacadoNumberInfo<
200  Number,
201  typename std::enable_if<
202  std::is_arithmetic<typename std::decay<Number>::type>::value>::type>
203  {
204  static const unsigned int n_supported_derivative_levels = 0;
205  };
206 
207 
212  template <typename ScalarType>
213  struct ADNumberInfoFromEnum<
214  ScalarType,
216  typename std::enable_if<
217  std::is_floating_point<ScalarType>::value>::type>
218  {
219  static const bool is_taped = false;
220  using real_type = Sacado::Fad::DFad<ScalarType>;
221  using derivative_type =
223  static const unsigned int n_supported_derivative_levels =
225  };
226 
227 
232  template <typename ScalarType>
233  struct ADNumberInfoFromEnum<
234  ScalarType,
236  typename std::enable_if<
237  std::is_floating_point<ScalarType>::value>::type>
238  {
239  static const bool is_taped = false;
240  using real_type = Sacado::Fad::DFad<Sacado::Fad::DFad<ScalarType>>;
241  using derivative_type =
243  static const unsigned int n_supported_derivative_levels =
245  };
246 
247 
252  template <typename ScalarType>
253  struct ADNumberInfoFromEnum<
254  ScalarType,
256  typename std::enable_if<
257  std::is_floating_point<ScalarType>::value>::type>
258  {
259  static const bool is_taped = false;
260  using real_type = Sacado::Rad::ADvar<ScalarType>;
261  using derivative_type =
263  static const unsigned int n_supported_derivative_levels =
265  };
266 
267 
272  template <typename ScalarType>
273  struct ADNumberInfoFromEnum<
274  ScalarType,
276  typename std::enable_if<
277  std::is_floating_point<ScalarType>::value>::type>
278  {
279  static const bool is_taped = false;
280  using real_type = Sacado::Rad::ADvar<Sacado::Fad::DFad<ScalarType>>;
281  using derivative_type =
283  static const unsigned int n_supported_derivative_levels =
285  };
286 
287 
292  template <typename NumberType>
293  struct Marking<Sacado::Fad::DFad<NumberType>>
294  {
295  using ad_type =
297  using derivative_type = typename SacadoNumberInfo<
298  Sacado::Fad::DFad<NumberType>>::derivative_type;
299  using scalar_type =
300  typename SacadoNumberInfo<Sacado::Fad::DFad<NumberType>>::scalar_type;
301 
302  /*
303  * Initialize the state of an independent variable.
304  */
305  static void
306  independent_variable(const scalar_type &in,
307  const unsigned int index,
308  const unsigned int n_independent_variables,
309  ad_type & out)
310  {
311  // It is required that we first initialise the outer number before
312  // any of the nested ones.
313  out = ad_type(n_independent_variables, index, in);
314 
315  // Initialize potential nested directional derivatives
317  in, index, n_independent_variables, out.val());
318  }
319 
320  /*
321  * Initialize the state of a dependent variable.
322  */
323  static void
324  dependent_variable(ad_type &out, const ad_type &func)
325  {
326  out = func;
327  }
328  };
329 
330 
335  template <typename NumberType>
336  struct Marking<Sacado::Rad::ADvar<NumberType>>
337  {
338  using ad_type =
340  using derivative_type = typename SacadoNumberInfo<
341  Sacado::Rad::ADvar<NumberType>>::derivative_type;
342  using scalar_type = typename SacadoNumberInfo<
343  Sacado::Rad::ADvar<NumberType>>::scalar_type;
344 
345  /*
346  * Initialize the state of an independent variable.
347  */
348  static void
349  independent_variable(const scalar_type &in,
350  const unsigned int index,
351  const unsigned int n_independent_variables,
352  ad_type & out)
353  {
354  // For Sacado::Rad::ADvar numbers, we have to initialize the
355  // ADNumber with an already fully-configured value. This means
356  // that if this nests another ADNumber then the nested number
357  // must already be setup and ready for use.
358 
359  // Initialize potential nested directional derivatives
360  derivative_type derivative_initializer;
362  in, index, n_independent_variables, derivative_initializer);
363 
364  // Initialize the outer ad_type
365  out = derivative_initializer;
366  }
367 
368  /*
369  * Initialize the state of a dependent variable.
370  */
371  static void
372  dependent_variable(ad_type &out, const ad_type &func)
373  {
374  out = func;
375  }
376  };
377 
378 
386  template <typename NumberType>
387  struct ExtractData<Sacado::Fad::DFad<NumberType>>
388  {
389  using derivative_type = typename SacadoNumberInfo<
390  Sacado::Fad::DFad<NumberType>>::derivative_type;
391  using scalar_type =
392  typename SacadoNumberInfo<Sacado::Fad::DFad<NumberType>>::scalar_type;
393  using value_type =
395 
399  static scalar_type
400  value(const Sacado::Fad::DFad<NumberType> &x)
401  {
402  return ExtractData<value_type>::value(x.val());
403  }
404 
405 
409  static unsigned int
410  n_directional_derivatives(const Sacado::Fad::DFad<NumberType> &x)
411  {
412  return x.size();
413  }
414 
415 
419  static derivative_type
420  directional_derivative(const Sacado::Fad::DFad<NumberType> &x,
421  const unsigned int direction)
422  {
423  if (x.hasFastAccess())
424  return x.fastAccessDx(direction);
425  else
426  return x.dx(direction);
427  }
428  };
429 
430 
438  template <typename NumberType>
439  struct ExtractData<Sacado::Rad::ADvar<NumberType>>
440  {
441  using derivative_type = typename SacadoNumberInfo<
442  Sacado::Rad::ADvar<NumberType>>::derivative_type;
443  using scalar_type = typename SacadoNumberInfo<
444  Sacado::Rad::ADvar<NumberType>>::scalar_type;
445  using value_type =
447 
451  static scalar_type
452  value(const Sacado::Rad::ADvar<NumberType> &x)
453  {
454  return ExtractData<value_type>::value(x.val());
455  }
456 
457 
461  static unsigned int
462  n_directional_derivatives(const Sacado::Rad::ADvar<NumberType> &)
463  {
464  // There are as many directional derivatives as there are
465  // independent variables, but each independent variable can
466  // only return one directional derivative.
467  return 1;
468  }
469 
470 
478  static derivative_type
479  directional_derivative(const Sacado::Rad::ADvar<NumberType> &x,
480  const unsigned int)
481  {
482  return x.adj();
483  }
484  };
485 
486  } // namespace internal
487 
488 
489  /* -------------- NumberTypes::sacado_dfad -------------- */
490 
491 
498  template <typename ADNumberType>
499  struct ADNumberTraits<
500  ADNumberType,
501  typename std::enable_if<std::is_same<
502  ADNumberType,
503  Sacado::Fad::DFad<typename ADNumberType::scalar_type>>::value>::type>
504  : NumberTraits<typename ADNumberType::scalar_type,
505  NumberTypes::sacado_dfad>
506  {};
507 
508 
515  template <typename ADNumberType>
516  struct ADNumberTraits<
517  ADNumberType,
518  typename std::enable_if<std::is_same<
519  ADNumberType,
520  std::complex<Sacado::Fad::DFad<
521  typename ADNumberType::value_type::scalar_type>>>::value>::type>
522  : NumberTraits<
523  std::complex<typename ADNumberType::value_type::scalar_type>,
524  NumberTypes::sacado_dfad>
525  {};
526 
527 
532  template <>
533  struct NumberTraits<Sacado::Fad::DFad<float>, NumberTypes::sacado_dfad>
534  : NumberTraits<
535  typename ADNumberTraits<Sacado::Fad::DFad<float>>::scalar_type,
536  NumberTypes::sacado_dfad>
537  {};
538 
539 
544  template <>
545  struct NumberTraits<std::complex<Sacado::Fad::DFad<float>>,
547  : NumberTraits<typename ADNumberTraits<
548  std::complex<Sacado::Fad::DFad<float>>>::scalar_type,
549  NumberTypes::sacado_dfad>
550  {};
551 
552 
557  template <>
558  struct NumberTraits<Sacado::Fad::DFad<double>, NumberTypes::sacado_dfad>
559  : NumberTraits<
560  typename ADNumberTraits<Sacado::Fad::DFad<double>>::scalar_type,
561  NumberTypes::sacado_dfad>
562  {};
563 
564 
569  template <>
570  struct NumberTraits<std::complex<Sacado::Fad::DFad<double>>,
572  : NumberTraits<typename ADNumberTraits<
573  std::complex<Sacado::Fad::DFad<double>>>::scalar_type,
574  NumberTypes::sacado_dfad>
575  {};
576 
577 
578  /* -------------- NumberTypes::sacado_rad -------------- */
579 
580 
587  template <typename ADNumberType>
588  struct ADNumberTraits<
589  ADNumberType,
590  typename std::enable_if<std::is_same<
591  ADNumberType,
592  Sacado::Rad::ADvar<typename ADNumberType::ADVari::scalar_type>>::
593  value>::type>
594  : NumberTraits<typename ADNumberType::ADVari::scalar_type,
595  NumberTypes::sacado_rad>
596  {};
597 
598 
603  template <>
604  struct NumberTraits<Sacado::Rad::ADvar<float>, NumberTypes::sacado_rad>
605  : NumberTraits<
606  typename ADNumberTraits<Sacado::Rad::ADvar<float>>::scalar_type,
607  NumberTypes::sacado_rad>
608  {};
609 
610 
615  template <>
616  struct NumberTraits<Sacado::Rad::ADvar<double>, NumberTypes::sacado_rad>
617  : NumberTraits<
618  typename ADNumberTraits<Sacado::Rad::ADvar<double>>::scalar_type,
619  NumberTypes::sacado_rad>
620  {};
621 
622 
623 # ifdef DEAL_II_TRILINOS_CXX_SUPPORTS_SACADO_COMPLEX_RAD
624 
625 
632  template <typename ADNumberType>
633  struct ADNumberTraits<
634  ADNumberType,
635  typename std::enable_if<std::is_same<
636  ADNumberType,
637  std::complex<Sacado::Rad::ADvar<typename ADNumberType::value_type::
638  ADVari::scalar_type>>>::value>::type>
639  : NumberTraits<
640  std::complex<typename ADNumberType::value_type::ADVari::scalar_type>,
641  NumberTypes::sacado_rad>
642  {};
643 
644 
649  template <>
650  struct NumberTraits<std::complex<Sacado::Rad::ADvar<float>>,
652  : NumberTraits<typename ADNumberTraits<
653  std::complex<Sacado::Rad::ADvar<float>>>::scalar_type,
654  NumberTypes::sacado_rad>
655  {};
656 
657 
662  template <>
663  struct NumberTraits<std::complex<Sacado::Rad::ADvar<double>>,
665  : NumberTraits<typename ADNumberTraits<
666  std::complex<Sacado::Rad::ADvar<double>>>::scalar_type,
667  NumberTypes::sacado_rad>
668  {};
669 
670 
671 # endif
672 
673 
674  /* -------------- NumberTypes::sacado_dfad_dfad -------------- */
675 
683  template <typename ADNumberType>
684  struct ADNumberTraits<
685  ADNumberType,
686  typename std::enable_if<
687  std::is_same<ADNumberType,
688  Sacado::Fad::DFad<Sacado::Fad::DFad<
689  typename ADNumberType::scalar_type>>>::value>::type>
690  : NumberTraits<typename ADNumberType::scalar_type,
691  NumberTypes::sacado_dfad_dfad>
692  {};
693 
694 
702  template <typename ADNumberType>
703  struct ADNumberTraits<
704  ADNumberType,
705  typename std::enable_if<std::is_same<
706  ADNumberType,
707  std::complex<Sacado::Fad::DFad<Sacado::Fad::DFad<
708  typename ADNumberType::value_type::scalar_type>>>>::value>::type>
709  : NumberTraits<
710  std::complex<typename ADNumberType::value_type::scalar_type>,
711  NumberTypes::sacado_dfad_dfad>
712  {};
713 
714 
719  template <>
720  struct NumberTraits<Sacado::Fad::DFad<Sacado::Fad::DFad<float>>,
722  : NumberTraits<typename ADNumberTraits<Sacado::Fad::DFad<
723  Sacado::Fad::DFad<float>>>::scalar_type,
724  NumberTypes::sacado_dfad_dfad>
725  {};
726 
727 
732  template <>
733  struct NumberTraits<
734  std::complex<Sacado::Fad::DFad<Sacado::Fad::DFad<float>>>,
736  : NumberTraits<typename ADNumberTraits<std::complex<Sacado::Fad::DFad<
737  Sacado::Fad::DFad<float>>>>::scalar_type,
738  NumberTypes::sacado_dfad_dfad>
739  {};
740 
741 
746  template <>
747  struct NumberTraits<Sacado::Fad::DFad<Sacado::Fad::DFad<double>>,
749  : NumberTraits<typename ADNumberTraits<Sacado::Fad::DFad<
750  Sacado::Fad::DFad<double>>>::scalar_type,
751  NumberTypes::sacado_dfad_dfad>
752  {};
753 
754 
759  template <>
760  struct NumberTraits<
761  std::complex<Sacado::Fad::DFad<Sacado::Fad::DFad<double>>>,
763  : NumberTraits<typename ADNumberTraits<std::complex<Sacado::Fad::DFad<
764  Sacado::Fad::DFad<double>>>>::scalar_type,
765  NumberTypes::sacado_dfad_dfad>
766  {};
767 
768 
769  /* -------------- NumberTypes::sacado_rad_dfad -------------- */
770 
778  template <typename ADNumberType>
779  struct ADNumberTraits<
780  ADNumberType,
781  typename std::enable_if<std::is_same<
782  ADNumberType,
783  Sacado::Rad::ADvar<Sacado::Fad::DFad<
784  typename ADNumberType::ADVari::scalar_type>>>::value>::type>
785  : NumberTraits<typename ADNumberType::ADVari::scalar_type,
786  NumberTypes::sacado_rad_dfad>
787  {};
788 
789 
794  template <>
795  struct NumberTraits<Sacado::Rad::ADvar<Sacado::Fad::DFad<float>>,
797  : NumberTraits<typename ADNumberTraits<Sacado::Rad::ADvar<
798  Sacado::Fad::DFad<float>>>::scalar_type,
799  NumberTypes::sacado_rad_dfad>
800  {};
801 
802 
807  template <>
808  struct NumberTraits<Sacado::Rad::ADvar<Sacado::Fad::DFad<double>>,
810  : NumberTraits<typename ADNumberTraits<Sacado::Rad::ADvar<
811  Sacado::Fad::DFad<double>>>::scalar_type,
812  NumberTypes::sacado_rad_dfad>
813  {};
814 
815 
816 # ifdef DEAL_II_TRILINOS_CXX_SUPPORTS_SACADO_COMPLEX_RAD
817 
818 
826  template <typename ADNumberType>
827  struct ADNumberTraits<
828  ADNumberType,
829  typename std::enable_if<std::is_same<
830  ADNumberType,
831  std::complex<Sacado::Rad::ADvar<Sacado::Fad::DFad<
832  typename ADNumberType::value_type::ADVari::scalar_type>>>>::value>::
833  type>
834  : NumberTraits<
835  std::complex<typename ADNumberType::value_type::ADVari::scalar_type>,
836  NumberTypes::sacado_rad_dfad>
837  {};
838 
839 
844  template <>
845  struct NumberTraits<
846  std::complex<Sacado::Rad::ADvar<Sacado::Fad::DFad<float>>>,
848  : NumberTraits<typename ADNumberTraits<std::complex<Sacado::Rad::ADvar<
849  Sacado::Fad::DFad<float>>>>::scalar_type,
850  NumberTypes::sacado_rad_dfad>
851  {};
852 
853 
858  template <>
859  struct NumberTraits<
860  std::complex<Sacado::Rad::ADvar<Sacado::Fad::DFad<double>>>,
862  : NumberTraits<typename ADNumberTraits<std::complex<Sacado::Rad::ADvar<
863  Sacado::Fad::DFad<double>>>>::scalar_type,
864  NumberTypes::sacado_rad_dfad>
865  {};
866 
867 
868 # endif
869 
870 
871  /* -------------- Additional type traits -------------- */
872 
873 
874  template <typename NumberType>
875  struct is_sacado_dfad_number<
876  NumberType,
877  typename std::enable_if<
878  ADNumberTraits<typename std::decay<NumberType>::type>::type_code ==
879  NumberTypes::sacado_dfad ||
880  ADNumberTraits<typename std::decay<NumberType>::type>::type_code ==
881  NumberTypes::sacado_dfad_dfad>::type> : std::true_type
882  {};
883 
884 
885  template <typename NumberType>
886  struct is_sacado_dfad_number<
887  NumberType,
888  typename std::enable_if<std::is_same<
889  NumberType,
890  Sacado::Fad::Expr<typename NumberType::value_type>>::value>::type>
891  : std::true_type
892  {};
893 
894 
895  template <typename NumberType>
896  struct is_sacado_rad_number<
897  NumberType,
898  typename std::enable_if<
899  ADNumberTraits<typename std::decay<NumberType>::type>::type_code ==
900  NumberTypes::sacado_rad ||
901  ADNumberTraits<typename std::decay<NumberType>::type>::type_code ==
902  NumberTypes::sacado_rad_dfad>::type> : std::true_type
903  {};
904 
905 
906  template <typename NumberType>
907  struct is_sacado_rad_number<
908  NumberType,
909  typename std::enable_if<std::is_same<
910  NumberType,
911  Sacado::Rad::ADvari<Sacado::Fad::DFad<
912  typename NumberType::ADVari::scalar_type>>>::value>::type>
913  : std::true_type
914  {};
915 
916 
917  template <typename NumberType>
918  struct is_sacado_number<
919  NumberType,
920  typename std::enable_if<is_sacado_dfad_number<NumberType>::value ||
921  is_sacado_rad_number<NumberType>::value>::type>
922  : std::true_type
923  {};
924 
925  } // namespace AD
926 } // namespace Differentiation
927 
928 
929 namespace numbers
930 {
931  template <typename NumberType>
932  struct is_cuda_compatible<
933  NumberType,
934  typename std::enable_if<::Differentiation::AD::is_sacado_rad_number<
935  NumberType>::value>::type> : std::false_type
936  {};
937 } // namespace numbers
938 
939 # endif // DOXYGEN
940 
941 
942 
944 
945 
946 #endif // DEAL_II_TRILINOS_WITH_SACADO
947 
948 #endif
STL namespace.
#define DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
Definition: config.h:372
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
#define DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
Definition: config.h:409
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
static const bool value