Reference documentation for deal.II version Git a285ab78cc 2020-08-04 19:09:24 +0200
\(\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  {
36  template <typename NumberType, typename = void>
37  struct is_sacado_number : std::false_type
38  {};
39 
40 
46  template <typename NumberType, typename = void>
47  struct is_sacado_dfad_number : std::false_type
48  {};
49 
50 
56  template <typename NumberType, typename = void>
57  struct is_sacado_rad_number : std::false_type
58  {};
59 
60  } // namespace AD
61 } // namespace Differentiation
62 
63 
65 
66 
67 
68 #ifdef DEAL_II_TRILINOS_WITH_SACADO
69 
71 # include <Sacado.hpp>
72 // It appears that some versions of Trilinos do not directly or indirectly
73 // include all the headers for all forward and reverse Sacado AD types.
74 // So we directly include these both here as a precaution.
75 // Standard forward AD classes (templated)
76 # include <Sacado_Fad_DFad.hpp>
77 // Reverse AD classes (templated)
78 # include <Sacado_trad.hpp>
80 
81 # include <deal.II/base/exceptions.h>
82 # include <deal.II/base/numbers.h>
83 
86 
87 # include <complex>
88 
90 
91 
92 namespace Differentiation
93 {
94  namespace AD
95  {
96  namespace internal
97  {
104  template <typename SacadoNumber, typename = void>
106 
107  } // namespace internal
108 
109 
110 
111  } // namespace AD
112 } // namespace Differentiation
113 
114 
115 /* ----------- inline and template functions and specializations ----------- */
116 
117 
118 # ifndef DOXYGEN
119 
120 namespace Differentiation
121 {
122  namespace AD
123  {
124  namespace internal
125  {
126  // The documentation on Sacado numbers is pretty sparse and/or hard to
127  // navigate. As a point of reference, see
128  // https://trilinos.org/docs/dev/packages/sacado/doc/html/classSacado_1_1Fad_1_1SimpleFad.html
129  // for semi-applicable documentation for the Sacado::Fad::Dfad class.
130  // and the examples in
131  // https://github.com/trilinos/Trilinos/tree/master/packages/sacado/example
132  //
133  // If one dares to venture there, the relevant files for the classes
134  // supported here are:
135  //
136  // Forward-mode auto-differentiable types:
137  // https://github.com/trilinos/Trilinos/blob/master/packages/sacado/src/sacado_dfad_DFad.hpp
138  // https://github.com/trilinos/Trilinos/blob/master/packages/sacado/src/sacado_dfad_GeneralFad.hpp
139  //
140  // Reverse-mode auto-differentiable types:
141  // https://github.com/trilinos/Trilinos/blob/master/packages/sacado/src/Sacado_trad.hpp
142 
143 
147  template <typename SacadoNumber>
148  struct SacadoNumberInfo<
149  SacadoNumber,
150  typename std::enable_if<std::is_same<
151  SacadoNumber,
152  Sacado::Fad::DFad<typename SacadoNumber::value_type>>::value>::type>
153  {
154  using ad_type = SacadoNumber;
155  using scalar_type = typename ad_type::scalar_type;
156  using value_type = typename ad_type::value_type;
157  using derivative_type = typename ad_type::value_type;
158 
159  static const unsigned int n_supported_derivative_levels =
161  };
162 
163 
167  template <typename SacadoNumber>
168  struct SacadoNumberInfo<
169  SacadoNumber,
170  typename std::enable_if<std::is_same<
171  SacadoNumber,
172  Sacado::Rad::ADvar<typename SacadoNumber::value_type>>::value>::type>
173  {
174  using ad_type = SacadoNumber;
175  using scalar_type = typename ad_type::ADVari::scalar_type;
176  using value_type = typename ad_type::ADVari::value_type;
177  using derivative_type = typename ad_type::ADVari::value_type;
178 
179  static const unsigned int n_supported_derivative_levels =
181  };
182 
183 
190  template <typename Number>
191  struct SacadoNumberInfo<
192  Number,
193  typename std::enable_if<
194  std::is_arithmetic<typename std::decay<Number>::type>::value>::type>
195  {
196  static const unsigned int n_supported_derivative_levels = 0;
197  };
198 
199 
204  template <typename ScalarType>
205  struct ADNumberInfoFromEnum<
206  ScalarType,
208  typename std::enable_if<
209  std::is_floating_point<ScalarType>::value>::type>
210  {
211  static const bool is_taped = false;
212  using real_type = Sacado::Fad::DFad<ScalarType>;
213  using derivative_type =
215  static const unsigned int n_supported_derivative_levels =
217  };
218 
219 
224  template <typename ScalarType>
225  struct ADNumberInfoFromEnum<
226  ScalarType,
228  typename std::enable_if<
229  std::is_floating_point<ScalarType>::value>::type>
230  {
231  static const bool is_taped = false;
232  using real_type = Sacado::Fad::DFad<Sacado::Fad::DFad<ScalarType>>;
233  using derivative_type =
235  static const unsigned int n_supported_derivative_levels =
237  };
238 
239 
244  template <typename ScalarType>
245  struct ADNumberInfoFromEnum<
246  ScalarType,
248  typename std::enable_if<
249  std::is_floating_point<ScalarType>::value>::type>
250  {
251  static const bool is_taped = false;
252  using real_type = Sacado::Rad::ADvar<ScalarType>;
253  using derivative_type =
255  static const unsigned int n_supported_derivative_levels =
257  };
258 
259 
264  template <typename ScalarType>
265  struct ADNumberInfoFromEnum<
266  ScalarType,
268  typename std::enable_if<
269  std::is_floating_point<ScalarType>::value>::type>
270  {
271  static const bool is_taped = false;
272  using real_type = Sacado::Rad::ADvar<Sacado::Fad::DFad<ScalarType>>;
273  using derivative_type =
275  static const unsigned int n_supported_derivative_levels =
277  };
278 
279 
284  template <typename NumberType>
285  struct Marking<Sacado::Fad::DFad<NumberType>>
286  {
287  using ad_type =
289  using derivative_type = typename SacadoNumberInfo<
290  Sacado::Fad::DFad<NumberType>>::derivative_type;
291  using scalar_type =
292  typename SacadoNumberInfo<Sacado::Fad::DFad<NumberType>>::scalar_type;
293 
294  /*
295  * Initialize the state of an independent variable.
296  */
297  static void
298  independent_variable(const scalar_type &in,
299  const unsigned int index,
300  const unsigned int n_independent_variables,
301  ad_type & out)
302  {
303  // It is required that we first initialise the outer number before
304  // any of the nested ones.
305  out = ad_type(n_independent_variables, index, in);
306 
307  // Initialize potential nested directional derivatives
309  in, index, n_independent_variables, out.val());
310  }
311 
312  /*
313  * Initialize the state of a dependent variable.
314  */
315  static void
316  dependent_variable(ad_type &out, const ad_type &func)
317  {
318  out = func;
319  }
320  };
321 
322 
327  template <typename NumberType>
328  struct Marking<Sacado::Rad::ADvar<NumberType>>
329  {
330  using ad_type =
332  using derivative_type = typename SacadoNumberInfo<
333  Sacado::Rad::ADvar<NumberType>>::derivative_type;
334  using scalar_type = typename SacadoNumberInfo<
335  Sacado::Rad::ADvar<NumberType>>::scalar_type;
336 
337  /*
338  * Initialize the state of an independent variable.
339  */
340  static void
341  independent_variable(const scalar_type &in,
342  const unsigned int index,
343  const unsigned int n_independent_variables,
344  ad_type & out)
345  {
346  // For Sacado::Rad::ADvar numbers, we have to initialize the
347  // ADNumber with an already fully-configured value. This means
348  // that if this nests another ADNumber then the nested number
349  // must already be setup and ready for use.
350 
351  // Initialize potential nested directional derivatives
352  derivative_type derivative_initializer;
354  in, index, n_independent_variables, derivative_initializer);
355 
356  // Initialize the outer ad_type
357  out = derivative_initializer;
358  }
359 
360  /*
361  * Initialize the state of a dependent variable.
362  */
363  static void
364  dependent_variable(ad_type &out, const ad_type &func)
365  {
366  out = func;
367  }
368  };
369 
370 
376  template <typename NumberType>
377  struct ExtractData<Sacado::Fad::DFad<NumberType>>
378  {
379  using derivative_type = typename SacadoNumberInfo<
380  Sacado::Fad::DFad<NumberType>>::derivative_type;
381  using scalar_type =
382  typename SacadoNumberInfo<Sacado::Fad::DFad<NumberType>>::scalar_type;
383  using value_type =
385 
389  static scalar_type
390  value(const Sacado::Fad::DFad<NumberType> &x)
391  {
392  return ExtractData<value_type>::value(x.val());
393  }
394 
395 
399  static unsigned int
400  n_directional_derivatives(const Sacado::Fad::DFad<NumberType> &x)
401  {
402  return x.size();
403  }
404 
405 
409  static derivative_type
410  directional_derivative(const Sacado::Fad::DFad<NumberType> &x,
411  const unsigned int direction)
412  {
413  if (x.hasFastAccess())
414  return x.fastAccessDx(direction);
415  else
416  return x.dx(direction);
417  }
418  };
419 
420 
426  template <typename NumberType>
427  struct ExtractData<Sacado::Rad::ADvar<NumberType>>
428  {
429  using derivative_type = typename SacadoNumberInfo<
430  Sacado::Rad::ADvar<NumberType>>::derivative_type;
431  using scalar_type = typename SacadoNumberInfo<
432  Sacado::Rad::ADvar<NumberType>>::scalar_type;
433  using value_type =
435 
439  static scalar_type
440  value(const Sacado::Rad::ADvar<NumberType> &x)
441  {
442  return ExtractData<value_type>::value(x.val());
443  }
444 
445 
449  static unsigned int
450  n_directional_derivatives(const Sacado::Rad::ADvar<NumberType> &)
451  {
452  // There are as many directional derivatives as there are
453  // independent variables, but each independent variable can
454  // only return one directional derivative.
455  return 1;
456  }
457 
458 
466  static derivative_type
467  directional_derivative(const Sacado::Rad::ADvar<NumberType> &x,
468  const unsigned int)
469  {
470  return x.adj();
471  }
472  };
473 
474  } // namespace internal
475 
476 
477  /* -------------- NumberTypes::sacado_dfad -------------- */
478 
479 
486  template <typename ADNumberType>
487  struct ADNumberTraits<
488  ADNumberType,
489  typename std::enable_if<std::is_same<
490  ADNumberType,
491  Sacado::Fad::DFad<typename ADNumberType::scalar_type>>::value>::type>
492  : NumberTraits<typename ADNumberType::scalar_type,
493  NumberTypes::sacado_dfad>
494  {};
495 
496 
503  template <typename ADNumberType>
504  struct ADNumberTraits<
505  ADNumberType,
506  typename std::enable_if<std::is_same<
507  ADNumberType,
508  std::complex<Sacado::Fad::DFad<
509  typename ADNumberType::value_type::scalar_type>>>::value>::type>
510  : NumberTraits<
511  std::complex<typename ADNumberType::value_type::scalar_type>,
512  NumberTypes::sacado_dfad>
513  {};
514 
515 
520  template <>
521  struct NumberTraits<Sacado::Fad::DFad<float>, NumberTypes::sacado_dfad>
522  : NumberTraits<
523  typename ADNumberTraits<Sacado::Fad::DFad<float>>::scalar_type,
524  NumberTypes::sacado_dfad>
525  {};
526 
527 
532  template <>
533  struct NumberTraits<std::complex<Sacado::Fad::DFad<float>>,
535  : NumberTraits<typename ADNumberTraits<
536  std::complex<Sacado::Fad::DFad<float>>>::scalar_type,
537  NumberTypes::sacado_dfad>
538  {};
539 
540 
545  template <>
546  struct NumberTraits<Sacado::Fad::DFad<double>, NumberTypes::sacado_dfad>
547  : NumberTraits<
548  typename ADNumberTraits<Sacado::Fad::DFad<double>>::scalar_type,
549  NumberTypes::sacado_dfad>
550  {};
551 
552 
557  template <>
558  struct NumberTraits<std::complex<Sacado::Fad::DFad<double>>,
560  : NumberTraits<typename ADNumberTraits<
561  std::complex<Sacado::Fad::DFad<double>>>::scalar_type,
562  NumberTypes::sacado_dfad>
563  {};
564 
565 
566  /* -------------- NumberTypes::sacado_rad -------------- */
567 
568 
575  template <typename ADNumberType>
576  struct ADNumberTraits<
577  ADNumberType,
578  typename std::enable_if<std::is_same<
579  ADNumberType,
580  Sacado::Rad::ADvar<typename ADNumberType::ADVari::scalar_type>>::
581  value>::type>
582  : NumberTraits<typename ADNumberType::ADVari::scalar_type,
583  NumberTypes::sacado_rad>
584  {};
585 
586 
591  template <>
592  struct NumberTraits<Sacado::Rad::ADvar<float>, NumberTypes::sacado_rad>
593  : NumberTraits<
594  typename ADNumberTraits<Sacado::Rad::ADvar<float>>::scalar_type,
595  NumberTypes::sacado_rad>
596  {};
597 
598 
603  template <>
604  struct NumberTraits<Sacado::Rad::ADvar<double>, NumberTypes::sacado_rad>
605  : NumberTraits<
606  typename ADNumberTraits<Sacado::Rad::ADvar<double>>::scalar_type,
607  NumberTypes::sacado_rad>
608  {};
609 
610 
611 # ifdef DEAL_II_TRILINOS_CXX_SUPPORTS_SACADO_COMPLEX_RAD
612 
613 
620  template <typename ADNumberType>
621  struct ADNumberTraits<
622  ADNumberType,
623  typename std::enable_if<std::is_same<
624  ADNumberType,
625  std::complex<Sacado::Rad::ADvar<typename ADNumberType::value_type::
626  ADVari::scalar_type>>>::value>::type>
627  : NumberTraits<
628  std::complex<typename ADNumberType::value_type::ADVari::scalar_type>,
629  NumberTypes::sacado_rad>
630  {};
631 
632 
637  template <>
638  struct NumberTraits<std::complex<Sacado::Rad::ADvar<float>>,
640  : NumberTraits<typename ADNumberTraits<
641  std::complex<Sacado::Rad::ADvar<float>>>::scalar_type,
642  NumberTypes::sacado_rad>
643  {};
644 
645 
650  template <>
651  struct NumberTraits<std::complex<Sacado::Rad::ADvar<double>>,
653  : NumberTraits<typename ADNumberTraits<
654  std::complex<Sacado::Rad::ADvar<double>>>::scalar_type,
655  NumberTypes::sacado_rad>
656  {};
657 
658 
659 # endif
660 
661 
662  /* -------------- NumberTypes::sacado_dfad_dfad -------------- */
663 
671  template <typename ADNumberType>
672  struct ADNumberTraits<
673  ADNumberType,
674  typename std::enable_if<
675  std::is_same<ADNumberType,
676  Sacado::Fad::DFad<Sacado::Fad::DFad<
677  typename ADNumberType::scalar_type>>>::value>::type>
678  : NumberTraits<typename ADNumberType::scalar_type,
679  NumberTypes::sacado_dfad_dfad>
680  {};
681 
682 
690  template <typename ADNumberType>
691  struct ADNumberTraits<
692  ADNumberType,
693  typename std::enable_if<std::is_same<
694  ADNumberType,
695  std::complex<Sacado::Fad::DFad<Sacado::Fad::DFad<
696  typename ADNumberType::value_type::scalar_type>>>>::value>::type>
697  : NumberTraits<
698  std::complex<typename ADNumberType::value_type::scalar_type>,
699  NumberTypes::sacado_dfad_dfad>
700  {};
701 
702 
707  template <>
708  struct NumberTraits<Sacado::Fad::DFad<Sacado::Fad::DFad<float>>,
710  : NumberTraits<typename ADNumberTraits<Sacado::Fad::DFad<
711  Sacado::Fad::DFad<float>>>::scalar_type,
712  NumberTypes::sacado_dfad_dfad>
713  {};
714 
715 
720  template <>
721  struct NumberTraits<
722  std::complex<Sacado::Fad::DFad<Sacado::Fad::DFad<float>>>,
724  : NumberTraits<typename ADNumberTraits<std::complex<Sacado::Fad::DFad<
725  Sacado::Fad::DFad<float>>>>::scalar_type,
726  NumberTypes::sacado_dfad_dfad>
727  {};
728 
729 
734  template <>
735  struct NumberTraits<Sacado::Fad::DFad<Sacado::Fad::DFad<double>>,
737  : NumberTraits<typename ADNumberTraits<Sacado::Fad::DFad<
738  Sacado::Fad::DFad<double>>>::scalar_type,
739  NumberTypes::sacado_dfad_dfad>
740  {};
741 
742 
747  template <>
748  struct NumberTraits<
749  std::complex<Sacado::Fad::DFad<Sacado::Fad::DFad<double>>>,
751  : NumberTraits<typename ADNumberTraits<std::complex<Sacado::Fad::DFad<
752  Sacado::Fad::DFad<double>>>>::scalar_type,
753  NumberTypes::sacado_dfad_dfad>
754  {};
755 
756 
757  /* -------------- NumberTypes::sacado_rad_dfad -------------- */
758 
766  template <typename ADNumberType>
767  struct ADNumberTraits<
768  ADNumberType,
769  typename std::enable_if<std::is_same<
770  ADNumberType,
771  Sacado::Rad::ADvar<Sacado::Fad::DFad<
772  typename ADNumberType::ADVari::scalar_type>>>::value>::type>
773  : NumberTraits<typename ADNumberType::ADVari::scalar_type,
774  NumberTypes::sacado_rad_dfad>
775  {};
776 
777 
782  template <>
783  struct NumberTraits<Sacado::Rad::ADvar<Sacado::Fad::DFad<float>>,
785  : NumberTraits<typename ADNumberTraits<Sacado::Rad::ADvar<
786  Sacado::Fad::DFad<float>>>::scalar_type,
787  NumberTypes::sacado_rad_dfad>
788  {};
789 
790 
795  template <>
796  struct NumberTraits<Sacado::Rad::ADvar<Sacado::Fad::DFad<double>>,
798  : NumberTraits<typename ADNumberTraits<Sacado::Rad::ADvar<
799  Sacado::Fad::DFad<double>>>::scalar_type,
800  NumberTypes::sacado_rad_dfad>
801  {};
802 
803 
804 # ifdef DEAL_II_TRILINOS_CXX_SUPPORTS_SACADO_COMPLEX_RAD
805 
806 
814  template <typename ADNumberType>
815  struct ADNumberTraits<
816  ADNumberType,
817  typename std::enable_if<std::is_same<
818  ADNumberType,
819  std::complex<Sacado::Rad::ADvar<Sacado::Fad::DFad<
820  typename ADNumberType::value_type::ADVari::scalar_type>>>>::value>::
821  type>
822  : NumberTraits<
823  std::complex<typename ADNumberType::value_type::ADVari::scalar_type>,
824  NumberTypes::sacado_rad_dfad>
825  {};
826 
827 
832  template <>
833  struct NumberTraits<
834  std::complex<Sacado::Rad::ADvar<Sacado::Fad::DFad<float>>>,
836  : NumberTraits<typename ADNumberTraits<std::complex<Sacado::Rad::ADvar<
837  Sacado::Fad::DFad<float>>>>::scalar_type,
838  NumberTypes::sacado_rad_dfad>
839  {};
840 
841 
846  template <>
847  struct NumberTraits<
848  std::complex<Sacado::Rad::ADvar<Sacado::Fad::DFad<double>>>,
850  : NumberTraits<typename ADNumberTraits<std::complex<Sacado::Rad::ADvar<
851  Sacado::Fad::DFad<double>>>>::scalar_type,
852  NumberTypes::sacado_rad_dfad>
853  {};
854 
855 
856 # endif
857 
858 
859  /* -------------- Additional type traits -------------- */
860 
861 
862  template <typename NumberType>
863  struct is_sacado_dfad_number<
864  NumberType,
865  typename std::enable_if<
866  ADNumberTraits<typename std::decay<NumberType>::type>::type_code ==
867  NumberTypes::sacado_dfad ||
868  ADNumberTraits<typename std::decay<NumberType>::type>::type_code ==
869  NumberTypes::sacado_dfad_dfad>::type> : std::true_type
870  {};
871 
872 
873  template <typename NumberType>
874  struct is_sacado_dfad_number<
875  NumberType,
876  typename std::enable_if<std::is_same<
877  NumberType,
878  Sacado::Fad::Expr<typename NumberType::value_type>>::value>::type>
879  : std::true_type
880  {};
881 
882 
883  template <typename NumberType>
884  struct is_sacado_rad_number<
885  NumberType,
886  typename std::enable_if<
887  ADNumberTraits<typename std::decay<NumberType>::type>::type_code ==
888  NumberTypes::sacado_rad ||
889  ADNumberTraits<typename std::decay<NumberType>::type>::type_code ==
890  NumberTypes::sacado_rad_dfad>::type> : std::true_type
891  {};
892 
893 
894  template <typename NumberType>
895  struct is_sacado_rad_number<
896  NumberType,
897  typename std::enable_if<std::is_same<
898  NumberType,
899  Sacado::Rad::ADvari<Sacado::Fad::DFad<
900  typename NumberType::ADVari::scalar_type>>>::value>::type>
901  : std::true_type
902  {};
903 
904 
905  template <typename NumberType>
906  struct is_sacado_number<
907  NumberType,
908  typename std::enable_if<is_sacado_dfad_number<NumberType>::value ||
909  is_sacado_rad_number<NumberType>::value>::type>
910  : std::true_type
911  {};
912 
913  } // namespace AD
914 } // namespace Differentiation
915 
916 
917 namespace numbers
918 {
919  template <typename NumberType>
920  struct is_cuda_compatible<
921  NumberType,
922  typename std::enable_if<::Differentiation::AD::is_sacado_rad_number<
923  NumberType>::value>::type> : std::false_type
924  {};
925 } // namespace numbers
926 
927 # endif // DOXYGEN
928 
929 
930 
932 
933 
934 #endif // DEAL_II_TRILINOS_WITH_SACADO
935 
936 #endif
STL namespace.
#define DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
Definition: config.h:376
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
#define DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
Definition: config.h:413
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362