Reference documentation for deal.II version 9.4.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\}}\)
ad_helpers.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2016 - 2022 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 at
12// the top level of the deal.II distribution.
13//
14// ---------------------------------------------------------------------
15
16#ifndef dealii_differentiation_ad_ad_helpers_h
17#define dealii_differentiation_ad_ad_helpers_h
18
19#include <deal.II/base/config.h>
20
21#if defined(DEAL_II_WITH_ADOLC) || defined(DEAL_II_TRILINOS_WITH_SACADO)
22
23# include <deal.II/base/numbers.h>
25# include <deal.II/base/tensor.h>
26
33
35
37# include <deal.II/lac/vector.h>
38
39# include <algorithm>
40# include <iostream>
41# include <iterator>
42# include <numeric>
43# include <set>
44
46
47namespace Differentiation
48{
49 namespace AD
50 {
168 template <enum AD::NumberTypes ADNumberTypeCode,
169 typename ScalarType = double>
171 {
172 public:
179
184 using ad_type =
186
191
206 HelperBase(const unsigned int n_independent_variables,
207 const unsigned int n_dependent_variables);
208
212 virtual ~HelperBase() = default;
213
215
220
225 std::size_t
227
232 std::size_t
233 n_dependent_variables() const;
234
243 void
244 print(std::ostream &stream) const;
245
252 void
253 print_values(std::ostream &stream) const;
254
265 void
266 print_tape_stats(const typename Types<ad_type>::tape_index tape_index,
267 std::ostream & stream) const;
268
270
275
294 static void
296 const bool ensure_persistent_setting = true);
297
299
304
340 virtual void
341 reset(const unsigned int n_independent_variables =
343 const unsigned int n_dependent_variables =
345 const bool clear_registered_tapes = true);
346
351 bool
352 is_recording() const;
353
359 active_tape_index() const;
360
365 bool
367 const typename Types<ad_type>::tape_index tape_index) const;
368
394 void
396 const typename Types<ad_type>::tape_buffer_sizes obufsize = 64 * 1024 *
397 1024,
398 const typename Types<ad_type>::tape_buffer_sizes lbufsize = 64 * 1024 *
399 1024,
400 const typename Types<ad_type>::tape_buffer_sizes vbufsize = 64 * 1024 *
401 1024,
402 const typename Types<ad_type>::tape_buffer_sizes tbufsize = 64 * 1024 *
403 1024);
404
450 bool
452 const typename Types<ad_type>::tape_index tape_index,
453 const bool overwrite_tape = false,
454 const bool keep_independent_values = true);
455
467 void
468 stop_recording_operations(const bool write_tapes_to_file = false);
469
479 void
481 const typename Types<ad_type>::tape_index tape_index);
482
522 bool
524 const typename Types<ad_type>::tape_index tape_index) const;
525
565 bool
567
575 void
577
579
580 protected:
585
593
601
618 void
619 activate_tape(const typename Types<ad_type>::tape_index tape_index,
620 const bool read_mode);
621
623
628
636 mutable std::vector<scalar_type> independent_variable_values;
637
646 mutable std::vector<ad_type> independent_variables;
647
653
659
665 void
667
674 void
675 set_sensitivity_value(const unsigned int index, const scalar_type &value);
676
695 void
696 mark_independent_variable(const unsigned int index, ad_type &out) const;
697
708 void
710
720 void
721 initialize_non_sensitive_independent_variable(const unsigned int index,
722 ad_type &out) const;
723
728 unsigned int
730
732
737
750 std::vector<ad_type> dependent_variables;
751
756
763 void
764 reset_registered_dependent_variables(const bool flag = false);
765
769 unsigned int
771
783 void
784 register_dependent_variable(const unsigned int index,
785 const ad_type & func);
786
788
789 }; // class HelperBase
790
791
792
834 template <enum AD::NumberTypes ADNumberTypeCode,
835 typename ScalarType = double>
836 class CellLevelBase : public HelperBase<ADNumberTypeCode, ScalarType>
837 {
838 public:
845
850 using ad_type =
852
857
872 CellLevelBase(const unsigned int n_independent_variables,
873 const unsigned int n_dependent_variables);
874
878 virtual ~CellLevelBase() = default;
879
881
886
904 void
905 register_dof_values(const std::vector<scalar_type> &dof_values);
906
925 template <typename VectorType>
926 void
928 const VectorType & values,
929 const std::vector<::types::global_dof_index> &local_dof_indices);
930
951 const std::vector<ad_type> &
953
955
960
980 void
981 set_dof_values(const std::vector<scalar_type> &dof_values);
982
1000 template <typename VectorType>
1001 void
1003 const VectorType & values,
1004 const std::vector<::types::global_dof_index> &local_dof_indices);
1005
1007
1012
1027 virtual void
1029
1047 virtual void
1049
1051
1052 }; // class CellLevelBase
1053
1054
1055
1216 template <enum AD::NumberTypes ADNumberTypeCode,
1217 typename ScalarType = double>
1218 class EnergyFunctional : public CellLevelBase<ADNumberTypeCode, ScalarType>
1219 {
1220 public:
1227
1232 using ad_type =
1234
1239
1255 EnergyFunctional(const unsigned int n_independent_variables);
1256
1260 virtual ~EnergyFunctional() = default;
1261
1263
1268
1282 void
1283 register_energy_functional(const ad_type &energy);
1284
1299 compute_energy() const;
1300
1320 void
1321 compute_residual(Vector<scalar_type> &residual) const override;
1322
1343 virtual void
1345 FullMatrix<scalar_type> &linearization) const override;
1346
1348
1349 }; // class EnergyFunctional
1350
1351
1352
1527 template <enum AD::NumberTypes ADNumberTypeCode,
1528 typename ScalarType = double>
1530 : public CellLevelBase<ADNumberTypeCode, ScalarType>
1531 {
1532 public:
1539
1544 using ad_type =
1546
1551
1567 const unsigned int n_dependent_variables);
1568
1572 virtual ~ResidualLinearization() = default;
1573
1575
1580
1594 void
1595 register_residual_vector(const std::vector<ad_type> &residual);
1596
1613 virtual void
1614 compute_residual(Vector<scalar_type> &residual) const override;
1615
1633 virtual void
1635 FullMatrix<scalar_type> &linearization) const override;
1636
1638
1639 }; // class ResidualLinearization
1640
1641
1642
1643 namespace internal
1644 {
1649 template <int dim, typename ExtractorType>
1651
1652
1658 template <int dim>
1659 struct Extractor<dim, FEValuesExtractors::Scalar>
1660 {
1664 static const unsigned int n_components = 1;
1665
1669 static const unsigned int rank = 0;
1670
1674 template <typename NumberType>
1676
1677 static_assert(
1679 "The number of components doesn't match that of the corresponding tensor type.");
1680 static_assert(
1682 "The rank doesn't match that of the corresponding tensor type.");
1683
1687 // Note: FEValuesViews::Scalar::tensor_type is double, so we can't
1688 // use it (FEValuesViews) in this context.
1689 // In fact, sadly, all FEValuesViews objects expect doubles as value
1690 // types.
1691 template <typename NumberType>
1693
1697 template <typename NumberType>
1699
1703 static inline unsigned int
1705 {
1706 return extractor.component;
1707 }
1708
1716 static bool
1717 symmetric_component(const unsigned int unrolled_index)
1718 {
1719 (void)unrolled_index;
1720 return false;
1721 }
1722
1730 template <typename IndexType = unsigned int, int rank_in>
1731 static IndexType
1733 const unsigned int column_offset)
1734 {
1735 Assert(column_offset <= rank_in, ExcInternalError());
1736 (void)table_indices;
1737 (void)column_offset;
1738 return 0;
1739 }
1740 };
1741
1742
1748 template <int dim>
1750 {
1754 static const unsigned int n_components = dim;
1755
1759 static const unsigned int rank = 1;
1760
1764 template <typename NumberType>
1766
1767 static_assert(
1769 "The number of components doesn't match that of the corresponding tensor type.");
1770 static_assert(
1772 "The rank doesn't match that of the corresponding tensor type.");
1773
1777 template <typename NumberType>
1779
1783 template <typename NumberType>
1785
1789 static inline unsigned int
1791 {
1792 return extractor.first_vector_component;
1793 }
1794
1803 static bool
1804 symmetric_component(const unsigned int unrolled_index)
1805 {
1806 (void)unrolled_index;
1807 return false;
1808 }
1809
1814 template <int rank_in>
1815 static TableIndices<rank>
1817 const unsigned int column_offset)
1818 {
1819 Assert(0 + column_offset < rank_in, ExcInternalError());
1820 return TableIndices<rank>(table_indices[column_offset]);
1821 }
1822
1838 template <typename IndexType = unsigned int, int rank_in>
1839 static IndexType
1841 const unsigned int column_offset)
1842 {
1843 static_assert(
1844 rank_in >= rank,
1845 "Cannot extract more table indices than the input table has!");
1846 using TensorType = tensor_type<double>;
1847 return TensorType::component_to_unrolled_index(
1848 table_index_view(table_indices, column_offset));
1849 }
1850 };
1851
1852
1858 template <int dim>
1860 {
1864 static const unsigned int n_components =
1866
1870 static const unsigned int rank = 1;
1871
1875 template <typename NumberType>
1877
1881 template <typename NumberType>
1883
1887 template <typename NumberType>
1889
1893 static inline unsigned int
1895 {
1896 return extractor.first_tensor_component;
1897 }
1898
1906 static bool
1907 symmetric_component(const unsigned int unrolled_index)
1908 {
1909 (void)unrolled_index;
1910 return false;
1911 }
1912
1917 template <int rank_in>
1918 static TableIndices<rank>
1920 const unsigned int column_offset)
1921 {
1922 Assert(column_offset < rank_in, ExcInternalError());
1923 return TableIndices<rank>(table_indices[column_offset]);
1924 }
1925
1941 template <typename IndexType = unsigned int, int rank_in>
1942 static IndexType
1944 const unsigned int column_offset)
1945 {
1946 static_assert(
1947 rank_in >= rank,
1948 "Cannot extract more table indices than the input table has!");
1949 using TensorType = tensor_type<double>;
1950 return TensorType::component_to_unrolled_index(
1951 table_index_view(table_indices, column_offset));
1952 }
1953 };
1954
1955
1961 template <int dim>
1963 {
1967 static const unsigned int n_components =
1969
1973 static const unsigned int rank = Tensor<2, dim>::rank;
1974
1978 template <typename NumberType>
1980
1984 template <typename NumberType>
1986
1990 template <typename NumberType>
1992
1996 static inline unsigned int
1998 {
1999 return extractor.first_tensor_component;
2000 }
2001
2009 static bool
2010 symmetric_component(const unsigned int unrolled_index)
2011 {
2012 (void)unrolled_index;
2013 return false;
2014 }
2015
2020 template <int rank_in>
2021 static TableIndices<rank>
2023 const unsigned int column_offset)
2024 {
2025 Assert(column_offset < rank_in, ExcInternalError());
2026 Assert(column_offset + 1 < rank_in, ExcInternalError());
2027 return TableIndices<rank>(table_indices[column_offset],
2028 table_indices[column_offset + 1]);
2029 }
2030
2046 template <typename IndexType = unsigned int, int rank_in>
2047 static IndexType
2049 const unsigned int column_offset)
2050 {
2051 static_assert(
2052 rank_in >= rank,
2053 "Cannot extract more table indices than the input table has!");
2054 using TensorType = tensor_type<double>;
2055 return TensorType::component_to_unrolled_index(
2056 table_index_view(table_indices, column_offset));
2057 }
2058 };
2059
2060
2066 template <int dim>
2068 {
2072 static const unsigned int n_components =
2074
2078 static const unsigned int rank = SymmetricTensor<2, dim>::rank;
2079
2083 template <typename NumberType>
2085
2089 template <typename NumberType>
2091
2095 template <typename NumberType>
2097
2101 static inline unsigned int
2103 {
2104 return extractor.first_tensor_component;
2105 }
2106
2115 static bool
2116 symmetric_component(const unsigned int unrolled_index)
2117 {
2118 const TableIndices<2> table_indices =
2120 return table_indices[0] != table_indices[1];
2121 }
2122
2127 template <int rank_in>
2128 static TableIndices<rank>
2130 const unsigned int column_offset)
2131 {
2132 Assert(column_offset < rank_in, ExcInternalError());
2133 Assert(column_offset + 1 < rank_in, ExcInternalError());
2134 return TableIndices<rank>(table_indices[column_offset],
2135 table_indices[column_offset + 1]);
2136 }
2137
2153 template <typename IndexType = unsigned int, int rank_in>
2154 static IndexType
2156 const unsigned int column_offset)
2157 {
2158 static_assert(
2159 rank_in >= rank,
2160 "Cannot extract more table indices than the input table has!");
2161 using TensorType = tensor_type<double>;
2162 return TensorType::component_to_unrolled_index(
2163 table_index_view(table_indices, column_offset));
2164 }
2165 };
2166
2167
2173 template <int dim, typename NumberType, typename ExtractorType>
2175 {
2180 using type =
2181 typename Extractor<dim,
2182 ExtractorType>::template tensor_type<NumberType>;
2183 };
2184
2185
2195 template <typename ExtractorType_Row, typename ExtractorType_Col>
2197 {
2209 template <int rank, int dim, typename NumberType>
2211 };
2212
2213
2222 template <>
2225 {
2233 template <int /*rank*/, int dim, typename NumberType>
2234 using type = SymmetricTensor<2 /*rank*/, dim, NumberType>;
2235 };
2236
2237
2246 template <>
2249 {
2258 template <int /*rank*/, int dim, typename NumberType>
2259 using type = SymmetricTensor<2 /*rank*/, dim, NumberType>;
2260 };
2261
2262
2270 template <>
2273 {
2282 template <int /*rank*/, int dim, typename NumberType>
2283 using type = SymmetricTensor<4 /*rank*/, dim, NumberType>;
2284 };
2285
2286
2296 template <int dim,
2297 typename NumberType,
2298 typename ExtractorType_Row,
2299 typename ExtractorType_Col>
2301 {
2307
2314 using type =
2317 };
2318
2319
2324 template <int dim, typename NumberType, typename ExtractorType_Field>
2327
2328
2338 template <int dim,
2339 typename NumberType,
2340 typename ExtractorType_Field,
2341 typename ExtractorType_Derivative>
2343 NumberType,
2344 ExtractorType_Field,
2345 ExtractorType_Derivative>;
2346
2347
2353 template <int dim,
2354 typename IndexType = unsigned int,
2355 typename ExtractorType>
2356 std::vector<IndexType>
2357 extract_field_component_indices(const ExtractorType &extractor,
2358 const bool ignore_symmetries = true)
2359 {
2360 (void)ignore_symmetries;
2361 const IndexType n_components =
2363 const IndexType comp_first =
2365 std::vector<IndexType> indices(n_components);
2366 std::iota(indices.begin(), indices.end(), comp_first);
2367 return indices;
2368 }
2369
2370
2379 template <int dim, typename IndexType = unsigned int>
2380 std::vector<IndexType>
2382 const FEValuesExtractors::SymmetricTensor<2> &extractor_symm_tensor,
2383 const bool ignore_symmetries = true)
2384 {
2385 using ExtractorType = FEValuesExtractors::SymmetricTensor<2>;
2386 const IndexType n_components =
2388 const IndexType comp_first =
2390 extractor_symm_tensor);
2391
2392 if (ignore_symmetries == true)
2393 {
2394 std::vector<IndexType> indices(n_components);
2395 std::iota(indices.begin(), indices.end(), comp_first);
2396 return indices;
2397 }
2398 else
2399 {
2400 // First get all of the indices of the non-symmetric tensor
2401 const FEValuesExtractors::Tensor<2> extractor_tensor(
2402 extractor_symm_tensor.first_tensor_component);
2403 std::vector<IndexType> indices =
2404 extract_field_component_indices<dim>(extractor_tensor, true);
2405
2406 // Then we overwrite any illegal entries with the equivalent indices
2407 // from the symmetric tensor
2408 for (unsigned int i = 0; i < indices.size(); ++i)
2409 {
2410 // The indices stored in the vector start with the extractor's
2411 // first_component_index. We need to account for this when
2412 // retrieving the tensor (local) index.
2413 const IndexType local_index_i = indices[i] - comp_first;
2414 if (local_index_i >= n_components)
2415 {
2416 const TableIndices<2> ti_tensor =
2418 local_index_i);
2419 const IndexType sti_new_index =
2421 ti_tensor);
2422 indices[i] = comp_first + sti_new_index;
2423 }
2424 }
2425
2426 return indices;
2427 }
2428 }
2429
2430
2435 template <typename TensorType, typename NumberType>
2436 inline void
2437 set_tensor_entry(TensorType & t,
2438 const unsigned int unrolled_index,
2439 const NumberType & value)
2440 {
2441 // Where possible, set values using TableIndices
2442 AssertIndexRange(unrolled_index, t.n_independent_components);
2443 t[TensorType::unrolled_to_component_indices(unrolled_index)] = value;
2444 }
2445
2446
2451 template <int dim, typename NumberType>
2452 inline void
2454 const unsigned int unrolled_index,
2455 const NumberType & value)
2456 {
2457 AssertIndexRange(unrolled_index, 1);
2458 (void)unrolled_index;
2459 t = value;
2460 }
2461
2462
2468 template <typename NumberType>
2469 inline void
2471 const unsigned int unrolled_index,
2472 const NumberType & value)
2473 {
2474 AssertIndexRange(unrolled_index, 1);
2475 (void)unrolled_index;
2476 t = value;
2477 }
2478
2479
2485 template <int dim, typename NumberType>
2486 inline void
2488 const unsigned int unrolled_index_row,
2489 const unsigned int unrolled_index_col,
2490 const NumberType & value)
2491 {
2492 // Fourth order symmetric tensors require a specialized interface
2493 // to extract values.
2494 using SubTensorType = SymmetricTensor<2, dim, NumberType>;
2495 AssertIndexRange(unrolled_index_row,
2496 SubTensorType::n_independent_components);
2497 AssertIndexRange(unrolled_index_col,
2498 SubTensorType::n_independent_components);
2499 const TableIndices<2> indices_row =
2500 SubTensorType::unrolled_to_component_indices(unrolled_index_row);
2501 const TableIndices<2> indices_col =
2502 SubTensorType::unrolled_to_component_indices(unrolled_index_col);
2503 t[indices_row[0]][indices_row[1]][indices_col[0]][indices_col[1]] =
2504 value;
2505 }
2506
2507
2512 template <int rank,
2513 int dim,
2514 typename NumberType,
2515 template <int, int, typename>
2516 class TensorType>
2517 inline NumberType
2518 get_tensor_entry(const TensorType<rank, dim, NumberType> &t,
2519 const unsigned int unrolled_index)
2520 {
2521 // Where possible, get values using TableIndices
2522 AssertIndexRange(unrolled_index, t.n_independent_components);
2523 return t[TensorType<rank, dim, NumberType>::
2524 unrolled_to_component_indices(unrolled_index)];
2525 }
2526
2527
2532 template <int dim,
2533 typename NumberType,
2534 template <int, int, typename>
2535 class TensorType>
2536 inline NumberType
2537 get_tensor_entry(const TensorType<0, dim, NumberType> &t,
2538 const unsigned int unrolled_index)
2539 {
2540 AssertIndexRange(unrolled_index, 1);
2541 (void)unrolled_index;
2542 return t;
2543 }
2544
2545
2551 template <typename NumberType>
2552 inline const NumberType &
2553 get_tensor_entry(const NumberType &t, const unsigned int unrolled_index)
2554 {
2555 AssertIndexRange(unrolled_index, 1);
2556 (void)unrolled_index;
2557 return t;
2558 }
2559
2560
2565 template <int rank,
2566 int dim,
2567 typename NumberType,
2568 template <int, int, typename>
2569 class TensorType>
2570 inline NumberType &
2571 get_tensor_entry(TensorType<rank, dim, NumberType> &t,
2572 const unsigned int unrolled_index)
2573 {
2574 // Where possible, get values using TableIndices
2575 AssertIndexRange(unrolled_index, t.n_independent_components);
2576 return t[TensorType<rank, dim, NumberType>::
2577 unrolled_to_component_indices(unrolled_index)];
2578 }
2579
2580
2585 template <int dim,
2586 typename NumberType,
2587 template <int, int, typename>
2588 class TensorType>
2589 NumberType &
2590 get_tensor_entry(TensorType<0, dim, NumberType> &t,
2591 const unsigned int index)
2592 {
2593 AssertIndexRange(index, 1);
2594 (void)index;
2595 return t;
2596 }
2597
2598
2604 template <typename NumberType>
2605 inline NumberType &
2606 get_tensor_entry(NumberType &t, const unsigned int index)
2607 {
2608 AssertIndexRange(index, 1);
2609 (void)index;
2610 return t;
2611 }
2612
2613 } // namespace internal
2614
2615
2616
2638 template <int dim,
2639 enum AD::NumberTypes ADNumberTypeCode,
2640 typename ScalarType = double>
2642 : public HelperBase<ADNumberTypeCode, ScalarType>
2643 {
2644 public:
2649 static constexpr unsigned int dimension = dim;
2650
2657
2662 using ad_type =
2664
2669
2685 const unsigned int n_dependent_variables);
2686
2690 virtual ~PointLevelFunctionsBase() = default;
2691
2693
2698
2702 virtual void
2703 reset(const unsigned int n_independent_variables =
2705 const unsigned int n_dependent_variables =
2707 const bool clear_registered_tapes = true) override;
2708
2724 void
2725 register_independent_variables(const std::vector<scalar_type> &values);
2726
2755 template <typename ValueType, typename ExtractorType>
2756 void
2757 register_independent_variable(const ValueType & value,
2758 const ExtractorType &extractor);
2759
2780 const std::vector<ad_type> &
2782
2783 /*
2784 * Extract a subset of the independent variables as represented by
2785 * auto-differentiable numbers. These are the independent
2786 * variables @f$\mathbf{A} \subset \mathbf{X}@f$ at which the dependent values
2787 * are evaluated and differentiated.
2788 *
2789 * This function indicates to the AD library that implements the
2790 * auto-differentiable number type that operations performed on these
2791 * numbers are to be tracked so they are considered "sensitive"
2792 * variables. This is, therefore, the set of variables with which one
2793 * would then perform computations, and based on which one can then
2794 * extract both the value of the function and its derivatives with the
2795 * member functions below. The values of the components of the returned
2796 * object are initialized to the values set with
2797 * register_independent_variable().
2798 *
2799 * @param[in] extractor An extractor associated with the input field
2800 * variables. This effectively defines which components of the global set
2801 * of independent variables this field is associated with.
2802 * @return An object of auto-differentiable type numbers. The return type is
2803 * based on the input extractor, and will be either a scalar,
2804 * Tensor<1,dim>, Tensor<2,dim>, or SymmetricTensor<2,dim>.
2805 *
2806 * @note For taped AD numbers, this operation is only valid in recording mode.
2807 */
2808 template <typename ExtractorType>
2809 typename internal::Extractor<dim,
2810 ExtractorType>::template tensor_type<ad_type>
2811 get_sensitive_variables(const ExtractorType &extractor) const;
2812
2814
2819
2837 void
2838 set_independent_variables(const std::vector<scalar_type> &values);
2839
2866 template <typename ValueType, typename ExtractorType>
2867 void
2868 set_independent_variable(const ValueType & value,
2869 const ExtractorType &extractor);
2870
2872
2873 protected:
2878
2889 void
2890 set_sensitivity_value(const unsigned int index,
2891 const bool symmetric_component,
2892 const scalar_type &value);
2893
2899 bool
2900 is_symmetric_independent_variable(const unsigned int index) const;
2901
2906 unsigned int
2908
2910
2911 private:
2916
2922
2924
2925 }; // class PointLevelFunctionsBase
2926
2927
2928
3083 template <int dim,
3084 enum AD::NumberTypes ADNumberTypeCode,
3085 typename ScalarType = double>
3087 : public PointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>
3088 {
3089 public:
3096
3101 using ad_type =
3103
3108
3118 ScalarFunction(const unsigned int n_independent_variables);
3119
3123 virtual ~ScalarFunction() = default;
3124
3126
3131
3143 void
3145
3154 compute_value() const;
3155
3168 void
3169 compute_gradient(Vector<scalar_type> &gradient) const;
3170
3185 void
3187
3222 template <typename ExtractorType_Row>
3223 static typename internal::
3226 const ExtractorType_Row & extractor_row);
3227
3266 template <typename ExtractorType_Row, typename ExtractorType_Col>
3267 static typename internal::ScalarFieldHessian<dim,
3269 ExtractorType_Row,
3270 ExtractorType_Col>::type
3272 const ExtractorType_Row & extractor_row,
3273 const ExtractorType_Col & extractor_col);
3274
3290 const FullMatrix<scalar_type> & hessian,
3291 const FEValuesExtractors::Scalar &extractor_row,
3292 const FEValuesExtractors::Scalar &extractor_col);
3293
3307 const FullMatrix<scalar_type> & hessian,
3308 const FEValuesExtractors::SymmetricTensor<2> &extractor_row,
3309 const FEValuesExtractors::SymmetricTensor<2> &extractor_col);
3310
3312
3313 }; // class ScalarFunction
3314
3315
3316
3474 template <int dim,
3475 enum AD::NumberTypes ADNumberTypeCode,
3476 typename ScalarType = double>
3478 : public PointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>
3479 {
3480 public:
3487
3492 using ad_type =
3494
3499
3514 VectorFunction(const unsigned int n_independent_variables,
3515 const unsigned int n_dependent_variables);
3516
3520 virtual ~VectorFunction() = default;
3521
3523
3528
3541 void
3542 register_dependent_variables(const std::vector<ad_type> &funcs);
3543
3562 template <typename ValueType, typename ExtractorType>
3563 void
3564 register_dependent_variable(const ValueType & funcs,
3565 const ExtractorType &extractor);
3566
3575 void
3577
3592 void
3594
3595
3608 template <typename ExtractorType_Row>
3609 static typename internal::
3612 const ExtractorType_Row & extractor_row);
3613
3661 template <typename ExtractorType_Row, typename ExtractorType_Col>
3662 static typename internal::VectorFieldJacobian<dim,
3664 ExtractorType_Row,
3665 ExtractorType_Col>::type
3667 const ExtractorType_Row & extractor_row,
3668 const ExtractorType_Col & extractor_col);
3669
3687 const FullMatrix<scalar_type> & jacobian,
3688 const FEValuesExtractors::Scalar &extractor_row,
3689 const FEValuesExtractors::Scalar &extractor_col);
3690
3706 const FullMatrix<scalar_type> & jacobian,
3707 const FEValuesExtractors::SymmetricTensor<2> &extractor_row,
3708 const FEValuesExtractors::SymmetricTensor<2> &extractor_col);
3709
3711
3712 }; // class VectorFunction
3713
3714
3715 } // namespace AD
3716} // namespace Differentiation
3717
3718
3719/* ----------------- inline and template functions ----------------- */
3720
3721
3722# ifndef DOXYGEN
3723
3724namespace Differentiation
3725{
3726 namespace AD
3727 {
3728 /* ----------------- CellLevelBase ----------------- */
3729
3730
3731
3732 template <enum AD::NumberTypes ADNumberTypeCode, typename ScalarType>
3733 template <typename VectorType>
3734 void
3736 const VectorType & values,
3737 const std::vector<::types::global_dof_index> &local_dof_indices)
3738 {
3739 // This is actually the same thing the set_dof_values() function,
3740 // in the sense that we simply populate our array of independent values
3741 // with a meaningful number. However, in this case we need to double check
3742 // that we're not registering these variables twice
3743 Assert(
3744 local_dof_indices.size() == this->n_independent_variables(),
3745 ExcMessage(
3746 "Degree of freedom index vector size does not match number of independent variables"));
3747 for (unsigned int i = 0; i < this->n_independent_variables(); ++i)
3748 {
3749 Assert(this->registered_independent_variable_values[i] == false,
3750 ExcMessage("Independent variables already registered."));
3751 }
3752 set_dof_values(values, local_dof_indices);
3753 }
3754
3755
3756
3757 template <enum AD::NumberTypes ADNumberTypeCode, typename ScalarType>
3758 template <typename VectorType>
3759 void
3761 const VectorType & values,
3762 const std::vector<::types::global_dof_index> &local_dof_indices)
3763 {
3764 Assert(local_dof_indices.size() == this->n_independent_variables(),
3765 ExcMessage(
3766 "Vector size does not match number of independent variables"));
3767 for (unsigned int i = 0; i < this->n_independent_variables(); ++i)
3769 i, values[local_dof_indices[i]]);
3770 }
3771
3772
3773
3774 /* ----------------- PointLevelFunctionsBase ----------------- */
3775
3776
3777
3778 template <int dim,
3779 enum AD::NumberTypes ADNumberTypeCode,
3780 typename ScalarType>
3781 template <typename ValueType, typename ExtractorType>
3782 void
3784 register_independent_variable(const ValueType & value,
3785 const ExtractorType &extractor)
3786 {
3787 // This is actually the same thing as the set_independent_variable
3788 // function, in the sense that we simply populate our array of independent
3789 // values with a meaningful number. However, in this case we need to
3790 // double check that we're not registering these variables twice
3791# ifdef DEBUG
3792 const std::vector<unsigned int> index_set(
3793 internal::extract_field_component_indices<dim>(extractor));
3794 for (const unsigned int index : index_set)
3795 {
3796 Assert(
3797 this->registered_independent_variable_values[index] == false,
3798 ExcMessage(
3799 "Overlapping indices for independent variables. "
3800 "One or more indices associated with the field that "
3801 "is being registered as an independent variable have "
3802 "already been associated with another field. This suggests "
3803 "that the component offsets used to construct their counterpart "
3804 "extractors are incompatible with one another. Make sure that "
3805 "the first component for each extractor properly takes into "
3806 "account the dimensionality of the preceding fields."));
3807 }
3808# endif
3809 set_independent_variable(value, extractor);
3810 }
3811
3812
3813
3814 template <int dim,
3815 enum AD::NumberTypes ADNumberTypeCode,
3816 typename ScalarType>
3817 template <typename ValueType, typename ExtractorType>
3818 void
3820 set_independent_variable(const ValueType & value,
3821 const ExtractorType &extractor)
3822 {
3823 const std::vector<unsigned int> index_set(
3824 internal::extract_field_component_indices<dim>(extractor));
3825 for (unsigned int i = 0; i < index_set.size(); ++i)
3826 {
3827 set_sensitivity_value(
3828 index_set[i],
3829 internal::Extractor<dim, ExtractorType>::symmetric_component(i),
3830 internal::get_tensor_entry(value, i));
3831 }
3832 }
3833
3834
3835
3836 template <int dim,
3837 enum AD::NumberTypes ADNumberTypeCode,
3838 typename ScalarType>
3839 template <typename ExtractorType>
3840 typename internal::Extractor<dim, ExtractorType>::template tensor_type<
3841 typename PointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::
3842 ad_type>
3844 get_sensitive_variables(const ExtractorType &extractor) const
3845 {
3846 if (ADNumberTraits<ad_type>::is_taped == true)
3847 {
3848 Assert(this->active_tape_index() !=
3850 ExcMessage("Invalid tape index"));
3851 }
3852
3853 // If necessary, finalize the internally stored vector of
3854 // AD numbers that represents the independent variables
3855 this->finalize_sensitive_independent_variables();
3856 Assert(this->independent_variables.size() ==
3857 this->n_independent_variables(),
3858 ExcDimensionMismatch(this->independent_variables.size(),
3859 this->n_independent_variables()));
3860
3861 const std::vector<unsigned int> index_set(
3862 internal::extract_field_component_indices<dim>(extractor));
3863 typename internal::Extractor<dim,
3864 ExtractorType>::template tensor_type<ad_type>
3865 out;
3866
3867 for (unsigned int i = 0; i < index_set.size(); ++i)
3868 {
3869 const unsigned int index = index_set[i];
3870 Assert(index < this->n_independent_variables(), ExcInternalError());
3871 Assert(this->registered_independent_variable_values[index] == true,
3874 this->independent_variables[index];
3875 }
3876
3877 return out;
3878 }
3879
3880
3881
3882 /* ----------------- ScalarFunction ----------------- */
3883
3884
3885
3886 template <int dim,
3887 enum AD::NumberTypes ADNumberTypeCode,
3888 typename ScalarType>
3889 template <typename ExtractorType_Row>
3890 typename internal::ScalarFieldGradient<
3891 dim,
3893 ExtractorType_Row>::type
3896 const ExtractorType_Row & extractor_row)
3897 {
3898 // NOTE: The order of components must be consistently defined throughout
3899 // this class.
3900 typename internal::
3902
3903 // Get indexsets for the subblock from which we wish to extract the
3904 // gradient values
3905 const std::vector<unsigned int> row_index_set(
3906 internal::extract_field_component_indices<dim>(extractor_row));
3907 Assert(out.n_independent_components == row_index_set.size(),
3908 ExcMessage("Not all tensor components have been extracted!"));
3909 for (unsigned int r = 0; r < row_index_set.size(); ++r)
3910 internal::set_tensor_entry(out, r, gradient[row_index_set[r]]);
3911
3912 return out;
3913 }
3914
3915
3916
3917 template <int dim,
3918 enum AD::NumberTypes ADNumberTypeCode,
3919 typename ScalarType>
3920 template <typename ExtractorType_Row, typename ExtractorType_Col>
3921 typename internal::ScalarFieldHessian<
3922 dim,
3924 ExtractorType_Row,
3925 ExtractorType_Col>::type
3928 const ExtractorType_Row & extractor_row,
3929 const ExtractorType_Col & extractor_col)
3930 {
3931 using InternalHessian = internal::ScalarFieldHessian<dim,
3932 scalar_type,
3933 ExtractorType_Row,
3934 ExtractorType_Col>;
3935 using InternalExtractorRow = internal::Extractor<dim, ExtractorType_Row>;
3936 using InternalExtractorCol = internal::Extractor<dim, ExtractorType_Col>;
3937 using HessianType = typename InternalHessian::type;
3938
3939 // NOTE: The order of components must be consistently defined throughout
3940 // this class.
3941 HessianType out;
3942
3943 // Get indexsets for the subblocks from which we wish to extract the
3944 // Hessian values
3945 // NOTE: Here we have to do some clever accounting when the
3946 // one extractor is a symmetric Tensor and the other is not, e.g.
3947 // <SymmTensor,Vector>. In this scenario the return type is a
3948 // non-symmetric Tensor<3,dim> but we have to fetch information from a
3949 // SymmTensor row/column that has too few entries to fill the output
3950 // tensor. So we must duplicate the relevant entries in the row/column
3951 // indexset to fetch off-diagonal components that are otherwise
3952 // non-existent in a SymmTensor.
3953 const std::vector<unsigned int> row_index_set(
3954 internal::extract_field_component_indices<dim>(
3955 extractor_row, false /*ignore_symmetries*/));
3956 const std::vector<unsigned int> col_index_set(
3957 internal::extract_field_component_indices<dim>(
3958 extractor_col, false /*ignore_symmetries*/));
3959
3960 for (unsigned int index = 0;
3961 index < HessianType::n_independent_components;
3962 ++index)
3963 {
3964 const TableIndices<HessianType::rank> ti_out =
3965 HessianType::unrolled_to_component_indices(index);
3966 const unsigned int r =
3967 InternalExtractorRow::local_component(ti_out, 0);
3968 const unsigned int c =
3969 InternalExtractorCol::local_component(ti_out,
3970 InternalExtractorRow::rank);
3971
3973 out, index, hessian[row_index_set[r]][col_index_set[c]]);
3974 }
3975
3976 return out;
3977 }
3978
3979
3980
3981 /* ----------------- VectorFunction ----------------- */
3982
3983
3984
3985 template <int dim,
3986 enum AD::NumberTypes ADNumberTypeCode,
3987 typename ScalarType>
3988 template <typename ValueType, typename ExtractorType>
3989 void
3991 register_dependent_variable(const ValueType & funcs,
3992 const ExtractorType &extractor)
3993 {
3994 const std::vector<unsigned int> index_set(
3995 internal::extract_field_component_indices<dim>(extractor));
3996 for (unsigned int i = 0; i < index_set.size(); ++i)
3997 {
3998 Assert(this->registered_marked_dependent_variables[index_set[i]] ==
3999 false,
4000 ExcMessage("Overlapping indices for dependent variables."));
4002 index_set[i], internal::get_tensor_entry(funcs, i));
4003 }
4004 }
4005
4006
4007
4008 template <int dim,
4009 enum AD::NumberTypes ADNumberTypeCode,
4010 typename ScalarType>
4011 template <typename ExtractorType_Row>
4013 dim,
4015 ExtractorType_Row>::type
4018 const ExtractorType_Row & extractor_row)
4019 {
4020 // NOTE: The order of components must be consistently defined throughout
4021 // this class.
4022 typename internal::VectorFieldValue<dim, scalar_type, ExtractorType_Row>::
4023 type out;
4024
4025 // Get indexsets for the subblock from which we wish to extract the
4026 // gradient values
4027 const std::vector<unsigned int> row_index_set(
4028 internal::extract_field_component_indices<dim>(extractor_row));
4029 Assert(out.n_independent_components == row_index_set.size(),
4030 ExcMessage("Not all tensor components have been extracted!"));
4031 for (unsigned int r = 0; r < row_index_set.size(); ++r)
4032 internal::set_tensor_entry(out, r, values[row_index_set[r]]);
4033
4034 return out;
4035 }
4036
4037
4038
4039 template <int dim,
4040 enum AD::NumberTypes ADNumberTypeCode,
4041 typename ScalarType>
4042 template <typename ExtractorType_Row, typename ExtractorType_Col>
4044 dim,
4046 ExtractorType_Row,
4047 ExtractorType_Col>::type
4050 const ExtractorType_Row & extractor_row,
4051 const ExtractorType_Col & extractor_col)
4052 {
4053 using InternalJacobian = internal::VectorFieldJacobian<dim,
4054 scalar_type,
4055 ExtractorType_Row,
4056 ExtractorType_Col>;
4057 using InternalExtractorRow = internal::Extractor<dim, ExtractorType_Row>;
4058 using InternalExtractorCol = internal::Extractor<dim, ExtractorType_Col>;
4059 using JacobianType = typename InternalJacobian::type;
4060
4061 // NOTE: The order of components must be consistently defined throughout
4062 // this class.
4063 JacobianType out;
4064
4065 // Get indexsets for the subblocks from which we wish to extract the
4066 // Hessian values.
4067 // NOTE: Here we have to do some clever accounting when the
4068 // one extractor is a symmetric Tensor and the other is not, e.g.
4069 // <SymmTensor,Vector>. In this scenario the return type is a
4070 // non-symmetric Tensor<3,dim> but we have to fetch information from a
4071 // SymmTensor row/column that has too few entries to fill the output
4072 // tensor. So we must duplicate the relevant entries in the row/column
4073 // indexset to fetch off-diagonal components that are otherwise
4074 // non-existent in a SymmTensor.
4075 const std::vector<unsigned int> row_index_set(
4076 internal::extract_field_component_indices<dim>(
4077 extractor_row, false /*ignore_symmetries*/));
4078 const std::vector<unsigned int> col_index_set(
4079 internal::extract_field_component_indices<dim>(
4080 extractor_col, false /*ignore_symmetries*/));
4081
4082 for (unsigned int index = 0;
4083 index < JacobianType::n_independent_components;
4084 ++index)
4085 {
4087 JacobianType::unrolled_to_component_indices(index);
4088 const unsigned int r =
4089 InternalExtractorRow::local_component(ti_out, 0);
4090 const unsigned int c =
4091 InternalExtractorCol::local_component(ti_out,
4092 InternalExtractorRow::rank);
4093
4095 out, index, jacobian[row_index_set[r]][col_index_set[c]]);
4096 }
4097
4098 return out;
4099 }
4100
4101
4102 } // namespace AD
4103} // namespace Differentiation
4104
4105
4106# endif // DOXYGEN
4107
4108
4110
4111#endif // defined(DEAL_II_WITH_ADOLC) || defined(DEAL_II_TRILINOS_WITH_SACADO)
4112
4113#endif // dealii_differentiation_ad_ad_helpers_h
void register_dof_values(const VectorType &values, const std::vector<::types::global_dof_index > &local_dof_indices)
const std::vector< ad_type > & get_sensitive_dof_values() const
Definition: ad_helpers.cc:744
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
Definition: ad_helpers.h:851
void set_dof_values(const std::vector< scalar_type > &dof_values)
Definition: ad_helpers.cc:769
virtual void compute_residual(Vector< scalar_type > &residual) const =0
void register_dof_values(const std::vector< scalar_type > &dof_values)
Definition: ad_helpers.cc:721
virtual void compute_linearization(FullMatrix< scalar_type > &linearization) const =0
CellLevelBase(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
Definition: ad_helpers.cc:710
void set_dof_values(const VectorType &values, const std::vector<::types::global_dof_index > &local_dof_indices)
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
Definition: ad_helpers.h:844
void register_energy_functional(const ad_type &energy)
Definition: ad_helpers.cc:802
EnergyFunctional(const unsigned int n_independent_variables)
Definition: ad_helpers.cc:793
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
Definition: ad_helpers.h:1233
void compute_residual(Vector< scalar_type > &residual) const override
Definition: ad_helpers.cc:868
virtual void compute_linearization(FullMatrix< scalar_type > &linearization) const override
Definition: ad_helpers.cc:933
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
Definition: ad_helpers.h:1226
unsigned int n_registered_dependent_variables() const
Definition: ad_helpers.cc:250
TapedDrivers< ad_type, scalar_type > taped_driver
Definition: ad_helpers.h:592
void stop_recording_operations(const bool write_tapes_to_file=false)
Definition: ad_helpers.cc:643
std::vector< bool > registered_marked_independent_variables
Definition: ad_helpers.h:658
typename AD::NumberTraits< ScalarType, ADNumberTypeCode >::scalar_type scalar_type
Definition: ad_helpers.h:178
std::vector< bool > registered_independent_variable_values
Definition: ad_helpers.h:652
Types< ad_type >::tape_index active_tape_index() const
Definition: ad_helpers.cc:284
void print_tape_stats(const typename Types< ad_type >::tape_index tape_index, std::ostream &stream) const
Definition: ad_helpers.cc:373
std::size_t n_independent_variables() const
Definition: ad_helpers.cc:241
bool start_recording_operations(const typename Types< ad_type >::tape_index tape_index, const bool overwrite_tape=false, const bool keep_independent_values=true)
Definition: ad_helpers.cc:578
std::vector< scalar_type > independent_variable_values
Definition: ad_helpers.h:636
void mark_independent_variable(const unsigned int index, ad_type &out) const
Definition: ad_helpers.cc:142
bool is_registered_tape(const typename Types< ad_type >::tape_index tape_index) const
Definition: ad_helpers.cc:296
std::vector< bool > registered_marked_dependent_variables
Definition: ad_helpers.h:755
void reset_registered_dependent_variables(const bool flag=false)
Definition: ad_helpers.cc:94
void activate_tape(const typename Types< ad_type >::tape_index tape_index, const bool read_mode)
Definition: ad_helpers.cc:522
void set_tape_buffer_sizes(const typename Types< ad_type >::tape_buffer_sizes obufsize=64 *1024 *1024, const typename Types< ad_type >::tape_buffer_sizes lbufsize=64 *1024 *1024, const typename Types< ad_type >::tape_buffer_sizes vbufsize=64 *1024 *1024, const typename Types< ad_type >::tape_buffer_sizes tbufsize=64 *1024 *1024)
Definition: ad_helpers.cc:559
unsigned int n_registered_independent_variables() const
Definition: ad_helpers.cc:230
void print(std::ostream &stream) const
Definition: ad_helpers.cc:309
void set_sensitivity_value(const unsigned int index, const scalar_type &value)
Definition: ad_helpers.cc:105
bool active_tape_requires_retaping() const
Definition: ad_helpers.cc:497
std::vector< ad_type > independent_variables
Definition: ad_helpers.h:646
HelperBase(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
Definition: ad_helpers.cc:38
virtual void reset(const unsigned int n_independent_variables=::numbers::invalid_unsigned_int, const unsigned int n_dependent_variables=::numbers::invalid_unsigned_int, const bool clear_registered_tapes=true)
Definition: ad_helpers.cc:390
TapelessDrivers< ad_type, scalar_type > tapeless_driver
Definition: ad_helpers.h:600
bool recorded_tape_requires_retaping(const typename Types< ad_type >::tape_index tape_index) const
Definition: ad_helpers.cc:484
std::vector< ad_type > dependent_variables
Definition: ad_helpers.h:750
std::size_t n_dependent_variables() const
Definition: ad_helpers.cc:262
void print_values(std::ostream &stream) const
Definition: ad_helpers.cc:359
void register_dependent_variable(const unsigned int index, const ad_type &func)
Definition: ad_helpers.cc:679
void finalize_sensitive_independent_variables() const
Definition: ad_helpers.cc:181
static void configure_tapeless_mode(const unsigned int n_independent_variables, const bool ensure_persistent_setting=true)
Definition: ad_helpers.cc:448
void initialize_non_sensitive_independent_variable(const unsigned int index, ad_type &out) const
Definition: ad_helpers.cc:205
typename AD::NumberTraits< ScalarType, ADNumberTypeCode >::ad_type ad_type
Definition: ad_helpers.h:185
void activate_recorded_tape(const typename Types< ad_type >::tape_index tape_index)
Definition: ad_helpers.cc:474
void set_independent_variable(const ValueType &value, const ExtractorType &extractor)
bool is_symmetric_independent_variable(const unsigned int index) const
Definition: ad_helpers.cc:1192
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
Definition: ad_helpers.h:2663
PointLevelFunctionsBase(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
Definition: ad_helpers.cc:1155
static constexpr unsigned int dimension
Definition: ad_helpers.h:2649
void set_independent_variables(const std::vector< scalar_type > &values)
Definition: ad_helpers.cc:1296
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
Definition: ad_helpers.h:2656
std::vector< bool > symmetric_independent_variables
Definition: ad_helpers.h:2921
unsigned int n_symmetric_independent_variables() const
Definition: ad_helpers.cc:1206
void register_independent_variables(const std::vector< scalar_type > &values)
Definition: ad_helpers.cc:1220
void set_sensitivity_value(const unsigned int index, const bool symmetric_component, const scalar_type &value)
Definition: ad_helpers.cc:1274
void register_independent_variable(const ValueType &value, const ExtractorType &extractor)
internal::Extractor< dim, ExtractorType >::template tensor_type< ad_type > get_sensitive_variables(const ExtractorType &extractor) const
const std::vector< ad_type > & get_sensitive_variables() const
Definition: ad_helpers.cc:1246
virtual void reset(const unsigned int n_independent_variables=::numbers::invalid_unsigned_int, const unsigned int n_dependent_variables=::numbers::invalid_unsigned_int, const bool clear_registered_tapes=true) override
Definition: ad_helpers.cc:1168
ResidualLinearization(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
Definition: ad_helpers.cc:1006
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
Definition: ad_helpers.h:1538
virtual void compute_residual(Vector< scalar_type > &residual) const override
Definition: ad_helpers.cc:1032
virtual void compute_linearization(FullMatrix< scalar_type > &linearization) const override
Definition: ad_helpers.cc:1086
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
Definition: ad_helpers.h:1545
void register_residual_vector(const std::vector< ad_type > &residual)
Definition: ad_helpers.cc:1018
static internal::ScalarFieldHessian< dim, scalar_type, ExtractorType_Row, ExtractorType_Col >::type extract_hessian_component(const FullMatrix< scalar_type > &hessian, const ExtractorType_Row &extractor_row, const ExtractorType_Col &extractor_col)
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
Definition: ad_helpers.h:3102
void register_dependent_variable(const ad_type &func)
Definition: ad_helpers.cc:1335
static internal::ScalarFieldGradient< dim, scalar_type, ExtractorType_Row >::type extract_gradient_component(const Vector< scalar_type > &gradient, const ExtractorType_Row &extractor_row)
ScalarFunction(const unsigned int n_independent_variables)
Definition: ad_helpers.cc:1321
void compute_gradient(Vector< scalar_type > &gradient) const
Definition: ad_helpers.cc:1398
void compute_hessian(FullMatrix< scalar_type > &hessian) const
Definition: ad_helpers.cc:1472
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
Definition: ad_helpers.h:3095
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
Definition: ad_helpers.h:3493
void compute_values(Vector< scalar_type > &values) const
Definition: ad_helpers.cc:1677
void register_dependent_variables(const std::vector< ad_type > &funcs)
Definition: ad_helpers.cc:1661
void compute_jacobian(FullMatrix< scalar_type > &jacobian) const
Definition: ad_helpers.cc:1733
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
Definition: ad_helpers.h:3486
static internal::VectorFieldJacobian< dim, scalar_type, ExtractorType_Row, ExtractorType_Col >::type extract_jacobian_component(const FullMatrix< scalar_type > &jacobian, const ExtractorType_Row &extractor_row, const ExtractorType_Col &extractor_col)
static internal::VectorFieldValue< dim, scalar_type, ExtractorType_Row >::type extract_value_component(const Vector< scalar_type > &values, const ExtractorType_Row &extractor_row)
void register_dependent_variable(const ValueType &funcs, const ExtractorType &extractor)
VectorFunction(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
Definition: ad_helpers.cc:1646
static constexpr unsigned int component_to_unrolled_index(const TableIndices< rank_ > &indices)
Definition: tensor.h:503
static constexpr TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
Definition: vector.h:109
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
#define Assert(cond, exc)
Definition: exceptions.h:1473
#define AssertIndexRange(index, range)
Definition: exceptions.h:1732
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
ScalarFieldHessian< dim, NumberType, ExtractorType_Field, ExtractorType_Derivative > VectorFieldJacobian
Definition: ad_helpers.h:2345
NumberType get_tensor_entry(const TensorType< rank, dim, NumberType > &t, const unsigned int unrolled_index)
Definition: ad_helpers.h:2518
void set_tensor_entry(TensorType &t, const unsigned int unrolled_index, const NumberType &value)
Definition: ad_helpers.h:2437
std::vector< IndexType > extract_field_component_indices(const ExtractorType &extractor, const bool ignore_symmetries=true)
Definition: ad_helpers.h:2357
ScalarFieldGradient< dim, NumberType, ExtractorType_Field > VectorFieldValue
Definition: ad_helpers.h:2326
void set_dof_values(const DoFCellAccessor< dim, spacedim, lda > &cell, const Vector< number > &local_values, OutputVector &values, const bool perform_check)
static const unsigned int invalid_unsigned_int
Definition: types.h:201
static const Types< ADNumberType >::tape_index invalid_tape_index
Definition: ad_drivers.h:121
static unsigned int first_component(const FEValuesExtractors::Scalar &extractor)
Definition: ad_helpers.h:1704
static IndexType local_component(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:1732
static bool symmetric_component(const unsigned int unrolled_index)
Definition: ad_helpers.h:1717
static TableIndices< rank > table_index_view(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:2129
static unsigned int first_component(const FEValuesExtractors::SymmetricTensor< 2 > &extractor)
Definition: ad_helpers.h:2102
static IndexType local_component(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:2155
static TableIndices< rank > table_index_view(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:1919
static unsigned int first_component(const FEValuesExtractors::Tensor< 1 > &extractor)
Definition: ad_helpers.h:1894
static IndexType local_component(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:1943
static TableIndices< rank > table_index_view(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:2022
static unsigned int first_component(const FEValuesExtractors::Tensor< 2 > &extractor)
Definition: ad_helpers.h:1997
static IndexType local_component(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:2048
static bool symmetric_component(const unsigned int unrolled_index)
Definition: ad_helpers.h:1804
static IndexType local_component(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:1840
static unsigned int first_component(const FEValuesExtractors::Vector &extractor)
Definition: ad_helpers.h:1790
static TableIndices< rank > table_index_view(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:1816
typename Extractor< dim, ExtractorType >::template tensor_type< NumberType > type
Definition: ad_helpers.h:2182
typename HessianType< ExtractorType_Row, ExtractorType_Col >::template type< rank, dim, NumberType > type
Definition: ad_helpers.h:2316