Reference documentation for deal.II version 9.2.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
ad_helpers.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 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 
36 # include <deal.II/lac/full_matrix.h>
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 
47 namespace Differentiation
48 {
49  namespace AD
50  {
170  template <enum AD::NumberTypes ADNumberTypeCode,
171  typename ScalarType = double>
173  {
174  public:
179  using scalar_type =
181 
186  using ad_type =
188 
193 
208  HelperBase(const unsigned int n_independent_variables,
209  const unsigned int n_dependent_variables);
210 
214  virtual ~HelperBase() = default;
215 
217 
222 
227  std::size_t
228  n_independent_variables() const;
229 
234  std::size_t
235  n_dependent_variables() const;
236 
245  void
246  print(std::ostream &stream) const;
247 
254  void
255  print_values(std::ostream &stream) const;
256 
267  void
268  print_tape_stats(const typename Types<ad_type>::tape_index tape_index,
269  std::ostream & stream) const;
270 
272 
277 
296  static void
298  const bool ensure_persistent_setting = true);
299 
301 
306 
342  virtual void
343  reset(const unsigned int n_independent_variables =
345  const unsigned int n_dependent_variables =
347  const bool clear_registered_tapes = true);
348 
353  bool
354  is_recording() const;
355 
361  active_tape_index() const;
362 
367  bool
369  const typename Types<ad_type>::tape_index tape_index) const;
370 
396  void
398  const typename Types<ad_type>::tape_buffer_sizes obufsize = 64 * 1024 *
399  1024,
400  const typename Types<ad_type>::tape_buffer_sizes lbufsize = 64 * 1024 *
401  1024,
402  const typename Types<ad_type>::tape_buffer_sizes vbufsize = 64 * 1024 *
403  1024,
404  const typename Types<ad_type>::tape_buffer_sizes tbufsize = 64 * 1024 *
405  1024);
406 
452  bool
454  const typename Types<ad_type>::tape_index tape_index,
455  const bool overwrite_tape = false,
456  const bool keep_independent_values = true);
457 
469  void
470  stop_recording_operations(const bool write_tapes_to_file = false);
471 
481  void
483  const typename Types<ad_type>::tape_index tape_index);
484 
524  bool
526  const typename Types<ad_type>::tape_index tape_index) const;
527 
567  bool
569 
577  void
579 
581 
582  protected:
587 
595 
603 
620  void
621  activate_tape(const typename Types<ad_type>::tape_index tape_index,
622  const bool read_mode);
623 
625 
630 
638  mutable std::vector<scalar_type> independent_variable_values;
639 
648  mutable std::vector<ad_type> independent_variables;
649 
655 
660  mutable std::vector<bool> registered_marked_independent_variables;
661 
667  void
669 
676  void
677  set_sensitivity_value(const unsigned int index, const scalar_type &value);
678 
697  void
698  mark_independent_variable(const unsigned int index, ad_type &out) const;
699 
710  void
712 
722  void
723  initialize_non_sensitive_independent_variable(const unsigned int index,
724  ad_type &out) const;
725 
730  unsigned int
732 
734 
739 
752  std::vector<ad_type> dependent_variables;
753 
758 
765  void
766  reset_registered_dependent_variables(const bool flag = false);
767 
771  unsigned int
773 
785  void
786  register_dependent_variable(const unsigned int index,
787  const ad_type & func);
788 
790 
791  }; // class HelperBase
792 
793 
794 
838  template <enum AD::NumberTypes ADNumberTypeCode,
839  typename ScalarType = double>
840  class CellLevelBase : public HelperBase<ADNumberTypeCode, ScalarType>
841  {
842  public:
847  using scalar_type =
849 
854  using ad_type =
856 
861 
876  CellLevelBase(const unsigned int n_independent_variables,
877  const unsigned int n_dependent_variables);
878 
882  virtual ~CellLevelBase() = default;
883 
885 
890 
908  void
909  register_dof_values(const std::vector<scalar_type> &dof_values);
910 
929  template <typename VectorType>
930  void
932  const VectorType & values,
933  const std::vector<::types::global_dof_index> &local_dof_indices);
934 
955  const std::vector<ad_type> &
956  get_sensitive_dof_values() const;
957 
959 
964 
984  void
985  set_dof_values(const std::vector<scalar_type> &dof_values);
986 
1004  template <typename VectorType>
1005  void
1007  const VectorType & values,
1008  const std::vector<::types::global_dof_index> &local_dof_indices);
1009 
1011 
1016 
1031  virtual void
1032  compute_residual(Vector<scalar_type> &residual) const = 0;
1033 
1051  virtual void
1052  compute_linearization(FullMatrix<scalar_type> &linearization) const = 0;
1053 
1055 
1056  }; // class CellLevelBase
1057 
1058 
1059 
1222  template <enum AD::NumberTypes ADNumberTypeCode,
1223  typename ScalarType = double>
1224  class EnergyFunctional : public CellLevelBase<ADNumberTypeCode, ScalarType>
1225  {
1226  public:
1231  using scalar_type =
1233 
1238  using ad_type =
1240 
1245 
1261  EnergyFunctional(const unsigned int n_independent_variables);
1262 
1266  virtual ~EnergyFunctional() = default;
1267 
1269 
1274 
1288  void
1289  register_energy_functional(const ad_type &energy);
1290 
1304  scalar_type
1305  compute_energy() const;
1306 
1326  void
1327  compute_residual(Vector<scalar_type> &residual) const override;
1328 
1349  virtual void
1351  FullMatrix<scalar_type> &linearization) const override;
1352 
1354 
1355  }; // class EnergyFunctional
1356 
1357 
1358 
1537  template <enum AD::NumberTypes ADNumberTypeCode,
1538  typename ScalarType = double>
1540  : public CellLevelBase<ADNumberTypeCode, ScalarType>
1541  {
1542  public:
1547  using scalar_type =
1549 
1554  using ad_type =
1556 
1561 
1577  const unsigned int n_dependent_variables);
1578 
1582  virtual ~ResidualLinearization() = default;
1583 
1585 
1590 
1604  void
1605  register_residual_vector(const std::vector<ad_type> &residual);
1606 
1623  virtual void
1624  compute_residual(Vector<scalar_type> &residual) const override;
1625 
1643  virtual void
1645  FullMatrix<scalar_type> &linearization) const override;
1646 
1648 
1649  }; // class ResidualLinearization
1650 
1651 
1652 
1653  namespace internal
1654  {
1659  template <int dim, typename ExtractorType>
1660  struct Extractor;
1661 
1662 
1668  template <int dim>
1669  struct Extractor<dim, FEValuesExtractors::Scalar>
1670  {
1674  static const unsigned int n_components = 1;
1675 
1679  static const unsigned int rank = 0;
1680 
1684  template <typename NumberType>
1686 
1687  static_assert(
1689  "The number of components doesn't match that of the corresponding tensor type.");
1690  static_assert(
1691  rank == tensor_type<double>::rank,
1692  "The rank doesn't match that of the corresponding tensor type.");
1693 
1697  // Note: FEValuesViews::Scalar::tensor_type is double, so we can't
1698  // use it (FEValuesViews) in this context.
1699  // In fact, sadly, all FEValuesViews objects expect doubles as value
1700  // types.
1701  template <typename NumberType>
1703 
1707  template <typename NumberType>
1709 
1713  static inline unsigned int
1715  {
1716  return extractor.component;
1717  }
1718 
1726  static bool
1727  symmetric_component(const unsigned int unrolled_index)
1728  {
1729  (void)unrolled_index;
1730  return false;
1731  }
1732 
1740  template <typename IndexType = unsigned int, int rank_in>
1741  static IndexType
1743  const unsigned int column_offset)
1744  {
1745  Assert(column_offset <= rank_in, ExcInternalError());
1746  (void)table_indices;
1747  (void)column_offset;
1748  return 0;
1749  }
1750  };
1751 
1752 
1758  template <int dim>
1760  {
1764  static const unsigned int n_components = dim;
1765 
1769  static const unsigned int rank = 1;
1770 
1774  template <typename NumberType>
1776 
1777  static_assert(
1779  "The number of components doesn't match that of the corresponding tensor type.");
1780  static_assert(
1781  rank == tensor_type<double>::rank,
1782  "The rank doesn't match that of the corresponding tensor type.");
1783 
1787  template <typename NumberType>
1789 
1793  template <typename NumberType>
1795 
1799  static inline unsigned int
1801  {
1802  return extractor.first_vector_component;
1803  }
1804 
1813  static bool
1814  symmetric_component(const unsigned int unrolled_index)
1815  {
1816  (void)unrolled_index;
1817  return false;
1818  }
1819 
1824  template <int rank_in>
1825  static TableIndices<rank>
1827  const unsigned int column_offset)
1828  {
1829  Assert(0 + column_offset < rank_in, ExcInternalError());
1830  return TableIndices<rank>(table_indices[column_offset]);
1831  }
1832 
1848  template <typename IndexType = unsigned int, int rank_in>
1849  static IndexType
1851  const unsigned int column_offset)
1852  {
1853  static_assert(
1854  rank_in >= rank,
1855  "Cannot extract more table indices than the input table has!");
1856  using TensorType = tensor_type<double>;
1857  return TensorType::component_to_unrolled_index(
1858  table_index_view(table_indices, column_offset));
1859  }
1860  };
1861 
1862 
1868  template <int dim>
1870  {
1874  static const unsigned int n_components =
1876 
1880  static const unsigned int rank = 1;
1881 
1885  template <typename NumberType>
1887 
1891  template <typename NumberType>
1893 
1897  template <typename NumberType>
1899 
1903  static inline unsigned int
1905  {
1906  return extractor.first_tensor_component;
1907  }
1908 
1916  static bool
1917  symmetric_component(const unsigned int unrolled_index)
1918  {
1919  (void)unrolled_index;
1920  return false;
1921  }
1922 
1927  template <int rank_in>
1928  static TableIndices<rank>
1930  const unsigned int column_offset)
1931  {
1932  Assert(column_offset < rank_in, ExcInternalError());
1933  return TableIndices<rank>(table_indices[column_offset]);
1934  }
1935 
1951  template <typename IndexType = unsigned int, int rank_in>
1952  static IndexType
1954  const unsigned int column_offset)
1955  {
1956  static_assert(
1957  rank_in >= rank,
1958  "Cannot extract more table indices than the input table has!");
1959  using TensorType = tensor_type<double>;
1960  return TensorType::component_to_unrolled_index(
1961  table_index_view(table_indices, column_offset));
1962  }
1963  };
1964 
1965 
1971  template <int dim>
1973  {
1977  static const unsigned int n_components =
1979 
1983  static const unsigned int rank = Tensor<2, dim>::rank;
1984 
1988  template <typename NumberType>
1990 
1994  template <typename NumberType>
1996 
2000  template <typename NumberType>
2002 
2006  static inline unsigned int
2008  {
2009  return extractor.first_tensor_component;
2010  }
2011 
2019  static bool
2020  symmetric_component(const unsigned int unrolled_index)
2021  {
2022  (void)unrolled_index;
2023  return false;
2024  }
2025 
2030  template <int rank_in>
2031  static TableIndices<rank>
2033  const unsigned int column_offset)
2034  {
2035  Assert(column_offset < rank_in, ExcInternalError());
2036  Assert(column_offset + 1 < rank_in, ExcInternalError());
2037  return TableIndices<rank>(table_indices[column_offset],
2038  table_indices[column_offset + 1]);
2039  }
2040 
2056  template <typename IndexType = unsigned int, int rank_in>
2057  static IndexType
2059  const unsigned int column_offset)
2060  {
2061  static_assert(
2062  rank_in >= rank,
2063  "Cannot extract more table indices than the input table has!");
2064  using TensorType = tensor_type<double>;
2065  return TensorType::component_to_unrolled_index(
2066  table_index_view(table_indices, column_offset));
2067  }
2068  };
2069 
2070 
2076  template <int dim>
2078  {
2082  static const unsigned int n_components =
2084 
2088  static const unsigned int rank = SymmetricTensor<2, dim>::rank;
2089 
2093  template <typename NumberType>
2095 
2099  template <typename NumberType>
2101 
2105  template <typename NumberType>
2107 
2111  static inline unsigned int
2113  {
2114  return extractor.first_tensor_component;
2115  }
2116 
2125  static bool
2126  symmetric_component(const unsigned int unrolled_index)
2127  {
2128  const TableIndices<2> table_indices =
2130  return table_indices[0] != table_indices[1];
2131  }
2132 
2137  template <int rank_in>
2138  static TableIndices<rank>
2140  const unsigned int column_offset)
2141  {
2142  Assert(column_offset < rank_in, ExcInternalError());
2143  Assert(column_offset + 1 < rank_in, ExcInternalError());
2144  return TableIndices<rank>(table_indices[column_offset],
2145  table_indices[column_offset + 1]);
2146  }
2147 
2163  template <typename IndexType = unsigned int, int rank_in>
2164  static IndexType
2166  const unsigned int column_offset)
2167  {
2168  static_assert(
2169  rank_in >= rank,
2170  "Cannot extract more table indices than the input table has!");
2171  using TensorType = tensor_type<double>;
2172  return TensorType::component_to_unrolled_index(
2173  table_index_view(table_indices, column_offset));
2174  }
2175  };
2176 
2177 
2183  template <int dim, typename NumberType, typename ExtractorType>
2185  {
2190  using type =
2191  typename Extractor<dim,
2192  ExtractorType>::template tensor_type<NumberType>;
2193  };
2194 
2195 
2205  template <typename ExtractorType_Row, typename ExtractorType_Col>
2207  {
2219  template <int rank, int dim, typename NumberType>
2221  };
2222 
2223 
2232  template <>
2235  {
2243  template <int /*rank*/, int dim, typename NumberType>
2244  using type = SymmetricTensor<2 /*rank*/, dim, NumberType>;
2245  };
2246 
2247 
2256  template <>
2259  {
2268  template <int /*rank*/, int dim, typename NumberType>
2269  using type = SymmetricTensor<2 /*rank*/, dim, NumberType>;
2270  };
2271 
2272 
2280  template <>
2283  {
2292  template <int /*rank*/, int dim, typename NumberType>
2293  using type = SymmetricTensor<4 /*rank*/, dim, NumberType>;
2294  };
2295 
2296 
2306  template <int dim,
2307  typename NumberType,
2308  typename ExtractorType_Row,
2309  typename ExtractorType_Col>
2311  {
2317 
2324  using type =
2327  };
2328 
2329 
2334  template <int dim, typename NumberType, typename ExtractorType_Field>
2335  using VectorFieldValue =
2337 
2338 
2348  template <int dim,
2349  typename NumberType,
2350  typename ExtractorType_Field,
2351  typename ExtractorType_Derivative>
2353  NumberType,
2354  ExtractorType_Field,
2355  ExtractorType_Derivative>;
2356 
2357 
2363  template <int dim,
2364  typename IndexType = unsigned int,
2365  typename ExtractorType>
2366  std::vector<IndexType>
2367  extract_field_component_indices(const ExtractorType &extractor,
2368  const bool ignore_symmetries = true)
2369  {
2370  (void)ignore_symmetries;
2371  const IndexType n_components =
2373  const IndexType comp_first =
2375  std::vector<IndexType> indices(n_components);
2376  std::iota(indices.begin(), indices.end(), comp_first);
2377  return indices;
2378  }
2379 
2380 
2389  template <int dim, typename IndexType = unsigned int>
2390  std::vector<IndexType>
2392  const FEValuesExtractors::SymmetricTensor<2> &extractor_symm_tensor,
2393  const bool ignore_symmetries = true)
2394  {
2395  using ExtractorType = FEValuesExtractors::SymmetricTensor<2>;
2396  const IndexType n_components =
2398  if (ignore_symmetries == true)
2399  {
2400  const IndexType comp_first =
2402  extractor_symm_tensor);
2403  std::vector<IndexType> indices(n_components);
2404  std::iota(indices.begin(), indices.end(), comp_first);
2405  return indices;
2406  }
2407  else
2408  {
2409  // First get all of the indices of the non-symmetric tensor
2410  const FEValuesExtractors::Tensor<2> extractor_tensor(
2411  extractor_symm_tensor.first_tensor_component);
2412  std::vector<IndexType> indices =
2413  extract_field_component_indices<dim>(extractor_tensor, true);
2414 
2415  // Then we overwrite any illegal entries with the equivalent indices
2416  // from the symmetric tensor
2417  for (unsigned int i = 0; i < indices.size(); ++i)
2418  {
2419  if (indices[i] >= n_components)
2420  {
2421  const TableIndices<2> ti_tensor =
2423  const IndexType sti_new_index =
2425  ti_tensor);
2426  indices[i] = sti_new_index;
2427  }
2428  }
2429 
2430  return indices;
2431  }
2432  }
2433 
2434 
2439  template <typename TensorType, typename NumberType>
2440  inline void
2441  set_tensor_entry(TensorType & t,
2442  const unsigned int unrolled_index,
2443  const NumberType & value)
2444  {
2445  // Where possible, set values using TableIndices
2446  AssertIndexRange(unrolled_index, t.n_independent_components);
2447  t[TensorType::unrolled_to_component_indices(unrolled_index)] = value;
2448  }
2449 
2450 
2455  template <int dim, typename NumberType>
2457  const unsigned int unrolled_index,
2458  const NumberType & value)
2459  {
2460  AssertIndexRange(unrolled_index, 1);
2461  (void)unrolled_index;
2462  t = value;
2463  }
2464 
2465 
2471  template <typename NumberType>
2472  inline void
2474  const unsigned int unrolled_index,
2475  const NumberType & value)
2476  {
2477  AssertIndexRange(unrolled_index, 1);
2478  (void)unrolled_index;
2479  t = value;
2480  }
2481 
2482 
2488  template <int dim, typename NumberType>
2490  const unsigned int unrolled_index_row,
2491  const unsigned int unrolled_index_col,
2492  const NumberType & value)
2493  {
2494  // Fourth order symmetric tensors require a specialized interface
2495  // to extract values.
2496  using SubTensorType = SymmetricTensor<2, dim, NumberType>;
2497  AssertIndexRange(unrolled_index_row,
2498  SubTensorType::n_independent_components);
2499  AssertIndexRange(unrolled_index_col,
2500  SubTensorType::n_independent_components);
2501  const TableIndices<2> indices_row =
2502  SubTensorType::unrolled_to_component_indices(unrolled_index_row);
2503  const TableIndices<2> indices_col =
2504  SubTensorType::unrolled_to_component_indices(unrolled_index_col);
2505  t[indices_row[0]][indices_row[1]][indices_col[0]][indices_col[1]] =
2506  value;
2507  }
2508 
2509 
2514  template <int rank,
2515  int dim,
2516  typename NumberType,
2517  template <int, int, typename> class TensorType>
2518  inline NumberType
2519  get_tensor_entry(const TensorType<rank, dim, NumberType> &t,
2520  const unsigned int unrolled_index)
2521  {
2522  // Where possible, get values using TableIndices
2523  AssertIndexRange(unrolled_index, t.n_independent_components);
2524  return t[TensorType<rank, dim, NumberType>::
2525  unrolled_to_component_indices(unrolled_index)];
2526  }
2527 
2528 
2533  template <int dim,
2534  typename NumberType,
2535  template <int, int, typename> 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> class TensorType>
2569  inline NumberType &
2570  get_tensor_entry(TensorType<rank, dim, NumberType> &t,
2571  const unsigned int unrolled_index)
2572  {
2573  // Where possible, get values using TableIndices
2574  AssertIndexRange(unrolled_index, t.n_independent_components);
2575  return t[TensorType<rank, dim, NumberType>::
2576  unrolled_to_component_indices(unrolled_index)];
2577  }
2578 
2579 
2584  template <int dim,
2585  typename NumberType,
2586  template <int, int, typename> class TensorType>
2587  NumberType &get_tensor_entry(TensorType<0, dim, NumberType> &t,
2588  const unsigned int index)
2589  {
2590  AssertIndexRange(index, 1);
2591  (void)index;
2592  return t;
2593  }
2594 
2595 
2601  template <typename NumberType>
2602  inline NumberType &
2603  get_tensor_entry(NumberType &t, const unsigned int index)
2604  {
2605  AssertIndexRange(index, 1);
2606  (void)index;
2607  return t;
2608  }
2609 
2610  } // namespace internal
2611 
2612 
2613 
2637  template <int dim,
2638  enum AD::NumberTypes ADNumberTypeCode,
2639  typename ScalarType = double>
2641  : public HelperBase<ADNumberTypeCode, ScalarType>
2642  {
2643  public:
2648  static const unsigned int dimension = dim;
2649 
2654  using scalar_type =
2656 
2661  using ad_type =
2663 
2668 
2684  const unsigned int n_dependent_variables);
2685 
2689  virtual ~PointLevelFunctionsBase() = default;
2690 
2692 
2697 
2701  virtual void
2702  reset(const unsigned int n_independent_variables =
2704  const unsigned int n_dependent_variables =
2706  const bool clear_registered_tapes = true) override;
2707 
2723  void
2724  register_independent_variables(const std::vector<scalar_type> &values);
2725 
2754  template <typename ValueType, typename ExtractorType>
2755  void
2756  register_independent_variable(const ValueType & value,
2757  const ExtractorType &extractor);
2758 
2779  const std::vector<ad_type> &
2780  get_sensitive_variables() const;
2781 
2782  /*
2783  * Extract a subset of the independent variables as represented by
2784  * auto-differentiable numbers. These are the independent
2785  * variables @f$\mathbf{A} \subset \mathbf{X}@f$ at which the dependent values
2786  * are evaluated and differentiated.
2787  *
2788  * This function indicates to the AD library that implements the
2789  * auto-differentiable number type that operations performed on these
2790  * numbers are to be tracked so they are considered "sensitive"
2791  * variables. This is, therefore, the set of variables with which one
2792  * would then perform computations, and based on which one can then
2793  * extract both the value of the function and its derivatives with the
2794  * member functions below. The values of the components of the returned
2795  * object are initialized to the values set with
2796  * register_independent_variable().
2797  *
2798  * @param[in] extractor An extractor associated with the input field
2799  * variables. This effectively defines which components of the global set
2800  * of independent variables this field is associated with.
2801  * @return An object of auto-differentiable type numbers. The return type is
2802  * based on the input extractor, and will be either a scalar,
2803  * Tensor<1,dim>, Tensor<2,dim>, or SymmetricTensor<2,dim>.
2804  *
2805  * @note For taped AD numbers, this operation is only valid in recording mode.
2806  */
2807  template <typename ExtractorType>
2808  typename internal::Extractor<dim,
2809  ExtractorType>::template tensor_type<ad_type>
2810  get_sensitive_variables(const ExtractorType &extractor) const;
2811 
2813 
2818 
2836  void
2837  set_independent_variables(const std::vector<scalar_type> &values);
2838 
2865  template <typename ValueType, typename ExtractorType>
2866  void
2867  set_independent_variable(const ValueType & value,
2868  const ExtractorType &extractor);
2869 
2871 
2872  protected:
2877 
2888  void
2889  set_sensitivity_value(const unsigned int index,
2890  const bool symmetric_component,
2891  const scalar_type &value);
2892 
2898  bool
2899  is_symmetric_independent_variable(const unsigned int index) const;
2900 
2905  unsigned int
2907 
2909 
2910  private:
2915 
2921 
2923 
2924  }; // class PointLevelFunctionsBase
2925 
2926 
2927 
3084  template <int dim,
3085  enum AD::NumberTypes ADNumberTypeCode,
3086  typename ScalarType = double>
3088  : public PointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>
3089  {
3090  public:
3095  using scalar_type =
3097 
3102  using ad_type =
3104 
3109 
3119  ScalarFunction(const unsigned int n_independent_variables);
3120 
3124  virtual ~ScalarFunction() = default;
3125 
3127 
3132 
3144  void
3145  register_dependent_variable(const ad_type &func);
3146 
3154  scalar_type
3155  compute_value() const;
3156 
3169  void
3170  compute_gradient(Vector<scalar_type> &gradient) const;
3171 
3186  void
3187  compute_hessian(FullMatrix<scalar_type> &hessian) const;
3188 
3223  template <typename ExtractorType_Row>
3224  static typename internal::
3226  extract_gradient_component(const Vector<scalar_type> &gradient,
3227  const ExtractorType_Row & extractor_row);
3228 
3267  template <typename ExtractorType_Row, typename ExtractorType_Col>
3268  static typename internal::ScalarFieldHessian<dim,
3269  scalar_type,
3270  ExtractorType_Row,
3271  ExtractorType_Col>::type
3273  const ExtractorType_Row & extractor_row,
3274  const ExtractorType_Col & extractor_col);
3275 
3291  const FullMatrix<scalar_type> & hessian,
3292  const FEValuesExtractors::Scalar &extractor_row,
3293  const FEValuesExtractors::Scalar &extractor_col);
3294 
3308  const FullMatrix<scalar_type> & hessian,
3309  const FEValuesExtractors::SymmetricTensor<2> &extractor_row,
3310  const FEValuesExtractors::SymmetricTensor<2> &extractor_col);
3311 
3313 
3314  }; // class ScalarFunction
3315 
3316 
3317 
3477  template <int dim,
3478  enum AD::NumberTypes ADNumberTypeCode,
3479  typename ScalarType = double>
3481  : public PointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>
3482  {
3483  public:
3488  using scalar_type =
3490 
3495  using ad_type =
3497 
3502 
3517  VectorFunction(const unsigned int n_independent_variables,
3518  const unsigned int n_dependent_variables);
3519 
3523  virtual ~VectorFunction() = default;
3524 
3526 
3531 
3544  void
3545  register_dependent_variables(const std::vector<ad_type> &funcs);
3546 
3565  template <typename ValueType, typename ExtractorType>
3566  void
3567  register_dependent_variable(const ValueType & funcs,
3568  const ExtractorType &extractor);
3569 
3578  void
3579  compute_values(Vector<scalar_type> &values) const;
3580 
3595  void
3596  compute_jacobian(FullMatrix<scalar_type> &jacobian) const;
3597 
3598 
3611  template <typename ExtractorType_Row>
3612  static typename internal::
3614  extract_value_component(const Vector<scalar_type> &values,
3615  const ExtractorType_Row & extractor_row);
3616 
3664  template <typename ExtractorType_Row, typename ExtractorType_Col>
3665  static typename internal::VectorFieldJacobian<dim,
3666  scalar_type,
3667  ExtractorType_Row,
3668  ExtractorType_Col>::type
3670  const ExtractorType_Row & extractor_row,
3671  const ExtractorType_Col & extractor_col);
3672 
3690  const FullMatrix<scalar_type> & jacobian,
3691  const FEValuesExtractors::Scalar &extractor_row,
3692  const FEValuesExtractors::Scalar &extractor_col);
3693 
3709  const FullMatrix<scalar_type> & jacobian,
3710  const FEValuesExtractors::SymmetricTensor<2> &extractor_row,
3711  const FEValuesExtractors::SymmetricTensor<2> &extractor_col);
3712 
3714 
3715  }; // class VectorFunction
3716 
3717 
3718  } // namespace AD
3719 } // namespace Differentiation
3720 
3721 
3722 /* ----------------- inline and template functions ----------------- */
3723 
3724 
3725 # ifndef DOXYGEN
3726 
3727 namespace Differentiation
3728 {
3729  namespace AD
3730  {
3731  /* ----------------- CellLevelBase ----------------- */
3732 
3733 
3734 
3735  template <enum AD::NumberTypes ADNumberTypeCode, typename ScalarType>
3736  template <typename VectorType>
3737  void
3739  const VectorType & values,
3740  const std::vector<::types::global_dof_index> &local_dof_indices)
3741  {
3742  // This is actually the same thing the set_dof_values() function,
3743  // in the sense that we simply populate our array of independent values
3744  // with a meaningful number. However, in this case we need to double check
3745  // that we're not registering these variables twice
3746  Assert(
3747  local_dof_indices.size() == this->n_independent_variables(),
3748  ExcMessage(
3749  "Degree of freedom index vector size does not match number of independent variables"));
3750  for (unsigned int i = 0; i < this->n_independent_variables(); ++i)
3751  {
3752  Assert(this->registered_independent_variable_values[i] == false,
3753  ExcMessage("Independent variables already registered."));
3754  }
3755  set_dof_values(values, local_dof_indices);
3756  }
3757 
3758 
3759 
3760  template <enum AD::NumberTypes ADNumberTypeCode, typename ScalarType>
3761  template <typename VectorType>
3762  void
3764  const VectorType & values,
3765  const std::vector<::types::global_dof_index> &local_dof_indices)
3766  {
3767  Assert(local_dof_indices.size() == this->n_independent_variables(),
3768  ExcMessage(
3769  "Vector size does not match number of independent variables"));
3770  for (unsigned int i = 0; i < this->n_independent_variables(); ++i)
3772  i, values[local_dof_indices[i]]);
3773  }
3774 
3775 
3776 
3777  /* ----------------- PointLevelFunctionsBase ----------------- */
3778 
3779 
3780 
3781  template <int dim,
3782  enum AD::NumberTypes ADNumberTypeCode,
3783  typename ScalarType>
3784  template <typename ValueType, typename ExtractorType>
3785  void
3787  register_independent_variable(const ValueType & value,
3788  const ExtractorType &extractor)
3789  {
3790  // This is actually the same thing as the set_independent_variable
3791  // function, in the sense that we simply populate our array of independent
3792  // values with a meaningful number. However, in this case we need to
3793  // double check that we're not registering these variables twice
3794 # ifdef DEBUG
3795  const std::vector<unsigned int> index_set(
3796  internal::extract_field_component_indices<dim>(extractor));
3797  for (const unsigned int index : index_set)
3798  {
3799  Assert(
3800  this->registered_independent_variable_values[index] == false,
3801  ExcMessage(
3802  "Overlapping indices for independent variables. "
3803  "One or more indices associated with the field that "
3804  "is being registered as an independent variable have "
3805  "already been associated with another field. This suggests "
3806  "that the component offsets used to construct their counterpart "
3807  "extractors are incompatible with one another. Make sure that "
3808  "the first component for each extractor properly takes into "
3809  "account the dimensionality of the preceding fields."));
3810  }
3811 # endif
3812  set_independent_variable(value, extractor);
3813  }
3814 
3815 
3816 
3817  template <int dim,
3818  enum AD::NumberTypes ADNumberTypeCode,
3819  typename ScalarType>
3820  template <typename ValueType, typename ExtractorType>
3821  void
3823  set_independent_variable(const ValueType & value,
3824  const ExtractorType &extractor)
3825  {
3826  const std::vector<unsigned int> index_set(
3827  internal::extract_field_component_indices<dim>(extractor));
3828  for (unsigned int i = 0; i < index_set.size(); ++i)
3829  {
3830  set_sensitivity_value(
3831  index_set[i],
3832  internal::Extractor<dim, ExtractorType>::symmetric_component(i),
3834  }
3835  }
3836 
3837 
3838 
3839  template <int dim,
3840  enum AD::NumberTypes ADNumberTypeCode,
3841  typename ScalarType>
3842  template <typename ExtractorType>
3843  typename internal::Extractor<dim, ExtractorType>::template tensor_type<
3844  typename PointLevelFunctionsBase<dim, ADNumberTypeCode, ScalarType>::
3845  ad_type>
3847  get_sensitive_variables(const ExtractorType &extractor) const
3848  {
3849  if (ADNumberTraits<ad_type>::is_taped == true)
3850  {
3851  Assert(this->active_tape_index() !=
3853  ExcMessage("Invalid tape index"));
3854  }
3855 
3856  // If necessary, finalize the internally stored vector of
3857  // AD numbers that represents the independent variables
3858  this->finalize_sensitive_independent_variables();
3859  Assert(this->independent_variables.size() ==
3860  this->n_independent_variables(),
3861  ExcDimensionMismatch(this->independent_variables.size(),
3862  this->n_independent_variables()));
3863 
3864  const std::vector<unsigned int> index_set(
3865  internal::extract_field_component_indices<dim>(extractor));
3866  typename internal::Extractor<dim,
3867  ExtractorType>::template tensor_type<ad_type>
3868  out;
3869 
3870  for (unsigned int i = 0; i < index_set.size(); ++i)
3871  {
3872  const unsigned int index = index_set[i];
3873  Assert(index < this->n_independent_variables(), ExcInternalError());
3874  Assert(this->registered_independent_variable_values[index] == true,
3875  ExcInternalError());
3876  internal::get_tensor_entry(out, i) =
3877  this->independent_variables[index];
3878  }
3879 
3880  return out;
3881  }
3882 
3883 
3884 
3885  /* ----------------- ScalarFunction ----------------- */
3886 
3887 
3888 
3889  template <int dim,
3890  enum AD::NumberTypes ADNumberTypeCode,
3891  typename ScalarType>
3892  template <typename ExtractorType_Row>
3893  typename internal::ScalarFieldGradient<
3894  dim,
3896  ExtractorType_Row>::type
3898  extract_gradient_component(const Vector<scalar_type> &gradient,
3899  const ExtractorType_Row & extractor_row)
3900  {
3901  // NOTE: The order of components must be consistently defined throughout
3902  // this class.
3903  typename internal::
3905 
3906  // Get indexsets for the subblock from which we wish to extract the
3907  // gradient values
3908  const std::vector<unsigned int> row_index_set(
3909  internal::extract_field_component_indices<dim>(extractor_row));
3910  Assert(out.n_independent_components == row_index_set.size(),
3911  ExcMessage("Not all tensor components have been extracted!"));
3912  for (unsigned int r = 0; r < row_index_set.size(); ++r)
3913  internal::set_tensor_entry(out, r, gradient[row_index_set[r]]);
3914 
3915  return out;
3916  }
3917 
3918 
3919 
3920  template <int dim,
3921  enum AD::NumberTypes ADNumberTypeCode,
3922  typename ScalarType>
3923  template <typename ExtractorType_Row, typename ExtractorType_Col>
3924  typename internal::ScalarFieldHessian<
3925  dim,
3927  ExtractorType_Row,
3928  ExtractorType_Col>::type
3931  const ExtractorType_Row & extractor_row,
3932  const ExtractorType_Col & extractor_col)
3933  {
3934  using InternalHessian = internal::ScalarFieldHessian<dim,
3935  scalar_type,
3936  ExtractorType_Row,
3937  ExtractorType_Col>;
3938  using InternalExtractorRow = internal::Extractor<dim, ExtractorType_Row>;
3939  using InternalExtractorCol = internal::Extractor<dim, ExtractorType_Col>;
3940  using HessianType = typename InternalHessian::type;
3941 
3942  // NOTE: The order of components must be consistently defined throughout
3943  // this class.
3944  HessianType out;
3945 
3946  // Get indexsets for the subblocks from which we wish to extract the
3947  // Hessian values
3948  // NOTE: Here we have to do some clever accounting when the
3949  // one extractor is a symmetric Tensor and the other is not, e.g.
3950  // <SymmTensor,Vector>. In this scenario the return type is a
3951  // non-symmetric Tensor<3,dim> but we have to fetch information from a
3952  // SymmTensor row/column that has too few entries to fill the output
3953  // tensor. So we must duplicate the relevant entries in the row/column
3954  // indexset to fetch off-diagonal components that are Otherwise
3955  // non-existent in a SymmTensor.
3956  const std::vector<unsigned int> row_index_set(
3957  internal::extract_field_component_indices<dim>(
3958  extractor_row, false /*ignore_symmetries*/));
3959  const std::vector<unsigned int> col_index_set(
3960  internal::extract_field_component_indices<dim>(
3961  extractor_col, false /*ignore_symmetries*/));
3962 
3963  for (unsigned int index = 0;
3964  index < HessianType::n_independent_components;
3965  ++index)
3966  {
3967  const TableIndices<HessianType::rank> ti_out =
3968  HessianType::unrolled_to_component_indices(index);
3969  const unsigned int r =
3970  InternalExtractorRow::local_component(ti_out, 0);
3971  const unsigned int c =
3972  InternalExtractorCol::local_component(ti_out,
3973  InternalExtractorRow::rank);
3974 
3976  out, index, hessian[row_index_set[r]][col_index_set[c]]);
3977  }
3978 
3979  return out;
3980  }
3981 
3982 
3983 
3984  /* ----------------- VectorFunction ----------------- */
3985 
3986 
3987 
3988  template <int dim,
3989  enum AD::NumberTypes ADNumberTypeCode,
3990  typename ScalarType>
3991  template <typename ValueType, typename ExtractorType>
3992  void
3994  register_dependent_variable(const ValueType & funcs,
3995  const ExtractorType &extractor)
3996  {
3997  const std::vector<unsigned int> index_set(
3998  internal::extract_field_component_indices<dim>(extractor));
3999  for (unsigned int i = 0; i < index_set.size(); ++i)
4000  {
4001  Assert(this->registered_marked_dependent_variables[index_set[i]] ==
4002  false,
4003  ExcMessage("Overlapping indices for dependent variables."));
4005  index_set[i], internal::get_tensor_entry(funcs, i));
4006  }
4007  }
4008 
4009 
4010 
4011  template <int dim,
4012  enum AD::NumberTypes ADNumberTypeCode,
4013  typename ScalarType>
4014  template <typename ExtractorType_Row>
4015  typename internal::VectorFieldValue<
4016  dim,
4018  ExtractorType_Row>::type
4020  const Vector<scalar_type> &values,
4021  const ExtractorType_Row & extractor_row)
4022  {
4023  // NOTE: The order of components must be consistently defined throughout
4024  // this class.
4025  typename internal::VectorFieldValue<dim, scalar_type, ExtractorType_Row>::
4026  type out;
4027 
4028  // Get indexsets for the subblock from which we wish to extract the
4029  // gradient values
4030  const std::vector<unsigned int> row_index_set(
4031  internal::extract_field_component_indices<dim>(extractor_row));
4032  Assert(out.n_independent_components == row_index_set.size(),
4033  ExcMessage("Not all tensor components have been extracted!"));
4034  for (unsigned int r = 0; r < row_index_set.size(); ++r)
4035  internal::set_tensor_entry(out, r, values[row_index_set[r]]);
4036 
4037  return out;
4038  }
4039 
4040 
4041 
4042  template <int dim,
4043  enum AD::NumberTypes ADNumberTypeCode,
4044  typename ScalarType>
4045  template <typename ExtractorType_Row, typename ExtractorType_Col>
4047  dim,
4049  ExtractorType_Row,
4050  ExtractorType_Col>::type
4053  const ExtractorType_Row & extractor_row,
4054  const ExtractorType_Col & extractor_col)
4055  {
4056  using InternalJacobian = internal::VectorFieldJacobian<dim,
4057  scalar_type,
4058  ExtractorType_Row,
4059  ExtractorType_Col>;
4060  using InternalExtractorRow = internal::Extractor<dim, ExtractorType_Row>;
4061  using InternalExtractorCol = internal::Extractor<dim, ExtractorType_Col>;
4062  using JacobianType = typename InternalJacobian::type;
4063 
4064  // NOTE: The order of components must be consistently defined throughout
4065  // this class.
4066  JacobianType out;
4067 
4068  // Get indexsets for the subblocks from which we wish to extract the
4069  // Hessian values.
4070  // NOTE: Here we have to do some clever accounting when the
4071  // one extractor is a symmetric Tensor and the other is not, e.g.
4072  // <SymmTensor,Vector>. In this scenario the return type is a
4073  // non-symmetric Tensor<3,dim> but we have to fetch information from a
4074  // SymmTensor row/column that has too few entries to fill the output
4075  // tensor. So we must duplicate the relevant entries in the row/column
4076  // indexset to fetch off-diagonal components that are Otherwise
4077  // non-existent in a SymmTensor.
4078  const std::vector<unsigned int> row_index_set(
4079  internal::extract_field_component_indices<dim>(
4080  extractor_row, false /*ignore_symmetries*/));
4081  const std::vector<unsigned int> col_index_set(
4082  internal::extract_field_component_indices<dim>(
4083  extractor_col, false /*ignore_symmetries*/));
4084 
4085  for (unsigned int index = 0;
4086  index < JacobianType::n_independent_components;
4087  ++index)
4088  {
4089  const TableIndices<JacobianType::rank> ti_out =
4090  JacobianType::unrolled_to_component_indices(index);
4091  const unsigned int r =
4092  InternalExtractorRow::local_component(ti_out, 0);
4093  const unsigned int c =
4094  InternalExtractorCol::local_component(ti_out,
4095  InternalExtractorRow::rank);
4096 
4098  out, index, jacobian[row_index_set[r]][col_index_set[c]]);
4099  }
4100 
4101  return out;
4102  }
4103 
4104 
4105  } // namespace AD
4106 } // namespace Differentiation
4107 
4108 
4109 # endif // DOXYGEN
4110 
4111 
4113 
4114 #endif // defined(DEAL_II_WITH_ADOLC) || defined(DEAL_II_TRILINOS_WITH_SACADO)
4115 
4116 #endif // dealii_differentiation_ad_ad_helpers_h
FEValuesExtractors::SymmetricTensor::first_tensor_component
unsigned int first_tensor_component
Definition: fe_values_extractors.h:204
Differentiation::AD::ScalarFunction::scalar_type
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
Definition: ad_helpers.h:3096
Differentiation::AD::ResidualLinearization::~ResidualLinearization
virtual ~ResidualLinearization()=default
Differentiation::AD::CellLevelBase::set_dof_values
void set_dof_values(const std::vector< scalar_type > &dof_values)
Definition: ad_helpers.cc:774
FEValuesExtractors
Definition: fe_values_extractors.h:82
Differentiation::AD::internal::VectorFieldValue
ScalarFieldGradient< dim, NumberType, ExtractorType_Field > VectorFieldValue
Definition: ad_helpers.h:2336
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Vector >::first_component
static unsigned int first_component(const FEValuesExtractors::Vector &extractor)
Definition: ad_helpers.h:1800
Differentiation::AD::HelperBase::set_tape_buffer_sizes
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:564
LinearAlgebraDealII::Vector
Vector< double > Vector
Definition: generic_linear_algebra.h:43
Differentiation::AD::HelperBase::reset
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:395
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Tensor< 2 > >::table_index_view
static TableIndices< rank > table_index_view(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:2032
Differentiation::AD::ResidualLinearization
Definition: ad_helpers.h:1539
TableIndices
Definition: table_indices.h:45
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Tensor< 2 > >::symmetric_component
static bool symmetric_component(const unsigned int unrolled_index)
Definition: ad_helpers.h:2020
Differentiation::AD::PointLevelFunctionsBase::n_symmetric_independent_variables
unsigned int n_symmetric_independent_variables() const
Definition: ad_helpers.cc:1211
Differentiation::AD::VectorFunction::register_dependent_variable
void register_dependent_variable(const ValueType &funcs, const ExtractorType &extractor)
SymmetricTensor
Definition: symmetric_tensor.h:611
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Tensor< 2 > >::local_component
static IndexType local_component(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:2058
Differentiation::AD::VectorFunction::compute_jacobian
void compute_jacobian(FullMatrix< scalar_type > &jacobian) const
Definition: ad_helpers.cc:1738
Differentiation::AD::VectorFunction::extract_jacobian_component
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)
Differentiation::AD::PointLevelFunctionsBase::~PointLevelFunctionsBase
virtual ~PointLevelFunctionsBase()=default
Differentiation::AD::internal::HessianType
Definition: ad_helpers.h:2206
Differentiation::AD::EnergyFunctional::compute_energy
scalar_type compute_energy() const
Definition: ad_helpers.cc:819
Differentiation::AD::HelperBase::configure_tapeless_mode
static void configure_tapeless_mode(const unsigned int n_independent_variables, const bool ensure_persistent_setting=true)
Definition: ad_helpers.cc:453
Differentiation::AD::VectorFunction::extract_value_component
static internal::VectorFieldValue< dim, scalar_type, ExtractorType_Row >::type extract_value_component(const Vector< scalar_type > &values, const ExtractorType_Row &extractor_row)
Differentiation::AD::HelperBase::independent_variable_values
std::vector< scalar_type > independent_variable_values
Definition: ad_helpers.h:638
Differentiation::AD::HelperBase::taped_driver
TapedDrivers< ad_type, scalar_type > taped_driver
Definition: ad_helpers.h:594
Differentiation::AD::PointLevelFunctionsBase::get_sensitive_variables
const std::vector< ad_type > & get_sensitive_variables() const
Definition: ad_helpers.cc:1251
FEValuesExtractors::Scalar
Definition: fe_values_extractors.h:95
VectorType
Differentiation::AD::internal::get_tensor_entry
NumberType get_tensor_entry(const TensorType< rank, dim, NumberType > &t, const unsigned int unrolled_index)
Definition: ad_helpers.h:2519
Differentiation::AD::internal::set_tensor_entry
void set_tensor_entry(TensorType &t, const unsigned int unrolled_index, const NumberType &value)
Definition: ad_helpers.h:2441
Differentiation::AD::HelperBase::reset_registered_dependent_variables
void reset_registered_dependent_variables(const bool flag=false)
Definition: ad_helpers.cc:94
AssertIndexRange
#define AssertIndexRange(index, range)
Definition: exceptions.h:1649
sacado_product_types.h
Differentiation::AD::HelperBase::set_sensitivity_value
void set_sensitivity_value(const unsigned int index, const scalar_type &value)
Definition: ad_helpers.cc:105
Differentiation::AD::HelperBase::~HelperBase
virtual ~HelperBase()=default
adolc_number_types.h
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Vector >::table_index_view
static TableIndices< rank > table_index_view(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:1826
Differentiation::AD::HelperBase::active_tape_index
Types< ad_type >::tape_index active_tape_index() const
Definition: ad_helpers.cc:284
Differentiation::AD::HelperBase::registered_marked_independent_variables
std::vector< bool > registered_marked_independent_variables
Definition: ad_helpers.h:660
Differentiation::AD::EnergyFunctional::compute_residual
void compute_residual(Vector< scalar_type > &residual) const override
Definition: ad_helpers.cc:873
DoFTools::n_components
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)
Differentiation::AD::ScalarFunction::compute_hessian
void compute_hessian(FullMatrix< scalar_type > &hessian) const
Definition: ad_helpers.cc:1477
Differentiation::AD::ScalarFunction
Definition: ad_helpers.h:3087
Differentiation::AD::ScalarFunction::extract_gradient_component
static internal::ScalarFieldGradient< dim, scalar_type, ExtractorType_Row >::type extract_gradient_component(const Vector< scalar_type > &gradient, const ExtractorType_Row &extractor_row)
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::SymmetricTensor< 2 > >::table_index_view
static TableIndices< rank > table_index_view(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:2139
Differentiation::AD::HelperBase::print_values
void print_values(std::ostream &stream) const
Definition: ad_helpers.cc:364
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Vector >::local_component
static IndexType local_component(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:1850
Differentiation::AD::PointLevelFunctionsBase::register_independent_variable
void register_independent_variable(const ValueType &value, const ExtractorType &extractor)
Differentiation::AD::PointLevelFunctionsBase::PointLevelFunctionsBase
PointLevelFunctionsBase(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
Definition: ad_helpers.cc:1160
Differentiation::AD::HelperBase< ADNumberTypeCode, double >::ad_type
typename AD::NumberTraits< double, ADNumberTypeCode >::ad_type ad_type
Definition: ad_helpers.h:187
Differentiation::AD::HelperBase
Definition: ad_helpers.h:172
fe_values_extractors.h
Differentiation::AD::internal::NumberType
Definition: numbers.h:653
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Tensor< 1 > >::symmetric_component
static bool symmetric_component(const unsigned int unrolled_index)
Definition: ad_helpers.h:1917
FEValuesExtractors::SymmetricTensor
Definition: fe_values_extractors.h:199
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Scalar >::symmetric_component
static bool symmetric_component(const unsigned int unrolled_index)
Definition: ad_helpers.h:1727
Tensor::unrolled_to_component_indices
static constexpr TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
Differentiation::AD::EnergyFunctional::register_energy_functional
void register_energy_functional(const ad_type &energy)
Definition: ad_helpers.cc:807
Differentiation::AD::CellLevelBase::get_sensitive_dof_values
const std::vector< ad_type > & get_sensitive_dof_values() const
Definition: ad_helpers.cc:749
Differentiation::AD::HelperBase::n_independent_variables
std::size_t n_independent_variables() const
Definition: ad_helpers.cc:241
Differentiation::AD::HelperBase::start_recording_operations
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:583
Differentiation::AD::internal::ScalarFieldGradient
Definition: ad_helpers.h:2184
Differentiation::AD::VectorFunction::compute_values
void compute_values(Vector< scalar_type > &values) const
Definition: ad_helpers.cc:1682
FEValuesExtractors::Vector::first_vector_component
unsigned int first_vector_component
Definition: fe_values_extractors.h:155
tensor.h
Differentiation::AD::internal::Extractor
Definition: ad_helpers.h:1660
Differentiation::AD::CellLevelBase::compute_residual
virtual void compute_residual(Vector< scalar_type > &residual) const =0
Differentiation::AD::PointLevelFunctionsBase
Definition: ad_helpers.h:2640
Differentiation::AD::HelperBase::stop_recording_operations
void stop_recording_operations(const bool write_tapes_to_file=false)
Definition: ad_helpers.cc:648
FEValuesExtractors::Vector
Definition: fe_values_extractors.h:150
Differentiation::AD::CellLevelBase::scalar_type
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
Definition: ad_helpers.h:848
Differentiation::AD::NumberTypes
NumberTypes
Definition: ad_number_types.h:36
StandardExceptions::ExcMessage
static ::ExceptionBase & ExcMessage(std::string arg1)
Differentiation::AD::HelperBase::finalize_sensitive_independent_variables
void finalize_sensitive_independent_variables() const
Definition: ad_helpers.cc:181
Differentiation::AD::HelperBase::dependent_variables
std::vector< ad_type > dependent_variables
Definition: ad_helpers.h:752
Differentiation::AD::HelperBase::n_registered_dependent_variables
unsigned int n_registered_dependent_variables() const
Definition: ad_helpers.cc:250
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::SymmetricTensor< 2 > >::symmetric_component
static bool symmetric_component(const unsigned int unrolled_index)
Definition: ad_helpers.h:2126
Differentiation::AD::VectorFunction::scalar_type
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
Definition: ad_helpers.h:3489
Differentiation::AD::CellLevelBase::CellLevelBase
CellLevelBase(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
Definition: ad_helpers.cc:715
Tensor
Definition: tensor.h:450
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
Differentiation::AD::CellLevelBase::compute_linearization
virtual void compute_linearization(FullMatrix< scalar_type > &linearization) const =0
ad_drivers.h
Differentiation
Definition: numbers.h:645
Differentiation::AD::internal::extract_field_component_indices
std::vector< IndexType > extract_field_component_indices(const ExtractorType &extractor, const bool ignore_symmetries=true)
Definition: ad_helpers.h:2367
Differentiation::AD::PointLevelFunctionsBase::dimension
static const unsigned int dimension
Definition: ad_helpers.h:2648
Differentiation::AD::ResidualLinearization::scalar_type
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
Definition: ad_helpers.h:1548
Differentiation::AD::HelperBase::registered_marked_dependent_variables
std::vector< bool > registered_marked_dependent_variables
Definition: ad_helpers.h:757
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Tensor< 2 > >::first_component
static unsigned int first_component(const FEValuesExtractors::Tensor< 2 > &extractor)
Definition: ad_helpers.h:2007
Differentiation::AD::HelperBase::register_dependent_variable
void register_dependent_variable(const unsigned int index, const ad_type &func)
Definition: ad_helpers.cc:684
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Scalar >::local_component
static IndexType local_component(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:1742
Differentiation::AD::ResidualLinearization::compute_residual
virtual void compute_residual(Vector< scalar_type > &residual) const override
Definition: ad_helpers.cc:1037
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::SymmetricTensor< 2 > >::local_component
static IndexType local_component(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:2165
Differentiation::AD::PointLevelFunctionsBase::set_independent_variables
void set_independent_variables(const std::vector< scalar_type > &values)
Definition: ad_helpers.cc:1301
Differentiation::AD::HelperBase::active_tape_requires_retaping
bool active_tape_requires_retaping() const
Definition: ad_helpers.cc:502
Differentiation::AD::ResidualLinearization::compute_linearization
virtual void compute_linearization(FullMatrix< scalar_type > &linearization) const override
Definition: ad_helpers.cc:1091
Differentiation::AD::EnergyFunctional::scalar_type
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
Definition: ad_helpers.h:1232
scalar_type
Differentiation::AD::ResidualLinearization::ResidualLinearization
ResidualLinearization(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
Definition: ad_helpers.cc:1011
Differentiation::AD::PointLevelFunctionsBase::symmetric_independent_variables
std::vector< bool > symmetric_independent_variables
Definition: ad_helpers.h:2920
Differentiation::AD::VectorFunction::~VectorFunction
virtual ~VectorFunction()=default
Differentiation::AD::HelperBase::scalar_type
typename AD::NumberTraits< ScalarType, ADNumberTypeCode >::scalar_type scalar_type
Definition: ad_helpers.h:180
Differentiation::AD::HelperBase::initialize_non_sensitive_independent_variable
void initialize_non_sensitive_independent_variable(const unsigned int index, ad_type &out) const
Definition: ad_helpers.cc:205
Differentiation::AD::EnergyFunctional::ad_type
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
Definition: ad_helpers.h:1239
Differentiation::AD::HelperBase::reset_registered_independent_variables
void reset_registered_independent_variables()
Definition: ad_helpers.cc:82
Differentiation::AD::PointLevelFunctionsBase::register_independent_variables
void register_independent_variables(const std::vector< scalar_type > &values)
Definition: ad_helpers.cc:1225
Differentiation::AD::Numbers::invalid_tape_index
static const Types< ADNumberType >::tape_index invalid_tape_index
Definition: ad_drivers.h:121
unsigned int
value
static const bool value
Definition: dof_tools_constraints.cc:433
StandardExceptions::ExcInternalError
static ::ExceptionBase & ExcInternalError()
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Tensor< 1 > >::first_component
static unsigned int first_component(const FEValuesExtractors::Tensor< 1 > &extractor)
Definition: ad_helpers.h:1904
Differentiation::AD::PointLevelFunctionsBase::set_sensitivity_value
void set_sensitivity_value(const unsigned int index, const bool symmetric_component, const scalar_type &value)
Definition: ad_helpers.cc:1279
Differentiation::AD::ScalarFunction::register_dependent_variable
void register_dependent_variable(const ad_type &func)
Definition: ad_helpers.cc:1340
Differentiation::AD::VectorFunction::VectorFunction
VectorFunction(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
Definition: ad_helpers.cc:1651
sacado_number_types.h
Differentiation::AD::ScalarFunction::ad_type
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
Definition: ad_helpers.h:3103
Assert
#define Assert(cond, exc)
Definition: exceptions.h:1419
Differentiation::AD::HelperBase::clear_active_tape
void clear_active_tape()
Definition: ad_helpers.cc:515
ad_number_traits.h
symmetric_tensor.h
Differentiation::AD::HelperBase::activate_recorded_tape
void activate_recorded_tape(const typename Types< ad_type >::tape_index tape_index)
Definition: ad_helpers.cc:479
Differentiation::AD::CellLevelBase
Definition: ad_helpers.h:840
vector.h
Differentiation::AD::HelperBase::recorded_tape_requires_retaping
bool recorded_tape_requires_retaping(const typename Types< ad_type >::tape_index tape_index) const
Definition: ad_helpers.cc:489
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Scalar >::first_component
static unsigned int first_component(const FEValuesExtractors::Scalar &extractor)
Definition: ad_helpers.h:1714
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Tensor< 1 > >::local_component
static IndexType local_component(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:1953
Differentiation::AD::HelperBase::print
void print(std::ostream &stream) const
Definition: ad_helpers.cc:309
Differentiation::AD::internal::ScalarFieldHessian::rank
static const int rank
Definition: ad_helpers.h:2315
Differentiation::AD::ScalarFunction::extract_hessian_component
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)
Differentiation::AD::ScalarFunction::~ScalarFunction
virtual ~ScalarFunction()=default
numbers::invalid_unsigned_int
static const unsigned int invalid_unsigned_int
Definition: types.h:191
Differentiation::AD::EnergyFunctional::compute_linearization
virtual void compute_linearization(FullMatrix< scalar_type > &linearization) const override
Definition: ad_helpers.cc:938
FEValuesExtractors::Tensor::first_tensor_component
unsigned int first_tensor_component
Definition: fe_values_extractors.h:253
Differentiation::AD::VectorFunction::register_dependent_variables
void register_dependent_variables(const std::vector< ad_type > &funcs)
Definition: ad_helpers.cc:1666
Differentiation::AD::HelperBase::print_tape_stats
void print_tape_stats(const typename Types< ad_type >::tape_index tape_index, std::ostream &stream) const
Definition: ad_helpers.cc:378
Differentiation::AD::HelperBase::is_registered_tape
bool is_registered_tape(const typename Types< ad_type >::tape_index tape_index) const
Definition: ad_helpers.cc:296
StandardExceptions::ExcDimensionMismatch
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
Differentiation::AD::PointLevelFunctionsBase::is_symmetric_independent_variable
bool is_symmetric_independent_variable(const unsigned int index) const
Definition: ad_helpers.cc:1197
Differentiation::AD::HelperBase::registered_independent_variable_values
std::vector< bool > registered_independent_variable_values
Definition: ad_helpers.h:654
FEValuesExtractors::Scalar::component
unsigned int component
Definition: fe_values_extractors.h:100
Differentiation::AD::HelperBase::n_registered_independent_variables
unsigned int n_registered_independent_variables() const
Definition: ad_helpers.cc:230
Differentiation::AD::VectorFunction::ad_type
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
Definition: ad_helpers.h:3496
Differentiation::AD::HelperBase::independent_variables
std::vector< ad_type > independent_variables
Definition: ad_helpers.h:648
Differentiation::AD::internal::ScalarFieldHessian
Definition: ad_helpers.h:2310
Differentiation::AD::HelperBase::n_dependent_variables
std::size_t n_dependent_variables() const
Definition: ad_helpers.cc:262
config.h
Differentiation::AD::ResidualLinearization::ad_type
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
Definition: ad_helpers.h:1555
Differentiation::AD::TapelessDrivers< ad_type, scalar_type >
FullMatrix
Definition: full_matrix.h:71
FEValuesExtractors::Tensor
Definition: fe_values_extractors.h:248
internal
Definition: aligned_vector.h:369
Differentiation::AD::internal::ScalarFieldGradient::type
typename Extractor< dim, ExtractorType >::template tensor_type< NumberType > type
Definition: ad_helpers.h:2192
Differentiation::AD::EnergyFunctional::EnergyFunctional
EnergyFunctional(const unsigned int n_independent_variables)
Definition: ad_helpers.cc:798
adolc_product_types.h
Differentiation::AD::PointLevelFunctionsBase::reset
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:1173
Differentiation::AD::EnergyFunctional::~EnergyFunctional
virtual ~EnergyFunctional()=default
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
Differentiation::AD::NumberTraits
Definition: ad_number_traits.h:53
Differentiation::AD::ScalarFunction::compute_value
scalar_type compute_value() const
Definition: ad_helpers.cc:1353
Differentiation::AD::HelperBase::tapeless_driver
TapelessDrivers< ad_type, scalar_type > tapeless_driver
Definition: ad_helpers.h:602
Differentiation::AD::ScalarFunction::compute_gradient
void compute_gradient(Vector< scalar_type > &gradient) const
Definition: ad_helpers.cc:1403
Differentiation::AD::EnergyFunctional
Definition: ad_helpers.h:1224
Differentiation::AD::CellLevelBase::register_dof_values
void register_dof_values(const std::vector< scalar_type > &dof_values)
Definition: ad_helpers.cc:726
Differentiation::AD::ScalarFunction::ScalarFunction
ScalarFunction(const unsigned int n_independent_variables)
Definition: ad_helpers.cc:1326
Differentiation::AD::VectorFunction
Definition: ad_helpers.h:3480
Differentiation::AD::internal::ScalarFieldHessian::type
typename HessianType< ExtractorType_Row, ExtractorType_Col >::template type< rank, dim, NumberType > type
Definition: ad_helpers.h:2326
Differentiation::AD::CellLevelBase::~CellLevelBase
virtual ~CellLevelBase()=default
Differentiation::AD::TapedDrivers< ad_type, scalar_type >
Differentiation::AD::HelperBase::activate_tape
void activate_tape(const typename Types< ad_type >::tape_index tape_index, const bool read_mode)
Definition: ad_helpers.cc:527
Differentiation::AD::PointLevelFunctionsBase::set_independent_variable
void set_independent_variable(const ValueType &value, const ExtractorType &extractor)
numbers.h
full_matrix.h
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Tensor< 1 > >::table_index_view
static TableIndices< rank > table_index_view(const TableIndices< rank_in > &table_indices, const unsigned int column_offset)
Definition: ad_helpers.h:1929
Differentiation::AD::HelperBase::mark_independent_variable
void mark_independent_variable(const unsigned int index, ad_type &out) const
Definition: ad_helpers.cc:142
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::SymmetricTensor< 2 > >::first_component
static unsigned int first_component(const FEValuesExtractors::SymmetricTensor< 2 > &extractor)
Definition: ad_helpers.h:2112
Differentiation::AD::ResidualLinearization::register_residual_vector
void register_residual_vector(const std::vector< ad_type > &residual)
Definition: ad_helpers.cc:1023
Differentiation::AD::HelperBase::is_recording
bool is_recording() const
Definition: ad_helpers.cc:271
int
Differentiation::AD::internal::Extractor< dim, FEValuesExtractors::Vector >::symmetric_component
static bool symmetric_component(const unsigned int unrolled_index)
Definition: ad_helpers.h:1814
Differentiation::AD::HelperBase::HelperBase
HelperBase(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
Definition: ad_helpers.cc:38
SymmetricTensor::component_to_unrolled_index
static constexpr unsigned int component_to_unrolled_index(const TableIndices< rank_ > &indices)
Differentiation::AD::internal::VectorFieldJacobian
ScalarFieldHessian< dim, NumberType, ExtractorType_Field, ExtractorType_Derivative > VectorFieldJacobian
Definition: ad_helpers.h:2355