Reference documentation for deal.II version 9.5.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\}}\)
Loading...
Searching...
No Matches
sacado_number_types.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2016 - 2023 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
27namespace 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 <deal.II/base/numbers.h>
72
75
76// It appears that some versions of Trilinos do not directly or indirectly
77// include all the headers for all forward and reverse Sacado AD types
78// in Sacado.hpp, so we also directly include these here as a precaution:
79# include <Sacado.hpp>
80# include <Sacado_Fad_DFad.hpp>
81# include <Sacado_trad.hpp>
82
83# include <complex>
84
86
87
88namespace Differentiation
89{
90 namespace AD
91 {
92 namespace internal
93 {
100 template <typename SacadoNumber, typename = void>
102
103 } // namespace internal
104
105
106
107 } // namespace AD
108} // namespace Differentiation
109
110
111/* ----------- inline and template functions and specializations ----------- */
112
113
114# ifndef DOXYGEN
115
116namespace Differentiation
117{
118 namespace AD
119 {
120 namespace internal
121 {
122 // The documentation on Sacado numbers is pretty sparse and/or hard to
123 // navigate. As a point of reference, see
124 // https://trilinos.org/docs/dev/packages/sacado/doc/html/classSacado_1_1Fad_1_1SimpleFad.html
125 // for semi-applicable documentation for the Sacado::Fad::Dfad class.
126 // and the examples in
127 // https://github.com/trilinos/Trilinos/tree/master/packages/sacado/example
128 //
129 // If one dares to venture there, the relevant files for the classes
130 // supported here are:
131 //
132 // Forward-mode auto-differentiable types:
133 // https://github.com/trilinos/Trilinos/blob/master/packages/sacado/src/sacado_dfad_DFad.hpp
134 // https://github.com/trilinos/Trilinos/blob/master/packages/sacado/src/sacado_dfad_GeneralFad.hpp
135 //
136 // Reverse-mode auto-differentiable types:
137 // https://github.com/trilinos/Trilinos/blob/master/packages/sacado/src/Sacado_trad.hpp
138
139
143 template <typename SacadoNumber>
144 struct SacadoNumberInfo<
145 SacadoNumber,
146 std::enable_if_t<std::is_same<
147 SacadoNumber,
148 Sacado::Fad::DFad<typename SacadoNumber::value_type>>::value>>
149 {
150 using ad_type = SacadoNumber;
151 using scalar_type = typename ad_type::scalar_type;
152 using value_type = typename ad_type::value_type;
153 using derivative_type = typename ad_type::value_type;
154
155 static const unsigned int n_supported_derivative_levels =
156 1 + SacadoNumberInfo<derivative_type>::n_supported_derivative_levels;
157 };
158
159
163 template <typename SacadoNumber>
164 struct SacadoNumberInfo<
165 SacadoNumber,
166 std::enable_if_t<std::is_same<
167 SacadoNumber,
168 Sacado::Rad::ADvar<typename SacadoNumber::value_type>>::value>>
169 {
170 using ad_type = SacadoNumber;
171 using scalar_type = typename ad_type::ADVari::scalar_type;
172 using value_type = typename ad_type::ADVari::value_type;
173 using derivative_type = typename ad_type::ADVari::value_type;
174
175 static const unsigned int n_supported_derivative_levels =
176 1 + SacadoNumberInfo<derivative_type>::n_supported_derivative_levels;
177 };
178
179
186 template <typename Number>
187 struct SacadoNumberInfo<Number,
188 std::enable_if_t<std::is_arithmetic<
189 typename std::decay<Number>::type>::value>>
190 {
191 static const unsigned int n_supported_derivative_levels = 0;
192 };
193
194
199 template <typename ScalarType>
200 struct ADNumberInfoFromEnum<
201 ScalarType,
203 std::enable_if_t<std::is_floating_point<ScalarType>::value>>
204 {
205 static const bool is_taped = false;
206 using real_type = Sacado::Fad::DFad<ScalarType>;
207 using derivative_type =
208 typename SacadoNumberInfo<real_type>::derivative_type;
209 static const unsigned int n_supported_derivative_levels =
210 SacadoNumberInfo<real_type>::n_supported_derivative_levels;
211 };
212
213
218 template <typename ScalarType>
219 struct ADNumberInfoFromEnum<
220 ScalarType,
222 std::enable_if_t<std::is_floating_point<ScalarType>::value>>
223 {
224 static const bool is_taped = false;
225 using real_type = Sacado::Fad::DFad<Sacado::Fad::DFad<ScalarType>>;
226 using derivative_type =
227 typename SacadoNumberInfo<real_type>::derivative_type;
228 static const unsigned int n_supported_derivative_levels =
229 SacadoNumberInfo<real_type>::n_supported_derivative_levels;
230 };
231
232
237 template <typename ScalarType>
238 struct ADNumberInfoFromEnum<
239 ScalarType,
241 std::enable_if_t<std::is_floating_point<ScalarType>::value>>
242 {
243 static const bool is_taped = false;
244 using real_type = Sacado::Rad::ADvar<ScalarType>;
245 using derivative_type =
246 typename SacadoNumberInfo<real_type>::derivative_type;
247 static const unsigned int n_supported_derivative_levels =
248 SacadoNumberInfo<real_type>::n_supported_derivative_levels;
249 };
250
251
256 template <typename ScalarType>
257 struct ADNumberInfoFromEnum<
258 ScalarType,
260 std::enable_if_t<std::is_floating_point<ScalarType>::value>>
261 {
262 static const bool is_taped = false;
263 using real_type = Sacado::Rad::ADvar<Sacado::Fad::DFad<ScalarType>>;
264 using derivative_type =
265 typename SacadoNumberInfo<real_type>::derivative_type;
266 static const unsigned int n_supported_derivative_levels =
267 SacadoNumberInfo<real_type>::n_supported_derivative_levels;
268 };
269
270
275 template <typename NumberType>
276 struct Marking<Sacado::Fad::DFad<NumberType>>
277 {
278 using ad_type =
279 typename SacadoNumberInfo<Sacado::Fad::DFad<NumberType>>::ad_type;
280 using derivative_type = typename SacadoNumberInfo<
281 Sacado::Fad::DFad<NumberType>>::derivative_type;
282 using scalar_type =
283 typename SacadoNumberInfo<Sacado::Fad::DFad<NumberType>>::scalar_type;
284
285 /*
286 * Initialize the state of an independent variable.
287 */
288 static void
289 independent_variable(const scalar_type &in,
290 const unsigned int index,
291 const unsigned int n_independent_variables,
292 ad_type & out)
293 {
294 // It is required that we first initialize the outer number before
295 // any of the nested ones.
296 out = ad_type(n_independent_variables, index, in);
297
298 // Initialize potential nested directional derivatives
299 Marking<derivative_type>::independent_variable(
300 in, index, n_independent_variables, out.val());
301 }
302
303 /*
304 * Initialize the state of a dependent variable.
305 */
306 static void
307 dependent_variable(ad_type &out, const ad_type &func)
308 {
309 out = func;
310 }
311 };
312
313
318 template <typename NumberType>
319 struct Marking<Sacado::Rad::ADvar<NumberType>>
320 {
321 using ad_type =
322 typename SacadoNumberInfo<Sacado::Rad::ADvar<NumberType>>::ad_type;
323 using derivative_type = typename SacadoNumberInfo<
324 Sacado::Rad::ADvar<NumberType>>::derivative_type;
325 using scalar_type = typename SacadoNumberInfo<
326 Sacado::Rad::ADvar<NumberType>>::scalar_type;
327
328 /*
329 * Initialize the state of an independent variable.
330 */
331 static void
332 independent_variable(const scalar_type &in,
333 const unsigned int index,
334 const unsigned int n_independent_variables,
335 ad_type & out)
336 {
337 // For Sacado::Rad::ADvar numbers, we have to initialize the
338 // ADNumber with an already fully-configured value. This means
339 // that if this nests another ADNumber then the nested number
340 // must already be set up and ready for use.
341
342 // Initialize potential nested directional derivatives
343 derivative_type derivative_initializer;
344 Marking<derivative_type>::independent_variable(
345 in, index, n_independent_variables, derivative_initializer);
346
347 // Initialize the outer ad_type
348 out = derivative_initializer;
349 }
350
351 /*
352 * Initialize the state of a dependent variable.
353 */
354 static void
355 dependent_variable(ad_type &out, const ad_type &func)
356 {
357 out = func;
358 }
359 };
360
361
367 template <typename NumberType>
368 struct ExtractData<Sacado::Fad::DFad<NumberType>>
369 {
370 using derivative_type = typename SacadoNumberInfo<
371 Sacado::Fad::DFad<NumberType>>::derivative_type;
372 using scalar_type =
373 typename SacadoNumberInfo<Sacado::Fad::DFad<NumberType>>::scalar_type;
374 using value_type =
375 typename SacadoNumberInfo<Sacado::Fad::DFad<NumberType>>::value_type;
376
380 static scalar_type
381 value(const Sacado::Fad::DFad<NumberType> &x)
382 {
383 return ExtractData<value_type>::value(x.val());
384 }
385
386
390 static unsigned int
391 n_directional_derivatives(const Sacado::Fad::DFad<NumberType> &x)
392 {
393 return x.size();
394 }
395
396
400 static derivative_type
401 directional_derivative(const Sacado::Fad::DFad<NumberType> &x,
402 const unsigned int direction)
403 {
404 if (x.hasFastAccess())
405 return x.fastAccessDx(direction);
406 else
407 return x.dx(direction);
408 }
409 };
410
411
417 template <typename NumberType>
418 struct ExtractData<Sacado::Rad::ADvar<NumberType>>
419 {
420 using derivative_type = typename SacadoNumberInfo<
421 Sacado::Rad::ADvar<NumberType>>::derivative_type;
422 using scalar_type = typename SacadoNumberInfo<
423 Sacado::Rad::ADvar<NumberType>>::scalar_type;
424 using value_type =
425 typename SacadoNumberInfo<Sacado::Rad::ADvar<NumberType>>::value_type;
426
430 static scalar_type
431 value(const Sacado::Rad::ADvar<NumberType> &x)
432 {
433 return ExtractData<value_type>::value(x.val());
434 }
435
436
440 static unsigned int
441 n_directional_derivatives(const Sacado::Rad::ADvar<NumberType> &)
442 {
443 // There are as many directional derivatives as there are
444 // independent variables, but each independent variable can
445 // only return one directional derivative.
446 return 1;
447 }
448
449
457 static derivative_type
458 directional_derivative(const Sacado::Rad::ADvar<NumberType> &x,
459 const unsigned int)
460 {
461 return x.adj();
462 }
463 };
464
465 } // namespace internal
466
467
468 /* -------------- NumberTypes::sacado_dfad -------------- */
469
470
477 template <typename ADNumberType>
478 struct ADNumberTraits<
479 ADNumberType,
480 std::enable_if_t<std::is_same<
481 ADNumberType,
482 Sacado::Fad::DFad<typename ADNumberType::scalar_type>>::value>>
483 : NumberTraits<typename ADNumberType::scalar_type,
484 NumberTypes::sacado_dfad>
485 {};
486
487
494 template <typename ADNumberType>
495 struct ADNumberTraits<
496 ADNumberType,
497 std::enable_if_t<std::is_same<
498 ADNumberType,
499 std::complex<Sacado::Fad::DFad<
500 typename ADNumberType::value_type::scalar_type>>>::value>>
501 : NumberTraits<
502 std::complex<typename ADNumberType::value_type::scalar_type>,
503 NumberTypes::sacado_dfad>
504 {};
505
506
511 template <>
512 struct NumberTraits<Sacado::Fad::DFad<float>, NumberTypes::sacado_dfad>
513 : NumberTraits<
514 typename ADNumberTraits<Sacado::Fad::DFad<float>>::scalar_type,
515 NumberTypes::sacado_dfad>
516 {};
517
518
523 template <>
524 struct NumberTraits<std::complex<Sacado::Fad::DFad<float>>,
526 : NumberTraits<typename ADNumberTraits<
527 std::complex<Sacado::Fad::DFad<float>>>::scalar_type,
528 NumberTypes::sacado_dfad>
529 {};
530
531
536 template <>
537 struct NumberTraits<Sacado::Fad::DFad<double>, NumberTypes::sacado_dfad>
538 : NumberTraits<
539 typename ADNumberTraits<Sacado::Fad::DFad<double>>::scalar_type,
540 NumberTypes::sacado_dfad>
541 {};
542
543
548 template <>
549 struct NumberTraits<std::complex<Sacado::Fad::DFad<double>>,
551 : NumberTraits<typename ADNumberTraits<
552 std::complex<Sacado::Fad::DFad<double>>>::scalar_type,
553 NumberTypes::sacado_dfad>
554 {};
555
556
557 /* -------------- NumberTypes::sacado_rad -------------- */
558
559
566 template <typename ADNumberType>
567 struct ADNumberTraits<
568 ADNumberType,
569 std::enable_if_t<std::is_same<
570 ADNumberType,
571 Sacado::Rad::ADvar<typename ADNumberType::ADVari::scalar_type>>::value>>
572 : NumberTraits<typename ADNumberType::ADVari::scalar_type,
573 NumberTypes::sacado_rad>
574 {};
575
576
581 template <>
582 struct NumberTraits<Sacado::Rad::ADvar<float>, NumberTypes::sacado_rad>
583 : NumberTraits<
584 typename ADNumberTraits<Sacado::Rad::ADvar<float>>::scalar_type,
585 NumberTypes::sacado_rad>
586 {};
587
588
593 template <>
594 struct NumberTraits<Sacado::Rad::ADvar<double>, NumberTypes::sacado_rad>
595 : NumberTraits<
596 typename ADNumberTraits<Sacado::Rad::ADvar<double>>::scalar_type,
597 NumberTypes::sacado_rad>
598 {};
599
600
601# ifdef DEAL_II_TRILINOS_CXX_SUPPORTS_SACADO_COMPLEX_RAD
602
603
610 template <typename ADNumberType>
611 struct ADNumberTraits<
612 ADNumberType,
613 std::enable_if_t<std::is_same<
614 ADNumberType,
615 std::complex<Sacado::Rad::ADvar<
616 typename ADNumberType::value_type::ADVari::scalar_type>>>::value>>
617 : NumberTraits<
618 std::complex<typename ADNumberType::value_type::ADVari::scalar_type>,
619 NumberTypes::sacado_rad>
620 {};
621
622
627 template <>
628 struct NumberTraits<std::complex<Sacado::Rad::ADvar<float>>,
630 : NumberTraits<typename ADNumberTraits<
631 std::complex<Sacado::Rad::ADvar<float>>>::scalar_type,
632 NumberTypes::sacado_rad>
633 {};
634
635
640 template <>
641 struct NumberTraits<std::complex<Sacado::Rad::ADvar<double>>,
643 : NumberTraits<typename ADNumberTraits<
644 std::complex<Sacado::Rad::ADvar<double>>>::scalar_type,
645 NumberTypes::sacado_rad>
646 {};
647
648
649# endif
650
651
652 /* -------------- NumberTypes::sacado_dfad_dfad -------------- */
653
661 template <typename ADNumberType>
662 struct ADNumberTraits<ADNumberType,
663 std::enable_if_t<std::is_same<
664 ADNumberType,
665 Sacado::Fad::DFad<Sacado::Fad::DFad<
666 typename ADNumberType::scalar_type>>>::value>>
667 : NumberTraits<typename ADNumberType::scalar_type,
668 NumberTypes::sacado_dfad_dfad>
669 {};
670
671
679 template <typename ADNumberType>
680 struct ADNumberTraits<
681 ADNumberType,
682 std::enable_if_t<std::is_same<
683 ADNumberType,
684 std::complex<Sacado::Fad::DFad<Sacado::Fad::DFad<
685 typename ADNumberType::value_type::scalar_type>>>>::value>>
686 : NumberTraits<
687 std::complex<typename ADNumberType::value_type::scalar_type>,
688 NumberTypes::sacado_dfad_dfad>
689 {};
690
691
696 template <>
697 struct NumberTraits<Sacado::Fad::DFad<Sacado::Fad::DFad<float>>,
699 : NumberTraits<typename ADNumberTraits<Sacado::Fad::DFad<
700 Sacado::Fad::DFad<float>>>::scalar_type,
701 NumberTypes::sacado_dfad_dfad>
702 {};
703
704
709 template <>
710 struct NumberTraits<
711 std::complex<Sacado::Fad::DFad<Sacado::Fad::DFad<float>>>,
713 : NumberTraits<typename ADNumberTraits<std::complex<Sacado::Fad::DFad<
714 Sacado::Fad::DFad<float>>>>::scalar_type,
715 NumberTypes::sacado_dfad_dfad>
716 {};
717
718
723 template <>
724 struct NumberTraits<Sacado::Fad::DFad<Sacado::Fad::DFad<double>>,
726 : NumberTraits<typename ADNumberTraits<Sacado::Fad::DFad<
727 Sacado::Fad::DFad<double>>>::scalar_type,
728 NumberTypes::sacado_dfad_dfad>
729 {};
730
731
736 template <>
737 struct NumberTraits<
738 std::complex<Sacado::Fad::DFad<Sacado::Fad::DFad<double>>>,
740 : NumberTraits<typename ADNumberTraits<std::complex<Sacado::Fad::DFad<
741 Sacado::Fad::DFad<double>>>>::scalar_type,
742 NumberTypes::sacado_dfad_dfad>
743 {};
744
745
746 /* -------------- NumberTypes::sacado_rad_dfad -------------- */
747
755 template <typename ADNumberType>
756 struct ADNumberTraits<
757 ADNumberType,
758 std::enable_if_t<
759 std::is_same<ADNumberType,
760 Sacado::Rad::ADvar<Sacado::Fad::DFad<
761 typename ADNumberType::ADVari::scalar_type>>>::value>>
762 : NumberTraits<typename ADNumberType::ADVari::scalar_type,
763 NumberTypes::sacado_rad_dfad>
764 {};
765
766
771 template <>
772 struct NumberTraits<Sacado::Rad::ADvar<Sacado::Fad::DFad<float>>,
774 : NumberTraits<typename ADNumberTraits<Sacado::Rad::ADvar<
775 Sacado::Fad::DFad<float>>>::scalar_type,
776 NumberTypes::sacado_rad_dfad>
777 {};
778
779
784 template <>
785 struct NumberTraits<Sacado::Rad::ADvar<Sacado::Fad::DFad<double>>,
787 : NumberTraits<typename ADNumberTraits<Sacado::Rad::ADvar<
788 Sacado::Fad::DFad<double>>>::scalar_type,
789 NumberTypes::sacado_rad_dfad>
790 {};
791
792
793# ifdef DEAL_II_TRILINOS_CXX_SUPPORTS_SACADO_COMPLEX_RAD
794
795
803 template <typename ADNumberType>
804 struct ADNumberTraits<
805 ADNumberType,
806 std::enable_if_t<std::is_same<
807 ADNumberType,
808 std::complex<Sacado::Rad::ADvar<Sacado::Fad::DFad<
809 typename ADNumberType::value_type::ADVari::scalar_type>>>>::value>>
810 : NumberTraits<
811 std::complex<typename ADNumberType::value_type::ADVari::scalar_type>,
812 NumberTypes::sacado_rad_dfad>
813 {};
814
815
820 template <>
821 struct NumberTraits<
822 std::complex<Sacado::Rad::ADvar<Sacado::Fad::DFad<float>>>,
824 : NumberTraits<typename ADNumberTraits<std::complex<Sacado::Rad::ADvar<
825 Sacado::Fad::DFad<float>>>>::scalar_type,
826 NumberTypes::sacado_rad_dfad>
827 {};
828
829
834 template <>
835 struct NumberTraits<
836 std::complex<Sacado::Rad::ADvar<Sacado::Fad::DFad<double>>>,
838 : NumberTraits<typename ADNumberTraits<std::complex<Sacado::Rad::ADvar<
839 Sacado::Fad::DFad<double>>>>::scalar_type,
840 NumberTypes::sacado_rad_dfad>
841 {};
842
843
844# endif
845
846
847 /* -------------- Additional type traits -------------- */
848
849
850 template <typename NumberType>
851 struct is_sacado_dfad_number<
852 NumberType,
853 std::enable_if_t<
854 ADNumberTraits<typename std::decay<NumberType>::type>::type_code ==
855 NumberTypes::sacado_dfad ||
856 ADNumberTraits<typename std::decay<NumberType>::type>::type_code ==
857 NumberTypes::sacado_dfad_dfad>> : std::true_type
858 {};
859
860
861 template <typename NumberType>
862 struct is_sacado_dfad_number<
863 NumberType,
864 std::enable_if_t<std::is_same<
865 NumberType,
866 Sacado::Fad::Expr<typename NumberType::value_type>>::value>>
867 : std::true_type
868 {};
869
870
871 template <typename NumberType>
872 struct is_sacado_rad_number<
873 NumberType,
874 std::enable_if_t<
875 ADNumberTraits<typename std::decay<NumberType>::type>::type_code ==
876 NumberTypes::sacado_rad ||
877 ADNumberTraits<typename std::decay<NumberType>::type>::type_code ==
878 NumberTypes::sacado_rad_dfad>> : std::true_type
879 {};
880
881
882 template <typename NumberType>
883 struct is_sacado_rad_number<
884 NumberType,
885 std::enable_if_t<
886 std::is_same<NumberType,
887 Sacado::Rad::ADvari<Sacado::Fad::DFad<
888 typename NumberType::ADVari::scalar_type>>>::value>>
889 : std::true_type
890 {};
891
892
893 template <typename NumberType>
894 struct is_sacado_number<
895 NumberType,
896 std::enable_if_t<is_sacado_dfad_number<NumberType>::value ||
897 is_sacado_rad_number<NumberType>::value>>
898 : std::true_type
899 {};
900
901 } // namespace AD
902} // namespace Differentiation
903
904
905# endif // DOXYGEN
906
907
908
910
911
912#endif // DEAL_II_TRILINOS_WITH_SACADO
913
914#endif
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
STL namespace.