Reference documentation for deal.II version GIT a2efd28e10 2023-02-01 14:45:02+00:00
\(\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\}}\)
fe_values.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 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.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_fe_values_h
17 #define dealii_fe_values_h
18 
19 
20 #include <deal.II/base/config.h>
21 
24 #include <deal.II/base/point.h>
29 
32 
33 #include <deal.II/fe/fe.h>
36 #include <deal.II/fe/mapping.h>
38 
39 #include <deal.II/grid/tria.h>
41 
43 
44 #include <algorithm>
45 #include <memory>
46 #include <type_traits>
47 
48 
49 // dummy include in order to have the
50 // definition of PetscScalar available
51 // without including other PETSc stuff
52 #ifdef DEAL_II_WITH_PETSC
53 # include <petsc.h>
54 #endif
55 
57 
58 // Forward declaration
59 #ifndef DOXYGEN
60 template <int dim, int spacedim = dim>
61 class FEValuesBase;
62 #endif
63 
64 namespace internal
65 {
70  template <int dim, class NumberType = double>
71  struct CurlType;
72 
79  template <class NumberType>
80  struct CurlType<1, NumberType>
81  {
83  };
84 
91  template <class NumberType>
92  struct CurlType<2, NumberType>
93  {
95  };
96 
103  template <class NumberType>
104  struct CurlType<3, NumberType>
105  {
107  };
108 } // namespace internal
109 
110 
111 
133 namespace FEValuesViews
134 {
146  template <int dim, int spacedim = dim>
147  class Scalar
148  {
149  public:
156 
163 
170 
177 
184  template <typename Number>
186 
193  template <typename Number>
196 
203  template <typename Number>
206 
213  template <typename Number>
216 
223  template <typename Number>
226 
233  template <typename Number>
235  {
240  using value_type =
241  typename ProductType<Number,
242  typename Scalar<dim, spacedim>::value_type>::type;
243 
248  using gradient_type = typename ProductType<
249  Number,
250  typename Scalar<dim, spacedim>::gradient_type>::type;
251 
257  typename ProductType<Number,
258  typename Scalar<dim, spacedim>::value_type>::type;
259 
264  using hessian_type = typename ProductType<
265  Number,
266  typename Scalar<dim, spacedim>::hessian_type>::type;
267 
273  Number,
275  };
276 
282  {
292 
301  unsigned int row_index;
302  };
303 
307  Scalar();
308 
314  Scalar(const FEValuesBase<dim, spacedim> &fe_values_base,
315  const unsigned int component);
316 
321  Scalar(const Scalar<dim, spacedim> &) = delete;
322 
326  // NOLINTNEXTLINE OSX does not compile with noexcept
328 
332  ~Scalar() = default;
333 
338  Scalar &
339  operator=(const Scalar<dim, spacedim> &) = delete;
340 
344  Scalar &
345  operator=(Scalar<dim, spacedim> &&) noexcept = default;
346 
360  value_type
361  value(const unsigned int shape_function, const unsigned int q_point) const;
362 
374  gradient(const unsigned int shape_function,
375  const unsigned int q_point) const;
376 
388  hessian(const unsigned int shape_function,
389  const unsigned int q_point) const;
390 
402  third_derivative(const unsigned int shape_function,
403  const unsigned int q_point) const;
404 
422  template <class InputVector>
423  void
425  const InputVector &fe_function,
426  std::vector<solution_value_type<typename InputVector::value_type>>
427  &values) const;
428 
463  template <class InputVector>
464  void
466  const InputVector &dof_values,
467  std::vector<solution_value_type<typename InputVector::value_type>>
468  &values) const;
469 
487  template <class InputVector>
488  void
490  const InputVector &fe_function,
491  std::vector<solution_gradient_type<typename InputVector::value_type>>
492  &gradients) const;
493 
500  template <class InputVector>
501  void
503  const InputVector &dof_values,
504  std::vector<solution_gradient_type<typename InputVector::value_type>>
505  &gradients) const;
506 
524  template <class InputVector>
525  void
527  const InputVector &fe_function,
528  std::vector<solution_hessian_type<typename InputVector::value_type>>
529  &hessians) const;
530 
537  template <class InputVector>
538  void
540  const InputVector &dof_values,
541  std::vector<solution_hessian_type<typename InputVector::value_type>>
542  &hessians) const;
543 
544 
563  template <class InputVector>
564  void
566  const InputVector &fe_function,
567  std::vector<solution_laplacian_type<typename InputVector::value_type>>
568  &laplacians) const;
569 
576  template <class InputVector>
577  void
579  const InputVector &dof_values,
580  std::vector<solution_laplacian_type<typename InputVector::value_type>>
581  &laplacians) const;
582 
583 
602  template <class InputVector>
603  void
605  const InputVector &fe_function,
606  std::vector<
607  solution_third_derivative_type<typename InputVector::value_type>>
608  &third_derivatives) const;
609 
616  template <class InputVector>
617  void
619  const InputVector &dof_values,
620  std::vector<
621  solution_third_derivative_type<typename InputVector::value_type>>
622  &third_derivatives) const;
623 
624 
625  private:
629  const SmartPointer<const FEValuesBase<dim, spacedim>> fe_values;
630 
635  const unsigned int component;
636 
641  };
642 
643 
644 
674  template <int dim, int spacedim = dim>
675  class Vector
676  {
677  public:
684 
694 
706 
713 
720  using curl_type = typename ::internal::CurlType<spacedim>::type;
721 
728 
735 
742  template <typename Number>
744 
751  template <typename Number>
754 
761  template <typename Number>
764 
771  template <typename Number>
774 
781  template <typename Number>
784 
791  template <typename Number>
793 
800  template <typename Number>
803 
810  template <typename Number>
813 
820  template <typename Number>
822  {
827  using value_type =
828  typename ProductType<Number,
829  typename Vector<dim, spacedim>::value_type>::type;
830 
835  using gradient_type = typename ProductType<
836  Number,
837  typename Vector<dim, spacedim>::gradient_type>::type;
838 
844  Number,
846 
851  using divergence_type = typename ProductType<
852  Number,
854 
860  typename ProductType<Number,
861  typename Vector<dim, spacedim>::value_type>::type;
862 
867  using curl_type =
868  typename ProductType<Number,
869  typename Vector<dim, spacedim>::curl_type>::type;
870 
875  using hessian_type = typename ProductType<
876  Number,
877  typename Vector<dim, spacedim>::hessian_type>::type;
878 
884  Number,
886  };
887 
893  {
903 
913  unsigned int row_index[spacedim];
914 
925  };
926 
930  Vector();
931 
940  Vector(const FEValuesBase<dim, spacedim> &fe_values_base,
941  const unsigned int first_vector_component);
942 
947  Vector(const Vector<dim, spacedim> &) = delete;
948 
952  // NOLINTNEXTLINE OSX does not compile with noexcept
954 
958  ~Vector() = default;
959 
964  Vector &
965  operator=(const Vector<dim, spacedim> &) = delete;
966 
970  // NOLINTNEXTLINE OSX does not compile with noexcept
971  Vector &
972  operator=(Vector<dim, spacedim> &&) = default; // NOLINT
973 
990  value_type
991  value(const unsigned int shape_function, const unsigned int q_point) const;
992 
1007  gradient(const unsigned int shape_function,
1008  const unsigned int q_point) const;
1009 
1026  symmetric_gradient(const unsigned int shape_function,
1027  const unsigned int q_point) const;
1028 
1040  divergence(const unsigned int shape_function,
1041  const unsigned int q_point) const;
1042 
1063  curl_type
1064  curl(const unsigned int shape_function, const unsigned int q_point) const;
1065 
1076  hessian_type
1077  hessian(const unsigned int shape_function,
1078  const unsigned int q_point) const;
1079 
1091  third_derivative(const unsigned int shape_function,
1092  const unsigned int q_point) const;
1093 
1111  template <class InputVector>
1112  void
1114  const InputVector &fe_function,
1116  &values) const;
1117 
1152  template <class InputVector>
1153  void
1155  const InputVector &dof_values,
1157  &values) const;
1158 
1176  template <class InputVector>
1177  void
1179  const InputVector &fe_function,
1181  &gradients) const;
1182 
1189  template <class InputVector>
1190  void
1192  const InputVector &dof_values,
1194  &gradients) const;
1195 
1219  template <class InputVector>
1220  void
1221  get_function_symmetric_gradients(
1222  const InputVector &fe_function,
1223  std::vector<
1225  &symmetric_gradients) const;
1226 
1233  template <class InputVector>
1234  void
1235  get_function_symmetric_gradients_from_local_dof_values(
1236  const InputVector &dof_values,
1237  std::vector<
1239  &symmetric_gradients) const;
1240 
1259  template <class InputVector>
1260  void
1261  get_function_divergences(
1262  const InputVector &fe_function,
1264  &divergences) const;
1265 
1272  template <class InputVector>
1273  void
1274  get_function_divergences_from_local_dof_values(
1275  const InputVector &dof_values,
1277  &divergences) const;
1278 
1297  template <class InputVector>
1298  void
1299  get_function_curls(
1300  const InputVector &fe_function,
1302  const;
1303 
1310  template <class InputVector>
1311  void
1312  get_function_curls_from_local_dof_values(
1313  const InputVector &dof_values,
1315  const;
1316 
1334  template <class InputVector>
1335  void
1337  const InputVector &fe_function,
1339  &hessians) const;
1340 
1347  template <class InputVector>
1348  void
1350  const InputVector &dof_values,
1352  &hessians) const;
1353 
1372  template <class InputVector>
1373  void
1375  const InputVector &fe_function,
1377  &laplacians) const;
1378 
1385  template <class InputVector>
1386  void
1388  const InputVector &dof_values,
1390  &laplacians) const;
1391 
1410  template <class InputVector>
1411  void
1413  const InputVector &fe_function,
1414  std::vector<
1416  &third_derivatives) const;
1417 
1424  template <class InputVector>
1425  void
1427  const InputVector &dof_values,
1428  std::vector<
1430  &third_derivatives) const;
1431 
1432  private:
1437 
1442  const unsigned int first_vector_component;
1443 
1447  std::vector<ShapeFunctionData> shape_function_data;
1448  };
1449 
1450 
1451  template <int rank, int dim, int spacedim = dim>
1453 
1476  template <int dim, int spacedim>
1477  class SymmetricTensor<2, dim, spacedim>
1478  {
1479  public:
1487 
1498 
1505  template <typename Number>
1507 
1514  template <typename Number>
1517 
1518 
1525  template <typename Number>
1526  struct DEAL_II_DEPRECATED OutputType
1527  {
1532  using value_type = typename ProductType<
1533  Number,
1535 
1540  using divergence_type = typename ProductType<
1541  Number,
1543  };
1544 
1549  struct ShapeFunctionData
1550  {
1559  bool is_nonzero_shape_function_component
1560  [value_type::n_independent_components];
1561 
1571  unsigned int row_index[value_type::n_independent_components];
1572 
1582 
1587  };
1588 
1592  SymmetricTensor();
1593 
1603  SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1604  const unsigned int first_tensor_component);
1605 
1611 
1615  // NOLINTNEXTLINE OSX does not compile with noexcept
1617 
1622  SymmetricTensor &
1624 
1628  SymmetricTensor &
1630 
1648  value_type
1649  value(const unsigned int shape_function, const unsigned int q_point) const;
1650 
1665  divergence(const unsigned int shape_function,
1666  const unsigned int q_point) const;
1667 
1685  template <class InputVector>
1686  void
1688  const InputVector &fe_function,
1689  std::vector<solution_value_type<typename InputVector::value_type>>
1690  &values) const;
1691 
1726  template <class InputVector>
1727  void
1729  const InputVector &dof_values,
1730  std::vector<solution_value_type<typename InputVector::value_type>>
1731  &values) const;
1732 
1754  template <class InputVector>
1755  void
1756  get_function_divergences(
1757  const InputVector &fe_function,
1758  std::vector<solution_divergence_type<typename InputVector::value_type>>
1759  &divergences) const;
1760 
1767  template <class InputVector>
1768  void
1769  get_function_divergences_from_local_dof_values(
1770  const InputVector &dof_values,
1771  std::vector<solution_divergence_type<typename InputVector::value_type>>
1772  &divergences) const;
1773 
1774  private:
1778  const SmartPointer<const FEValuesBase<dim, spacedim>> fe_values;
1779 
1784  const unsigned int first_tensor_component;
1785 
1789  std::vector<ShapeFunctionData> shape_function_data;
1790  };
1791 
1792 
1793  template <int rank, int dim, int spacedim = dim>
1794  class Tensor;
1795 
1814  template <int dim, int spacedim>
1815  class Tensor<2, dim, spacedim>
1816  {
1817  public:
1823 
1828 
1834 
1841  template <typename Number>
1843 
1850  template <typename Number>
1853 
1860  template <typename Number>
1863 
1864 
1871  template <typename Number>
1872  struct DEAL_II_DEPRECATED OutputType
1873  {
1878  using value_type = typename ProductType<
1879  Number,
1880  typename Tensor<2, dim, spacedim>::value_type>::type;
1881 
1886  using divergence_type = typename ProductType<
1887  Number,
1889 
1894  using gradient_type = typename ProductType<
1895  Number,
1897  };
1898 
1903  struct ShapeFunctionData
1904  {
1913  bool is_nonzero_shape_function_component
1914  [value_type::n_independent_components];
1915 
1925  unsigned int row_index[value_type::n_independent_components];
1926 
1936 
1941  };
1942 
1946  Tensor();
1947 
1953 
1957  // NOLINTNEXTLINE OSX does not compile with noexcept
1959 
1963  ~Tensor() = default;
1964 
1974  Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1975  const unsigned int first_tensor_component);
1976 
1977 
1982  Tensor &
1984 
1988  Tensor &
1989  operator=(Tensor<2, dim, spacedim> &&) = default; // NOLINT
1990 
2007  value_type
2008  value(const unsigned int shape_function, const unsigned int q_point) const;
2009 
2024  divergence(const unsigned int shape_function,
2025  const unsigned int q_point) const;
2026 
2041  gradient(const unsigned int shape_function,
2042  const unsigned int q_point) const;
2043 
2061  template <class InputVector>
2062  void
2064  const InputVector &fe_function,
2066  &values) const;
2067 
2102  template <class InputVector>
2103  void
2105  const InputVector &dof_values,
2107  &values) const;
2108 
2130  template <class InputVector>
2131  void
2132  get_function_divergences(
2133  const InputVector &fe_function,
2135  &divergences) const;
2136 
2143  template <class InputVector>
2144  void
2145  get_function_divergences_from_local_dof_values(
2146  const InputVector &dof_values,
2148  &divergences) const;
2149 
2166  template <class InputVector>
2167  void
2169  const InputVector &fe_function,
2171  &gradients) const;
2172 
2179  template <class InputVector>
2180  void
2182  const InputVector &dof_values,
2184  &gradients) const;
2185 
2186  private:
2191 
2196  const unsigned int first_tensor_component;
2197 
2201  std::vector<ShapeFunctionData> shape_function_data;
2202  };
2203 
2204 } // namespace FEValuesViews
2205 
2206 
2207 namespace internal
2208 {
2209  namespace FEValuesViews
2210  {
2215  template <int dim, int spacedim, typename Extractor>
2216  struct ViewType
2217  {};
2218 
2226  template <int dim, int spacedim>
2227  struct ViewType<dim, spacedim, FEValuesExtractors::Scalar>
2228  {
2229  using type = typename ::FEValuesViews::Scalar<dim, spacedim>;
2230  };
2231 
2239  template <int dim, int spacedim>
2240  struct ViewType<dim, spacedim, FEValuesExtractors::Vector>
2241  {
2242  using type = typename ::FEValuesViews::Vector<dim, spacedim>;
2243  };
2244 
2252  template <int dim, int spacedim, int rank>
2253  struct ViewType<dim, spacedim, FEValuesExtractors::Tensor<rank>>
2254  {
2255  using type = typename ::FEValuesViews::Tensor<rank, dim, spacedim>;
2256  };
2257 
2265  template <int dim, int spacedim, int rank>
2266  struct ViewType<dim, spacedim, FEValuesExtractors::SymmetricTensor<rank>>
2267  {
2268  using type =
2269  typename ::FEValuesViews::SymmetricTensor<rank, dim, spacedim>;
2270  };
2271 
2279  template <int dim, int spacedim>
2280  struct Cache
2281  {
2286  std::vector<::FEValuesViews::Scalar<dim, spacedim>> scalars;
2287  std::vector<::FEValuesViews::Vector<dim, spacedim>> vectors;
2288  std::vector<::FEValuesViews::SymmetricTensor<2, dim, spacedim>>
2290  std::vector<::FEValuesViews::Tensor<2, dim, spacedim>>
2292 
2297  };
2298  } // namespace FEValuesViews
2299 } // namespace internal
2300 
2301 namespace FEValuesViews
2302 {
2307  template <int dim, int spacedim, typename Extractor>
2308  using View = typename ::internal::FEValuesViews::
2309  ViewType<dim, spacedim, Extractor>::type;
2310 } // namespace FEValuesViews
2311 
2312 
2412 template <int dim, int spacedim>
2414 {
2415 public:
2419  static constexpr unsigned int dimension = dim;
2420 
2424  static constexpr unsigned int space_dimension = spacedim;
2425 
2433  const unsigned int n_quadrature_points;
2434 
2444  const unsigned int max_n_quadrature_points;
2445 
2451  const unsigned int dofs_per_cell;
2452 
2453 
2461  FEValuesBase(const unsigned int n_q_points,
2462  const unsigned int dofs_per_cell,
2463  const UpdateFlags update_flags,
2464  const Mapping<dim, spacedim> & mapping,
2465  const FiniteElement<dim, spacedim> &fe);
2466 
2471  FEValuesBase &
2472  operator=(const FEValuesBase &) = delete;
2473 
2478  FEValuesBase(const FEValuesBase &) = delete;
2479 
2483  virtual ~FEValuesBase() override;
2484 
2485 
2489 
2511  const double &
2512  shape_value(const unsigned int function_no,
2513  const unsigned int point_no) const;
2514 
2535  double
2536  shape_value_component(const unsigned int function_no,
2537  const unsigned int point_no,
2538  const unsigned int component) const;
2539 
2565  const Tensor<1, spacedim> &
2566  shape_grad(const unsigned int function_no,
2567  const unsigned int quadrature_point) const;
2568 
2586  shape_grad_component(const unsigned int function_no,
2587  const unsigned int point_no,
2588  const unsigned int component) const;
2589 
2609  const Tensor<2, spacedim> &
2610  shape_hessian(const unsigned int function_no,
2611  const unsigned int point_no) const;
2612 
2630  shape_hessian_component(const unsigned int function_no,
2631  const unsigned int point_no,
2632  const unsigned int component) const;
2633 
2653  const Tensor<3, spacedim> &
2654  shape_3rd_derivative(const unsigned int function_no,
2655  const unsigned int point_no) const;
2656 
2674  shape_3rd_derivative_component(const unsigned int function_no,
2675  const unsigned int point_no,
2676  const unsigned int component) const;
2677 
2680 
2718  template <class InputVector>
2719  void
2721  const InputVector & fe_function,
2722  std::vector<typename InputVector::value_type> &values) const;
2723 
2737  template <class InputVector>
2738  void
2740  const InputVector & fe_function,
2741  std::vector<Vector<typename InputVector::value_type>> &values) const;
2742 
2799  template <class InputVector>
2800  void
2802  const InputVector & fe_function,
2804  std::vector<typename InputVector::value_type> & values) const;
2805 
2814  template <class InputVector>
2815  void
2817  const InputVector & fe_function,
2819  std::vector<Vector<typename InputVector::value_type>> &values) const;
2820 
2821 
2843  template <class InputVector>
2844  void
2846  const InputVector & fe_function,
2848  ArrayView<std::vector<typename InputVector::value_type>> values,
2849  const bool quadrature_points_fastest) const;
2850 
2853 
2891  template <class InputVector>
2892  void
2894  const InputVector &fe_function,
2896  &gradients) const;
2897 
2914  template <class InputVector>
2915  void
2917  const InputVector &fe_function,
2918  std::vector<
2920  &gradients) const;
2921 
2930  template <class InputVector>
2931  void
2933  const InputVector & fe_function,
2936  &gradients) const;
2937 
2946  template <class InputVector>
2947  void
2949  const InputVector & fe_function,
2951  ArrayView<
2953  gradients,
2954  const bool quadrature_points_fastest = false) const;
2955 
2960 
2999  template <class InputVector>
3000  void
3002  const InputVector &fe_function,
3004  &hessians) const;
3005 
3023  template <class InputVector>
3024  void
3026  const InputVector &fe_function,
3027  std::vector<
3029  & hessians,
3030  const bool quadrature_points_fastest = false) const;
3031 
3040  template <class InputVector>
3041  void
3043  const InputVector & fe_function,
3046  &hessians) const;
3047 
3056  template <class InputVector>
3057  void
3059  const InputVector & fe_function,
3061  ArrayView<
3063  hessians,
3064  const bool quadrature_points_fastest = false) const;
3065 
3106  template <class InputVector>
3107  void
3109  const InputVector & fe_function,
3110  std::vector<typename InputVector::value_type> &laplacians) const;
3111 
3131  template <class InputVector>
3132  void
3134  const InputVector & fe_function,
3135  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
3136 
3145  template <class InputVector>
3146  void
3148  const InputVector & fe_function,
3150  std::vector<typename InputVector::value_type> & laplacians) const;
3151 
3160  template <class InputVector>
3161  void
3163  const InputVector & fe_function,
3165  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
3166 
3175  template <class InputVector>
3176  void
3178  const InputVector & fe_function,
3180  std::vector<std::vector<typename InputVector::value_type>> &laplacians,
3181  const bool quadrature_points_fastest = false) const;
3182 
3185 
3225  template <class InputVector>
3226  void
3228  const InputVector &fe_function,
3230  &third_derivatives) const;
3231 
3250  template <class InputVector>
3251  void
3253  const InputVector &fe_function,
3254  std::vector<
3256  & third_derivatives,
3257  const bool quadrature_points_fastest = false) const;
3258 
3267  template <class InputVector>
3268  void
3270  const InputVector & fe_function,
3273  &third_derivatives) const;
3274 
3283  template <class InputVector>
3284  void
3286  const InputVector & fe_function,
3288  ArrayView<
3290  third_derivatives,
3291  const bool quadrature_points_fastest = false) const;
3295 
3321  dof_indices() const;
3322 
3356  dof_indices_starting_at(const unsigned int start_dof_index) const;
3357 
3389  dof_indices_ending_at(const unsigned int end_dof_index) const;
3390 
3394 
3418 
3424  const Point<spacedim> &
3425  quadrature_point(const unsigned int q) const;
3426 
3432  const std::vector<Point<spacedim>> &
3434 
3450  double
3451  JxW(const unsigned int quadrature_point) const;
3452 
3456  const std::vector<double> &
3458 
3466  jacobian(const unsigned int quadrature_point) const;
3467 
3474  const std::vector<DerivativeForm<1, dim, spacedim>> &
3475  get_jacobians() const;
3476 
3485  jacobian_grad(const unsigned int quadrature_point) const;
3486 
3493  const std::vector<DerivativeForm<2, dim, spacedim>> &
3495 
3504  const Tensor<3, spacedim> &
3505  jacobian_pushed_forward_grad(const unsigned int quadrature_point) const;
3506 
3513  const std::vector<Tensor<3, spacedim>> &
3515 
3524  jacobian_2nd_derivative(const unsigned int quadrature_point) const;
3525 
3532  const std::vector<DerivativeForm<3, dim, spacedim>> &
3534 
3544  const Tensor<4, spacedim> &
3546  const unsigned int quadrature_point) const;
3547 
3554  const std::vector<Tensor<4, spacedim>> &
3556 
3566  jacobian_3rd_derivative(const unsigned int quadrature_point) const;
3567 
3574  const std::vector<DerivativeForm<4, dim, spacedim>> &
3576 
3586  const Tensor<5, spacedim> &
3588  const unsigned int quadrature_point) const;
3589 
3596  const std::vector<Tensor<5, spacedim>> &
3598 
3606  inverse_jacobian(const unsigned int quadrature_point) const;
3607 
3614  const std::vector<DerivativeForm<1, spacedim, dim>> &
3616 
3636  const Tensor<1, spacedim> &
3637  normal_vector(const unsigned int i) const;
3638 
3646  const std::vector<Tensor<1, spacedim>> &
3647  get_normal_vectors() const;
3648 
3652 
3664 
3675 
3687 
3688 
3699 
3703 
3708  const Mapping<dim, spacedim> &
3709  get_mapping() const;
3710 
3715  get_fe() const;
3716 
3720  UpdateFlags
3722 
3727  get_cell() const;
3728 
3735  get_cell_similarity() const;
3736 
3741  std::size_t
3742  memory_consumption() const;
3754  std::string,
3755  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3756  << "object for which this kind of information has not been computed. What "
3757  << "information these objects compute is determined by the update_* flags you "
3758  << "pass to the constructor. Here, the operation you are attempting requires "
3759  << "the <" << arg1
3760  << "> flag to be set, but it was apparently not specified "
3761  << "upon construction.");
3762 
3768  DeclExceptionMsg(ExcNotReinited,
3769  "FEValues object is not reinit'ed to any cell");
3770 
3778  ExcFEDontMatch,
3779  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3780  "to the DoFHandler that provided the cell iterator do not match.");
3786  DeclException1(ExcShapeFunctionNotPrimitive,
3787  int,
3788  << "The shape function with index " << arg1
3789  << " is not primitive, i.e. it is vector-valued and "
3790  << "has more than one non-zero vector component. This "
3791  << "function cannot be called for these shape functions. "
3792  << "Maybe you want to use the same function with the "
3793  << "_component suffix?");
3794 
3803  "The given FiniteElement is not a primitive element but the requested operation "
3804  "only works for those. See FiniteElement::is_primitive() for more information.");
3805 
3806 protected:
3814  {
3815  public:
3817  ExcNeedsDoFHandler,
3818  "You have previously called the FEValues::reinit() function with a "
3819  "cell iterator of type Triangulation<dim,spacedim>::cell_iterator. However, "
3820  "when you do this, you cannot call some functions in the FEValues "
3821  "class, such as the get_function_values/gradients/hessians/third_derivatives "
3822  "functions. If you need these functions, then you need to call "
3823  "FEValues::reinit() with an iterator type that allows to extract "
3824  "degrees of freedom, such as DoFHandler<dim,spacedim>::cell_iterator.");
3825 
3830 
3834  template <bool lda>
3837 
3842  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3843 
3847  bool
3848  is_initialized() const;
3849 
3856  operator typename Triangulation<dim, spacedim>::cell_iterator() const;
3857 
3863  n_dofs_for_dof_handler() const;
3864 
3869  template <typename VectorType>
3870  void
3871  get_interpolated_dof_values(
3872  const VectorType & in,
3874 
3879  void
3880  get_interpolated_dof_values(const IndexSet & in,
3881  Vector<IndexSet::value_type> &out) const;
3882 
3883  private:
3888  };
3889 
3896 
3904  boost::signals2::connection tria_listener_refinement;
3905 
3913  boost::signals2::connection tria_listener_mesh_transform;
3914 
3920  void
3921  invalidate_present_cell();
3922 
3932  void
3933  maybe_invalidate_previous_present_cell(
3934  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3935 
3941 
3947  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3949 
3956 
3957 
3965 
3971  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3973 
3979  spacedim>
3981 
3982 
3987 
3996  UpdateFlags
3997  compute_update_flags(const UpdateFlags update_flags) const;
3998 
4005 
4011  void
4012  check_cell_similarity(
4013  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
4014 
4015 private:
4020 
4021  // Make the view classes friends of this class, since they access internal
4022  // data.
4023  template <int, int>
4025  template <int, int>
4027  template <int, int, int>
4029  template <int, int, int>
4031 };
4032 
4033 
4034 
4044 template <int dim, int spacedim = dim>
4045 class FEValues : public FEValuesBase<dim, spacedim>
4046 {
4047 public:
4052  static constexpr unsigned int integral_dimension = dim;
4053 
4058  FEValues(const Mapping<dim, spacedim> & mapping,
4059  const FiniteElement<dim, spacedim> &fe,
4060  const Quadrature<dim> & quadrature,
4061  const UpdateFlags update_flags);
4062 
4069  FEValues(const Mapping<dim, spacedim> & mapping,
4070  const FiniteElement<dim, spacedim> &fe,
4071  const hp::QCollection<dim> & quadrature,
4072  const UpdateFlags update_flags);
4073 
4080  const Quadrature<dim> & quadrature,
4081  const UpdateFlags update_flags);
4082 
4090  const hp::QCollection<dim> & quadrature,
4091  const UpdateFlags update_flags);
4092 
4099  template <bool level_dof_access>
4100  void
4103 
4117  void
4119 
4124  const Quadrature<dim> &
4126 
4131  std::size_t
4132  memory_consumption() const;
4133 
4148  const FEValues<dim, spacedim> &
4150 
4151 private:
4156 
4160  void
4161  initialize(const UpdateFlags update_flags);
4162 
4169  void
4170  do_reinit();
4171 };
4172 
4173 
4183 template <int dim, int spacedim = dim>
4184 class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
4185 {
4186 public:
4191  static constexpr unsigned int integral_dimension = dim - 1;
4192 
4204  FEFaceValuesBase(const unsigned int dofs_per_cell,
4205  const UpdateFlags update_flags,
4206  const Mapping<dim, spacedim> & mapping,
4207  const FiniteElement<dim, spacedim> &fe,
4208  const Quadrature<dim - 1> & quadrature);
4209 
4216  FEFaceValuesBase(const unsigned int dofs_per_cell,
4217  const UpdateFlags update_flags,
4218  const Mapping<dim, spacedim> & mapping,
4219  const FiniteElement<dim, spacedim> &fe,
4220  const hp::QCollection<dim - 1> & quadrature);
4221 
4229  const Tensor<1, spacedim> &
4230  boundary_form(const unsigned int i) const;
4231 
4238  const std::vector<Tensor<1, spacedim>> &
4239  get_boundary_forms() const;
4240 
4245  unsigned int
4247 
4252  unsigned int
4254 
4259  const Quadrature<dim - 1> &
4261 
4266  std::size_t
4267  memory_consumption() const;
4268 
4269 protected:
4274  unsigned int present_face_no;
4275 
4280  unsigned int present_face_index;
4281 
4285  const hp::QCollection<dim - 1> quadrature;
4286 };
4287 
4288 
4289 
4303 template <int dim, int spacedim = dim>
4304 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
4305 {
4306 public:
4311  static constexpr unsigned int dimension = dim;
4312 
4313  static constexpr unsigned int space_dimension = spacedim;
4314 
4319  static constexpr unsigned int integral_dimension = dim - 1;
4320 
4325  FEFaceValues(const Mapping<dim, spacedim> & mapping,
4326  const FiniteElement<dim, spacedim> &fe,
4327  const Quadrature<dim - 1> & quadrature,
4328  const UpdateFlags update_flags);
4329 
4336  FEFaceValues(const Mapping<dim, spacedim> & mapping,
4337  const FiniteElement<dim, spacedim> &fe,
4338  const hp::QCollection<dim - 1> & quadrature,
4339  const UpdateFlags update_flags);
4340 
4347  const Quadrature<dim - 1> & quadrature,
4348  const UpdateFlags update_flags);
4349 
4357  const hp::QCollection<dim - 1> & quadrature,
4358  const UpdateFlags update_flags);
4359 
4364  template <bool level_dof_access>
4365  void
4368  const unsigned int face_no);
4369 
4376  template <bool level_dof_access>
4377  void
4380  const typename Triangulation<dim, spacedim>::face_iterator & face);
4381 
4395  void
4397  const unsigned int face_no);
4398 
4399  /*
4400  * Reinitialize the gradients, Jacobi determinants, etc for the given face
4401  * on a given cell of type "iterator into a Triangulation object", and the
4402  * given finite element. Since iterators into a triangulation alone only
4403  * convey information about the geometry of a cell, but not about degrees of
4404  * freedom possibly associated with this cell, you will not be able to call
4405  * some functions of this class if they need information about degrees of
4406  * freedom. These functions are, above all, the
4407  * <tt>get_function_value/gradients/hessians/third_derivatives</tt>
4408  * functions. If you want to call these functions, you have to call the @p
4409  * reinit variants that take iterators into DoFHandler or other DoF handler
4410  * type objects.
4411  *
4412  * @note @p face must be one of @p cell's face iterators.
4413  */
4414  void
4416  const typename Triangulation<dim, spacedim>::face_iterator &face);
4417 
4434 
4435 private:
4439  void
4440  initialize(const UpdateFlags update_flags);
4441 
4448  void
4449  do_reinit(const unsigned int face_no);
4450 };
4451 
4452 
4469 template <int dim, int spacedim = dim>
4470 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
4471 {
4472 public:
4476  static constexpr unsigned int dimension = dim;
4477 
4481  static constexpr unsigned int space_dimension = spacedim;
4482 
4487  static constexpr unsigned int integral_dimension = dim - 1;
4488 
4493  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
4494  const FiniteElement<dim, spacedim> &fe,
4495  const Quadrature<dim - 1> & face_quadrature,
4496  const UpdateFlags update_flags);
4497 
4504  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
4505  const FiniteElement<dim, spacedim> &fe,
4506  const hp::QCollection<dim - 1> & face_quadrature,
4507  const UpdateFlags update_flags);
4508 
4515  const Quadrature<dim - 1> & face_quadrature,
4516  const UpdateFlags update_flags);
4517 
4525  const hp::QCollection<dim - 1> & face_quadrature,
4526  const UpdateFlags update_flags);
4527 
4534  template <bool level_dof_access>
4535  void
4538  const unsigned int face_no,
4539  const unsigned int subface_no);
4540 
4545  template <bool level_dof_access>
4546  void
4549  const typename Triangulation<dim, spacedim>::face_iterator & face,
4550  const typename Triangulation<dim, spacedim>::face_iterator &subface);
4551 
4565  void
4567  const unsigned int face_no,
4568  const unsigned int subface_no);
4569 
4589  void
4591  const typename Triangulation<dim, spacedim>::face_iterator &face,
4592  const typename Triangulation<dim, spacedim>::face_iterator &subface);
4593 
4610 
4616  DeclException0(ExcReinitCalledWithBoundaryFace);
4617 
4623  DeclException0(ExcFaceHasNoSubfaces);
4624 
4625 private:
4629  void
4630  initialize(const UpdateFlags update_flags);
4631 
4638  void
4639  do_reinit(const unsigned int face_no, const unsigned int subface_no);
4640 };
4641 
4642 
4643 #ifndef DOXYGEN
4644 
4645 
4646 /*------------------------ Inline functions: namespace FEValuesViews --------*/
4647 
4648 namespace FEValuesViews
4649 {
4650  template <int dim, int spacedim>
4651  inline typename Scalar<dim, spacedim>::value_type
4652  Scalar<dim, spacedim>::value(const unsigned int shape_function,
4653  const unsigned int q_point) const
4654  {
4655  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4656  Assert(
4657  fe_values->update_flags & update_values,
4659  "update_values"))));
4660 
4661  // an adaptation of the FEValuesBase::shape_value_component function
4662  // except that here we know the component as fixed and we have
4663  // pre-computed and cached a bunch of information. See the comments there.
4664  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4665  return fe_values->finite_element_output.shape_values(
4666  shape_function_data[shape_function].row_index, q_point);
4667  else
4668  return 0;
4669  }
4670 
4671 
4672 
4673  template <int dim, int spacedim>
4674  inline typename Scalar<dim, spacedim>::gradient_type
4675  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
4676  const unsigned int q_point) const
4677  {
4678  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4679  Assert(fe_values->update_flags & update_gradients,
4681  "update_gradients")));
4682 
4683  // an adaptation of the FEValuesBase::shape_grad_component
4684  // function except that here we know the component as fixed and we have
4685  // pre-computed and cached a bunch of information. See the comments there.
4686  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4687  return fe_values->finite_element_output
4688  .shape_gradients[shape_function_data[shape_function].row_index]
4689  [q_point];
4690  else
4691  return gradient_type();
4692  }
4693 
4694 
4695 
4696  template <int dim, int spacedim>
4697  inline typename Scalar<dim, spacedim>::hessian_type
4698  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
4699  const unsigned int q_point) const
4700  {
4701  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4702  Assert(fe_values->update_flags & update_hessians,
4704  "update_hessians")));
4705 
4706  // an adaptation of the FEValuesBase::shape_hessian_component
4707  // function except that here we know the component as fixed and we have
4708  // pre-computed and cached a bunch of information. See the comments there.
4709  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4710  return fe_values->finite_element_output
4711  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
4712  else
4713  return hessian_type();
4714  }
4715 
4716 
4717 
4718  template <int dim, int spacedim>
4720  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
4721  const unsigned int q_point) const
4722  {
4723  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4724  Assert(fe_values->update_flags & update_3rd_derivatives,
4726  "update_3rd_derivatives")));
4727 
4728  // an adaptation of the FEValuesBase::shape_3rdderivative_component
4729  // function except that here we know the component as fixed and we have
4730  // pre-computed and cached a bunch of information. See the comments there.
4731  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4732  return fe_values->finite_element_output
4733  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
4734  [q_point];
4735  else
4736  return third_derivative_type();
4737  }
4738 
4739 
4740 
4741  template <int dim, int spacedim>
4742  inline typename Vector<dim, spacedim>::value_type
4743  Vector<dim, spacedim>::value(const unsigned int shape_function,
4744  const unsigned int q_point) const
4745  {
4746  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4747  Assert(fe_values->update_flags & update_values,
4749  "update_values")));
4750 
4751  // same as for the scalar case except that we have one more index
4752  const int snc =
4753  shape_function_data[shape_function].single_nonzero_component;
4754  if (snc == -2)
4755  return value_type();
4756  else if (snc != -1)
4757  {
4758  value_type return_value;
4759  return_value[shape_function_data[shape_function]
4760  .single_nonzero_component_index] =
4761  fe_values->finite_element_output.shape_values(snc, q_point);
4762  return return_value;
4763  }
4764  else
4765  {
4766  value_type return_value;
4767  for (unsigned int d = 0; d < dim; ++d)
4768  if (shape_function_data[shape_function]
4769  .is_nonzero_shape_function_component[d])
4770  return_value[d] = fe_values->finite_element_output.shape_values(
4771  shape_function_data[shape_function].row_index[d], q_point);
4772 
4773  return return_value;
4774  }
4775  }
4776 
4777 
4778 
4779  template <int dim, int spacedim>
4780  inline typename Vector<dim, spacedim>::gradient_type
4781  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4782  const unsigned int q_point) const
4783  {
4784  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4785  Assert(fe_values->update_flags & update_gradients,
4787  "update_gradients")));
4788 
4789  // same as for the scalar case except that we have one more index
4790  const int snc =
4791  shape_function_data[shape_function].single_nonzero_component;
4792  if (snc == -2)
4793  return gradient_type();
4794  else if (snc != -1)
4795  {
4796  gradient_type return_value;
4797  return_value[shape_function_data[shape_function]
4798  .single_nonzero_component_index] =
4799  fe_values->finite_element_output.shape_gradients[snc][q_point];
4800  return return_value;
4801  }
4802  else
4803  {
4804  gradient_type return_value;
4805  for (unsigned int d = 0; d < dim; ++d)
4806  if (shape_function_data[shape_function]
4807  .is_nonzero_shape_function_component[d])
4808  return_value[d] =
4809  fe_values->finite_element_output.shape_gradients
4810  [shape_function_data[shape_function].row_index[d]][q_point];
4811 
4812  return return_value;
4813  }
4814  }
4815 
4816 
4817 
4818  template <int dim, int spacedim>
4820  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4821  const unsigned int q_point) const
4822  {
4823  // this function works like in the case above
4824  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4825  Assert(fe_values->update_flags & update_gradients,
4827  "update_gradients")));
4828 
4829  // same as for the scalar case except that we have one more index
4830  const int snc =
4831  shape_function_data[shape_function].single_nonzero_component;
4832  if (snc == -2)
4833  return divergence_type();
4834  else if (snc != -1)
4835  return fe_values->finite_element_output
4836  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4837  .single_nonzero_component_index];
4838  else
4839  {
4840  divergence_type return_value = 0;
4841  for (unsigned int d = 0; d < dim; ++d)
4842  if (shape_function_data[shape_function]
4843  .is_nonzero_shape_function_component[d])
4844  return_value +=
4845  fe_values->finite_element_output.shape_gradients
4846  [shape_function_data[shape_function].row_index[d]][q_point][d];
4847 
4848  return return_value;
4849  }
4850  }
4851 
4852 
4853 
4854  template <int dim, int spacedim>
4855  inline typename Vector<dim, spacedim>::curl_type
4856  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4857  const unsigned int q_point) const
4858  {
4859  // this function works like in the case above
4860 
4861  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4862  Assert(fe_values->update_flags & update_gradients,
4864  "update_gradients")));
4865  // same as for the scalar case except that we have one more index
4866  const int snc =
4867  shape_function_data[shape_function].single_nonzero_component;
4868 
4869  if (snc == -2)
4870  return curl_type();
4871 
4872  else
4873  switch (dim)
4874  {
4875  case 1:
4876  {
4877  Assert(false,
4878  ExcMessage(
4879  "Computing the curl in 1d is not a useful operation"));
4880  return curl_type();
4881  }
4882 
4883  case 2:
4884  {
4885  if (snc != -1)
4886  {
4887  curl_type return_value;
4888 
4889  // the single nonzero component can only be zero or one in 2d
4890  if (shape_function_data[shape_function]
4891  .single_nonzero_component_index == 0)
4892  return_value[0] =
4893  -1.0 * fe_values->finite_element_output
4894  .shape_gradients[snc][q_point][1];
4895  else
4896  return_value[0] = fe_values->finite_element_output
4897  .shape_gradients[snc][q_point][0];
4898 
4899  return return_value;
4900  }
4901 
4902  else
4903  {
4904  curl_type return_value;
4905 
4906  return_value[0] = 0.0;
4907 
4908  if (shape_function_data[shape_function]
4909  .is_nonzero_shape_function_component[0])
4910  return_value[0] -=
4911  fe_values->finite_element_output
4912  .shape_gradients[shape_function_data[shape_function]
4913  .row_index[0]][q_point][1];
4914 
4915  if (shape_function_data[shape_function]
4916  .is_nonzero_shape_function_component[1])
4917  return_value[0] +=
4918  fe_values->finite_element_output
4919  .shape_gradients[shape_function_data[shape_function]
4920  .row_index[1]][q_point][0];
4921 
4922  return return_value;
4923  }
4924  }
4925 
4926  case 3:
4927  {
4928  if (snc != -1)
4929  {
4930  curl_type return_value;
4931 
4932  switch (shape_function_data[shape_function]
4933  .single_nonzero_component_index)
4934  {
4935  case 0:
4936  {
4937  return_value[0] = 0;
4938  return_value[1] = fe_values->finite_element_output
4939  .shape_gradients[snc][q_point][2];
4940  return_value[2] =
4941  -1.0 * fe_values->finite_element_output
4942  .shape_gradients[snc][q_point][1];
4943  return return_value;
4944  }
4945 
4946  case 1:
4947  {
4948  return_value[0] =
4949  -1.0 * fe_values->finite_element_output
4950  .shape_gradients[snc][q_point][2];
4951  return_value[1] = 0;
4952  return_value[2] = fe_values->finite_element_output
4953  .shape_gradients[snc][q_point][0];
4954  return return_value;
4955  }
4956 
4957  default:
4958  {
4959  return_value[0] = fe_values->finite_element_output
4960  .shape_gradients[snc][q_point][1];
4961  return_value[1] =
4962  -1.0 * fe_values->finite_element_output
4963  .shape_gradients[snc][q_point][0];
4964  return_value[2] = 0;
4965  return return_value;
4966  }
4967  }
4968  }
4969 
4970  else
4971  {
4972  curl_type return_value;
4973 
4974  for (unsigned int i = 0; i < dim; ++i)
4975  return_value[i] = 0.0;
4976 
4977  if (shape_function_data[shape_function]
4978  .is_nonzero_shape_function_component[0])
4979  {
4980  return_value[1] +=
4981  fe_values->finite_element_output
4982  .shape_gradients[shape_function_data[shape_function]
4983  .row_index[0]][q_point][2];
4984  return_value[2] -=
4985  fe_values->finite_element_output
4986  .shape_gradients[shape_function_data[shape_function]
4987  .row_index[0]][q_point][1];
4988  }
4989 
4990  if (shape_function_data[shape_function]
4991  .is_nonzero_shape_function_component[1])
4992  {
4993  return_value[0] -=
4994  fe_values->finite_element_output
4995  .shape_gradients[shape_function_data[shape_function]
4996  .row_index[1]][q_point][2];
4997  return_value[2] +=
4998  fe_values->finite_element_output
4999  .shape_gradients[shape_function_data[shape_function]
5000  .row_index[1]][q_point][0];
5001  }
5002 
5003  if (shape_function_data[shape_function]
5004  .is_nonzero_shape_function_component[2])
5005  {
5006  return_value[0] +=
5007  fe_values->finite_element_output
5008  .shape_gradients[shape_function_data[shape_function]
5009  .row_index[2]][q_point][1];
5010  return_value[1] -=
5011  fe_values->finite_element_output
5012  .shape_gradients[shape_function_data[shape_function]
5013  .row_index[2]][q_point][0];
5014  }
5015 
5016  return return_value;
5017  }
5018  }
5019  }
5020  // should not end up here
5021  Assert(false, ExcInternalError());
5022  return curl_type();
5023  }
5024 
5025 
5026 
5027  template <int dim, int spacedim>
5028  inline typename Vector<dim, spacedim>::hessian_type
5029  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
5030  const unsigned int q_point) const
5031  {
5032  // this function works like in the case above
5033  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5034  Assert(fe_values->update_flags & update_hessians,
5036  "update_hessians")));
5037 
5038  // same as for the scalar case except that we have one more index
5039  const int snc =
5040  shape_function_data[shape_function].single_nonzero_component;
5041  if (snc == -2)
5042  return hessian_type();
5043  else if (snc != -1)
5044  {
5045  hessian_type return_value;
5046  return_value[shape_function_data[shape_function]
5047  .single_nonzero_component_index] =
5048  fe_values->finite_element_output.shape_hessians[snc][q_point];
5049  return return_value;
5050  }
5051  else
5052  {
5053  hessian_type return_value;
5054  for (unsigned int d = 0; d < dim; ++d)
5055  if (shape_function_data[shape_function]
5056  .is_nonzero_shape_function_component[d])
5057  return_value[d] =
5058  fe_values->finite_element_output.shape_hessians
5059  [shape_function_data[shape_function].row_index[d]][q_point];
5060 
5061  return return_value;
5062  }
5063  }
5064 
5065 
5066 
5067  template <int dim, int spacedim>
5069  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
5070  const unsigned int q_point) const
5071  {
5072  // this function works like in the case above
5073  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5074  Assert(fe_values->update_flags & update_3rd_derivatives,
5076  "update_3rd_derivatives")));
5077 
5078  // same as for the scalar case except that we have one more index
5079  const int snc =
5080  shape_function_data[shape_function].single_nonzero_component;
5081  if (snc == -2)
5082  return third_derivative_type();
5083  else if (snc != -1)
5084  {
5085  third_derivative_type return_value;
5086  return_value[shape_function_data[shape_function]
5087  .single_nonzero_component_index] =
5088  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
5089  return return_value;
5090  }
5091  else
5092  {
5093  third_derivative_type return_value;
5094  for (unsigned int d = 0; d < dim; ++d)
5095  if (shape_function_data[shape_function]
5096  .is_nonzero_shape_function_component[d])
5097  return_value[d] =
5098  fe_values->finite_element_output.shape_3rd_derivatives
5099  [shape_function_data[shape_function].row_index[d]][q_point];
5100 
5101  return return_value;
5102  }
5103  }
5104 
5105 
5106 
5107  namespace internal
5108  {
5113  inline ::SymmetricTensor<2, 1>
5114  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
5115  {
5116  AssertIndexRange(n, 1);
5117  (void)n;
5118 
5119  return {{t[0]}};
5120  }
5121 
5122 
5123 
5124  inline ::SymmetricTensor<2, 2>
5125  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
5126  {
5127  switch (n)
5128  {
5129  case 0:
5130  {
5131  return {{t[0], 0, t[1] / 2}};
5132  }
5133  case 1:
5134  {
5135  return {{0, t[1], t[0] / 2}};
5136  }
5137  default:
5138  {
5139  AssertIndexRange(n, 2);
5140  return {};
5141  }
5142  }
5143  }
5144 
5145 
5146 
5147  inline ::SymmetricTensor<2, 3>
5148  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
5149  {
5150  switch (n)
5151  {
5152  case 0:
5153  {
5154  return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
5155  }
5156  case 1:
5157  {
5158  return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
5159  }
5160  case 2:
5161  {
5162  return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
5163  }
5164  default:
5165  {
5166  AssertIndexRange(n, 3);
5167  return {};
5168  }
5169  }
5170  }
5171  } // namespace internal
5172 
5173 
5174 
5175  template <int dim, int spacedim>
5177  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
5178  const unsigned int q_point) const
5179  {
5180  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5181  Assert(fe_values->update_flags & update_gradients,
5183  "update_gradients")));
5184 
5185  // same as for the scalar case except that we have one more index
5186  const int snc =
5187  shape_function_data[shape_function].single_nonzero_component;
5188  if (snc == -2)
5189  return symmetric_gradient_type();
5190  else if (snc != -1)
5191  return internal::symmetrize_single_row(
5192  shape_function_data[shape_function].single_nonzero_component_index,
5193  fe_values->finite_element_output.shape_gradients[snc][q_point]);
5194  else
5195  {
5196  gradient_type return_value;
5197  for (unsigned int d = 0; d < dim; ++d)
5198  if (shape_function_data[shape_function]
5199  .is_nonzero_shape_function_component[d])
5200  return_value[d] =
5201  fe_values->finite_element_output.shape_gradients
5202  [shape_function_data[shape_function].row_index[d]][q_point];
5203 
5204  return symmetrize(return_value);
5205  }
5206  }
5207 
5208 
5209 
5210  template <int dim, int spacedim>
5212  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
5213  const unsigned int q_point) const
5214  {
5215  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5216  Assert(fe_values->update_flags & update_values,
5218  "update_values")));
5219 
5220  // similar to the vector case where we have more then one index and we need
5221  // to convert between unrolled and component indexing for tensors
5222  const int snc =
5223  shape_function_data[shape_function].single_nonzero_component;
5224 
5225  if (snc == -2)
5226  {
5227  // shape function is zero for the selected components
5228  return value_type();
5229  }
5230  else if (snc != -1)
5231  {
5232  value_type return_value;
5233  const unsigned int comp =
5234  shape_function_data[shape_function].single_nonzero_component_index;
5235  return_value[value_type::unrolled_to_component_indices(comp)] =
5236  fe_values->finite_element_output.shape_values(snc, q_point);
5237  return return_value;
5238  }
5239  else
5240  {
5241  value_type return_value;
5242  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
5243  if (shape_function_data[shape_function]
5244  .is_nonzero_shape_function_component[d])
5245  return_value[value_type::unrolled_to_component_indices(d)] =
5246  fe_values->finite_element_output.shape_values(
5247  shape_function_data[shape_function].row_index[d], q_point);
5248  return return_value;
5249  }
5250  }
5251 
5252 
5253 
5254  template <int dim, int spacedim>
5257  const unsigned int shape_function,
5258  const unsigned int q_point) const
5259  {
5260  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5261  Assert(fe_values->update_flags & update_gradients,
5263  "update_gradients")));
5264 
5265  const int snc =
5266  shape_function_data[shape_function].single_nonzero_component;
5267 
5268  if (snc == -2)
5269  {
5270  // shape function is zero for the selected components
5271  return divergence_type();
5272  }
5273  else if (snc != -1)
5274  {
5275  // we have a single non-zero component when the symmetric tensor is
5276  // represented in unrolled form. this implies we potentially have
5277  // two non-zero components when represented in component form! we
5278  // will only have one non-zero entry if the non-zero component lies on
5279  // the diagonal of the tensor.
5280  //
5281  // the divergence of a second-order tensor is a first order tensor.
5282  //
5283  // assume the second-order tensor is A with components A_{ij}. then
5284  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
5285  // entries in the tensorial representation. define the
5286  // divergence as:
5287  // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
5288  // (which is incidentally also
5289  // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
5290  // In both cases, a sum is implied.
5291  //
5292  // Now, we know the nonzero component in unrolled form: it is indicated
5293  // by 'snc'. we can figure out which tensor components belong to this:
5294  const unsigned int comp =
5295  shape_function_data[shape_function].single_nonzero_component_index;
5296  const unsigned int ii =
5297  value_type::unrolled_to_component_indices(comp)[0];
5298  const unsigned int jj =
5299  value_type::unrolled_to_component_indices(comp)[1];
5300 
5301  // given the form of the divergence above, if ii=jj there is only a
5302  // single nonzero component of the full tensor and the gradient
5303  // equals
5304  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
5305  // all other entries of 'b' are zero
5306  //
5307  // on the other hand, if ii!=jj, then there are two nonzero entries in
5308  // the full tensor and
5309  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
5310  // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
5311  // again, all other entries of 'b' are zero
5312  const ::Tensor<1, spacedim> &phi_grad =
5313  fe_values->finite_element_output.shape_gradients[snc][q_point];
5314 
5315  divergence_type return_value;
5316  return_value[ii] = phi_grad[jj];
5317 
5318  if (ii != jj)
5319  return_value[jj] = phi_grad[ii];
5320 
5321  return return_value;
5322  }
5323  else
5324  {
5325  Assert(false, ExcNotImplemented());
5326  divergence_type return_value;
5327  return return_value;
5328  }
5329  }
5330 
5331 
5332 
5333  template <int dim, int spacedim>
5334  inline typename Tensor<2, dim, spacedim>::value_type
5335  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
5336  const unsigned int q_point) const
5337  {
5338  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5339  Assert(fe_values->update_flags & update_values,
5341  "update_values")));
5342 
5343  // similar to the vector case where we have more then one index and we need
5344  // to convert between unrolled and component indexing for tensors
5345  const int snc =
5346  shape_function_data[shape_function].single_nonzero_component;
5347 
5348  if (snc == -2)
5349  {
5350  // shape function is zero for the selected components
5351  return value_type();
5352  }
5353  else if (snc != -1)
5354  {
5355  value_type return_value;
5356  const unsigned int comp =
5357  shape_function_data[shape_function].single_nonzero_component_index;
5358  const TableIndices<2> indices =
5360  return_value[indices] =
5361  fe_values->finite_element_output.shape_values(snc, q_point);
5362  return return_value;
5363  }
5364  else
5365  {
5366  value_type return_value;
5367  for (unsigned int d = 0; d < dim * dim; ++d)
5368  if (shape_function_data[shape_function]
5369  .is_nonzero_shape_function_component[d])
5370  {
5371  const TableIndices<2> indices =
5373  return_value[indices] =
5374  fe_values->finite_element_output.shape_values(
5375  shape_function_data[shape_function].row_index[d], q_point);
5376  }
5377  return return_value;
5378  }
5379  }
5380 
5381 
5382 
5383  template <int dim, int spacedim>
5385  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
5386  const unsigned int q_point) const
5387  {
5388  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5389  Assert(fe_values->update_flags & update_gradients,
5391  "update_gradients")));
5392 
5393  const int snc =
5394  shape_function_data[shape_function].single_nonzero_component;
5395 
5396  if (snc == -2)
5397  {
5398  // shape function is zero for the selected components
5399  return divergence_type();
5400  }
5401  else if (snc != -1)
5402  {
5403  // we have a single non-zero component when the tensor is
5404  // represented in unrolled form.
5405  //
5406  // the divergence of a second-order tensor is a first order tensor.
5407  //
5408  // assume the second-order tensor is A with components A_{ij},
5409  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
5410  //
5411  // Now, we know the nonzero component in unrolled form: it is indicated
5412  // by 'snc'. we can figure out which tensor components belong to this:
5413  const unsigned int comp =
5414  shape_function_data[shape_function].single_nonzero_component_index;
5415  const TableIndices<2> indices =
5417  const unsigned int ii = indices[0];
5418  const unsigned int jj = indices[1];
5419 
5420  const ::Tensor<1, spacedim> &phi_grad =
5421  fe_values->finite_element_output.shape_gradients[snc][q_point];
5422 
5423  divergence_type return_value;
5424  // note that we contract \nabla from the right
5425  return_value[ii] = phi_grad[jj];
5426 
5427  return return_value;
5428  }
5429  else
5430  {
5431  Assert(false, ExcNotImplemented());
5432  divergence_type return_value;
5433  return return_value;
5434  }
5435  }
5436 
5437 
5438 
5439  template <int dim, int spacedim>
5441  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
5442  const unsigned int q_point) const
5443  {
5444  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5445  Assert(fe_values->update_flags & update_gradients,
5447  "update_gradients")));
5448 
5449  const int snc =
5450  shape_function_data[shape_function].single_nonzero_component;
5451 
5452  if (snc == -2)
5453  {
5454  // shape function is zero for the selected components
5455  return gradient_type();
5456  }
5457  else if (snc != -1)
5458  {
5459  // we have a single non-zero component when the tensor is
5460  // represented in unrolled form.
5461  //
5462  // the gradient of a second-order tensor is a third order tensor.
5463  //
5464  // assume the second-order tensor is A with components A_{ij},
5465  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
5466  //
5467  // Now, we know the nonzero component in unrolled form: it is indicated
5468  // by 'snc'. we can figure out which tensor components belong to this:
5469  const unsigned int comp =
5470  shape_function_data[shape_function].single_nonzero_component_index;
5471  const TableIndices<2> indices =
5473  const unsigned int ii = indices[0];
5474  const unsigned int jj = indices[1];
5475 
5476  const ::Tensor<1, spacedim> &phi_grad =
5477  fe_values->finite_element_output.shape_gradients[snc][q_point];
5478 
5479  gradient_type return_value;
5480  return_value[ii][jj] = phi_grad;
5481 
5482  return return_value;
5483  }
5484  else
5485  {
5486  Assert(false, ExcNotImplemented());
5487  gradient_type return_value;
5488  return return_value;
5489  }
5490  }
5491 
5492 } // namespace FEValuesViews
5493 
5494 
5495 
5496 /*---------------------- Inline functions: FEValuesBase ---------------------*/
5497 
5498 
5499 
5500 template <int dim, int spacedim>
5501 template <bool lda>
5505  : initialized(true)
5506  , cell(cell)
5507  , dof_handler(&cell->get_dof_handler())
5508  , level_dof_access(lda)
5509 {}
5510 
5511 
5512 
5513 template <int dim, int spacedim>
5516  const FEValuesExtractors::Scalar &scalar) const
5517 {
5518  AssertIndexRange(scalar.component, fe_values_views_cache.scalars.size());
5519 
5520  return fe_values_views_cache.scalars[scalar.component];
5521 }
5522 
5523 
5524 
5525 template <int dim, int spacedim>
5528  const FEValuesExtractors::Vector &vector) const
5529 {
5531  fe_values_views_cache.vectors.size());
5532 
5533  return fe_values_views_cache.vectors[vector.first_vector_component];
5534 }
5535 
5536 
5537 
5538 template <int dim, int spacedim>
5541  const FEValuesExtractors::SymmetricTensor<2> &tensor) const
5542 {
5543  Assert(
5544  tensor.first_tensor_component <
5545  fe_values_views_cache.symmetric_second_order_tensors.size(),
5547  0,
5548  fe_values_views_cache.symmetric_second_order_tensors.size()));
5549 
5550  return fe_values_views_cache
5551  .symmetric_second_order_tensors[tensor.first_tensor_component];
5552 }
5553 
5554 
5555 
5556 template <int dim, int spacedim>
5559  const FEValuesExtractors::Tensor<2> &tensor) const
5560 {
5562  fe_values_views_cache.second_order_tensors.size());
5563 
5564  return fe_values_views_cache
5565  .second_order_tensors[tensor.first_tensor_component];
5566 }
5567 
5568 
5569 
5570 template <int dim, int spacedim>
5571 inline const double &
5572 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
5573  const unsigned int j) const
5574 {
5575  AssertIndexRange(i, fe->n_dofs_per_cell());
5576  Assert(this->update_flags & update_values,
5577  ExcAccessToUninitializedField("update_values"));
5578  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5579  Assert(present_cell.is_initialized(), ExcNotReinited());
5580  // if the entire FE is primitive,
5581  // then we can take a short-cut:
5582  if (fe->is_primitive())
5583  return this->finite_element_output.shape_values(i, j);
5584  else
5585  {
5586  // otherwise, use the mapping
5587  // between shape function
5588  // numbers and rows. note that
5589  // by the assertions above, we
5590  // know that this particular
5591  // shape function is primitive,
5592  // so we can call
5593  // system_to_component_index
5594  const unsigned int row =
5595  this->finite_element_output
5596  .shape_function_to_row_table[i * fe->n_components() +
5597  fe->system_to_component_index(i).first];
5598  return this->finite_element_output.shape_values(row, j);
5599  }
5600 }
5601 
5602 
5603 
5604 template <int dim, int spacedim>
5605 inline double
5607  const unsigned int i,
5608  const unsigned int j,
5609  const unsigned int component) const
5610 {
5611  AssertIndexRange(i, fe->n_dofs_per_cell());
5612  Assert(this->update_flags & update_values,
5613  ExcAccessToUninitializedField("update_values"));
5614  AssertIndexRange(component, fe->n_components());
5615  Assert(present_cell.is_initialized(), ExcNotReinited());
5616 
5617  // check whether the shape function
5618  // is non-zero at all within
5619  // this component:
5620  if (fe->get_nonzero_components(i)[component] == false)
5621  return 0;
5622 
5623  // look up the right row in the
5624  // table and take the data from
5625  // there
5626  const unsigned int row =
5627  this->finite_element_output
5628  .shape_function_to_row_table[i * fe->n_components() + component];
5629  return this->finite_element_output.shape_values(row, j);
5630 }
5631 
5632 
5633 
5634 template <int dim, int spacedim>
5635 inline const Tensor<1, spacedim> &
5636 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
5637  const unsigned int j) const
5638 {
5639  AssertIndexRange(i, fe->n_dofs_per_cell());
5640  Assert(this->update_flags & update_gradients,
5641  ExcAccessToUninitializedField("update_gradients"));
5642  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5643  Assert(present_cell.is_initialized(), ExcNotReinited());
5644  // if the entire FE is primitive,
5645  // then we can take a short-cut:
5646  if (fe->is_primitive())
5647  return this->finite_element_output.shape_gradients[i][j];
5648  else
5649  {
5650  // otherwise, use the mapping
5651  // between shape function
5652  // numbers and rows. note that
5653  // by the assertions above, we
5654  // know that this particular
5655  // shape function is primitive,
5656  // so we can call
5657  // system_to_component_index
5658  const unsigned int row =
5659  this->finite_element_output
5660  .shape_function_to_row_table[i * fe->n_components() +
5661  fe->system_to_component_index(i).first];
5662  return this->finite_element_output.shape_gradients[row][j];
5663  }
5664 }
5665 
5666 
5667 
5668 template <int dim, int spacedim>
5669 inline Tensor<1, spacedim>
5671  const unsigned int i,
5672  const unsigned int j,
5673  const unsigned int component) const
5674 {
5675  AssertIndexRange(i, fe->n_dofs_per_cell());
5676  Assert(this->update_flags & update_gradients,
5677  ExcAccessToUninitializedField("update_gradients"));
5678  AssertIndexRange(component, fe->n_components());
5679  Assert(present_cell.is_initialized(), ExcNotReinited());
5680  // check whether the shape function
5681  // is non-zero at all within
5682  // this component:
5683  if (fe->get_nonzero_components(i)[component] == false)
5684  return Tensor<1, spacedim>();
5685 
5686  // look up the right row in the
5687  // table and take the data from
5688  // there
5689  const unsigned int row =
5690  this->finite_element_output
5691  .shape_function_to_row_table[i * fe->n_components() + component];
5692  return this->finite_element_output.shape_gradients[row][j];
5693 }
5694 
5695 
5696 
5697 template <int dim, int spacedim>
5698 inline const Tensor<2, spacedim> &
5699 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
5700  const unsigned int j) const
5701 {
5702  AssertIndexRange(i, fe->n_dofs_per_cell());
5703  Assert(this->update_flags & update_hessians,
5704  ExcAccessToUninitializedField("update_hessians"));
5705  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5706  Assert(present_cell.is_initialized(), ExcNotReinited());
5707  // if the entire FE is primitive,
5708  // then we can take a short-cut:
5709  if (fe->is_primitive())
5710  return this->finite_element_output.shape_hessians[i][j];
5711  else
5712  {
5713  // otherwise, use the mapping
5714  // between shape function
5715  // numbers and rows. note that
5716  // by the assertions above, we
5717  // know that this particular
5718  // shape function is primitive,
5719  // so we can call
5720  // system_to_component_index
5721  const unsigned int row =
5722  this->finite_element_output
5723  .shape_function_to_row_table[i * fe->n_components() +
5724  fe->system_to_component_index(i).first];
5725  return this->finite_element_output.shape_hessians[row][j];
5726  }
5727 }
5728 
5729 
5730 
5731 template <int dim, int spacedim>
5732 inline Tensor<2, spacedim>
5734  const unsigned int i,
5735  const unsigned int j,
5736  const unsigned int component) const
5737 {
5738  AssertIndexRange(i, fe->n_dofs_per_cell());
5739  Assert(this->update_flags & update_hessians,
5740  ExcAccessToUninitializedField("update_hessians"));
5741  AssertIndexRange(component, fe->n_components());
5742  Assert(present_cell.is_initialized(), ExcNotReinited());
5743  // check whether the shape function
5744  // is non-zero at all within
5745  // this component:
5746  if (fe->get_nonzero_components(i)[component] == false)
5747  return Tensor<2, spacedim>();
5748 
5749  // look up the right row in the
5750  // table and take the data from
5751  // there
5752  const unsigned int row =
5753  this->finite_element_output
5754  .shape_function_to_row_table[i * fe->n_components() + component];
5755  return this->finite_element_output.shape_hessians[row][j];
5756 }
5757 
5758 
5759 
5760 template <int dim, int spacedim>
5761 inline const Tensor<3, spacedim> &
5763  const unsigned int j) const
5764 {
5765  AssertIndexRange(i, fe->n_dofs_per_cell());
5766  Assert(this->update_flags & update_3rd_derivatives,
5767  ExcAccessToUninitializedField("update_3rd_derivatives"));
5768  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5769  Assert(present_cell.is_initialized(), ExcNotReinited());
5770  // if the entire FE is primitive,
5771  // then we can take a short-cut:
5772  if (fe->is_primitive())
5773  return this->finite_element_output.shape_3rd_derivatives[i][j];
5774  else
5775  {
5776  // otherwise, use the mapping
5777  // between shape function
5778  // numbers and rows. note that
5779  // by the assertions above, we
5780  // know that this particular
5781  // shape function is primitive,
5782  // so we can call
5783  // system_to_component_index
5784  const unsigned int row =
5785  this->finite_element_output
5786  .shape_function_to_row_table[i * fe->n_components() +
5787  fe->system_to_component_index(i).first];
5788  return this->finite_element_output.shape_3rd_derivatives[row][j];
5789  }
5790 }
5791 
5792 
5793 
5794 template <int dim, int spacedim>
5795 inline Tensor<3, spacedim>
5797  const unsigned int i,
5798  const unsigned int j,
5799  const unsigned int component) const
5800 {
5801  AssertIndexRange(i, fe->n_dofs_per_cell());
5802  Assert(this->update_flags & update_3rd_derivatives,
5803  ExcAccessToUninitializedField("update_3rd_derivatives"));
5804  AssertIndexRange(component, fe->n_components());
5805  Assert(present_cell.is_initialized(), ExcNotReinited());
5806  // check whether the shape function
5807  // is non-zero at all within
5808  // this component:
5809  if (fe->get_nonzero_components(i)[component] == false)
5810  return Tensor<3, spacedim>();
5811 
5812  // look up the right row in the
5813  // table and take the data from
5814  // there
5815  const unsigned int row =
5816  this->finite_element_output
5817  .shape_function_to_row_table[i * fe->n_components() + component];
5818  return this->finite_element_output.shape_3rd_derivatives[row][j];
5819 }
5820 
5821 
5822 
5823 template <int dim, int spacedim>
5824 inline const FiniteElement<dim, spacedim> &
5826 {
5827  return *fe;
5828 }
5829 
5830 
5831 
5832 template <int dim, int spacedim>
5833 inline const Mapping<dim, spacedim> &
5835 {
5836  return *mapping;
5837 }
5838 
5839 
5840 
5841 template <int dim, int spacedim>
5842 inline UpdateFlags
5844 {
5845  return this->update_flags;
5846 }
5847 
5848 
5849 
5850 template <int dim, int spacedim>
5851 inline const std::vector<Point<spacedim>> &
5853 {
5854  Assert(this->update_flags & update_quadrature_points,
5855  ExcAccessToUninitializedField("update_quadrature_points"));
5856  Assert(present_cell.is_initialized(), ExcNotReinited());
5857  return this->mapping_output.quadrature_points;
5858 }
5859 
5860 
5861 
5862 template <int dim, int spacedim>
5863 inline const std::vector<double> &
5865 {
5866  Assert(this->update_flags & update_JxW_values,
5867  ExcAccessToUninitializedField("update_JxW_values"));
5868  Assert(present_cell.is_initialized(), ExcNotReinited());
5869  return this->mapping_output.JxW_values;
5870 }
5871 
5872 
5873 
5874 template <int dim, int spacedim>
5875 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5877 {
5878  Assert(this->update_flags & update_jacobians,
5879  ExcAccessToUninitializedField("update_jacobians"));
5880  Assert(present_cell.is_initialized(), ExcNotReinited());
5881  return this->mapping_output.jacobians;
5882 }
5883 
5884 
5885 
5886 template <int dim, int spacedim>
5887 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5889 {
5890  Assert(this->update_flags & update_jacobian_grads,
5891  ExcAccessToUninitializedField("update_jacobians_grads"));
5892  Assert(present_cell.is_initialized(), ExcNotReinited());
5893  return this->mapping_output.jacobian_grads;
5894 }
5895 
5896 
5897 
5898 template <int dim, int spacedim>
5899 inline const Tensor<3, spacedim> &
5901  const unsigned int i) const
5902 {
5903  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5904  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5905  Assert(present_cell.is_initialized(), ExcNotReinited());
5906  return this->mapping_output.jacobian_pushed_forward_grads[i];
5907 }
5908 
5909 
5910 
5911 template <int dim, int spacedim>
5912 inline const std::vector<Tensor<3, spacedim>> &
5914 {
5915  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5916  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5917  Assert(present_cell.is_initialized(), ExcNotReinited());
5918  return this->mapping_output.jacobian_pushed_forward_grads;
5919 }
5920 
5921 
5922 
5923 template <int dim, int spacedim>
5924 inline const DerivativeForm<3, dim, spacedim> &
5925 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5926 {
5927  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5928  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5929  Assert(present_cell.is_initialized(), ExcNotReinited());
5930  return this->mapping_output.jacobian_2nd_derivatives[i];
5931 }
5932 
5933 
5934 
5935 template <int dim, int spacedim>
5936 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5938 {
5939  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5940  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5941  Assert(present_cell.is_initialized(), ExcNotReinited());
5942  return this->mapping_output.jacobian_2nd_derivatives;
5943 }
5944 
5945 
5946 
5947 template <int dim, int spacedim>
5948 inline const Tensor<4, spacedim> &
5950  const unsigned int i) const
5951 {
5954  "update_jacobian_pushed_forward_2nd_derivatives"));
5955  Assert(present_cell.is_initialized(), ExcNotReinited());
5956  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5957 }
5958 
5959 
5960 
5961 template <int dim, int spacedim>
5962 inline const std::vector<Tensor<4, spacedim>> &
5964 {
5967  "update_jacobian_pushed_forward_2nd_derivatives"));
5968  Assert(present_cell.is_initialized(), ExcNotReinited());
5969  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5970 }
5971 
5972 
5973 
5974 template <int dim, int spacedim>
5975 inline const DerivativeForm<4, dim, spacedim> &
5976 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5977 {
5978  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5979  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5980  Assert(present_cell.is_initialized(), ExcNotReinited());
5981  return this->mapping_output.jacobian_3rd_derivatives[i];
5982 }
5983 
5984 
5985 
5986 template <int dim, int spacedim>
5987 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5989 {
5990  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5991  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5992  Assert(present_cell.is_initialized(), ExcNotReinited());
5993  return this->mapping_output.jacobian_3rd_derivatives;
5994 }
5995 
5996 
5997 
5998 template <int dim, int spacedim>
5999 inline const Tensor<5, spacedim> &
6001  const unsigned int i) const
6002 {
6005  "update_jacobian_pushed_forward_3rd_derivatives"));
6006  Assert(present_cell.is_initialized(), ExcNotReinited());
6007  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
6008 }
6009 
6010 
6011 
6012 template <int dim, int spacedim>
6013 inline const std::vector<Tensor<5, spacedim>> &
6015 {
6018  "update_jacobian_pushed_forward_3rd_derivatives"));
6019  Assert(present_cell.is_initialized(), ExcNotReinited());
6020  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
6021 }
6022 
6023 
6024 
6025 template <int dim, int spacedim>
6026 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
6028 {
6029  Assert(this->update_flags & update_inverse_jacobians,
6030  ExcAccessToUninitializedField("update_inverse_jacobians"));
6031  Assert(present_cell.is_initialized(), ExcNotReinited());
6032  return this->mapping_output.inverse_jacobians;
6033 }
6034 
6035 
6036 
6037 template <int dim, int spacedim>
6040 {
6041  return {0U, dofs_per_cell};
6042 }
6043 
6044 
6045 
6046 template <int dim, int spacedim>
6049  const unsigned int start_dof_index) const
6050 {
6051  Assert(start_dof_index <= dofs_per_cell,
6052  ExcIndexRange(start_dof_index, 0, dofs_per_cell + 1));
6053  return {start_dof_index, dofs_per_cell};
6054 }
6055 
6056 
6057 
6058 template <int dim, int spacedim>
6061  const unsigned int end_dof_index) const
6062 {
6063  Assert(end_dof_index < dofs_per_cell,
6064  ExcIndexRange(end_dof_index, 0, dofs_per_cell));
6065  return {0U, end_dof_index + 1};
6066 }
6067 
6068 
6069 
6070 template <int dim, int spacedim>
6073 {
6074  return {0U, n_quadrature_points};
6075 }
6076 
6077 
6078 
6079 template <int dim, int spacedim>
6080 inline const Point<spacedim> &
6081 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
6082 {
6083  Assert(this->update_flags & update_quadrature_points,
6084  ExcAccessToUninitializedField("update_quadrature_points"));
6085  AssertIndexRange(i, this->mapping_output.quadrature_points.size());
6086  Assert(present_cell.is_initialized(), ExcNotReinited());
6087 
6088  return this->mapping_output.quadrature_points[i];
6089 }
6090 
6091 
6092 
6093 template <int dim, int spacedim>
6094 inline double
6095 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
6096 {
6097  Assert(this->update_flags & update_JxW_values,
6098  ExcAccessToUninitializedField("update_JxW_values"));
6099  AssertIndexRange(i, this->mapping_output.JxW_values.size());
6100  Assert(present_cell.is_initialized(), ExcNotReinited());
6101 
6102  return this->mapping_output.JxW_values[i];
6103 }
6104 
6105 
6106 
6107 template <int dim, int spacedim>
6108 inline const DerivativeForm<1, dim, spacedim> &
6109 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
6110 {
6111  Assert(this->update_flags & update_jacobians,
6112  ExcAccessToUninitializedField("update_jacobians"));
6113  AssertIndexRange(i, this->mapping_output.jacobians.size());
6114  Assert(present_cell.is_initialized(), ExcNotReinited());
6115 
6116  return this->mapping_output.jacobians[i];
6117 }
6118 
6119 
6120 
6121 template <int dim, int spacedim>
6122 inline const DerivativeForm<2, dim, spacedim> &
6123 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
6124 {
6125  Assert(this->update_flags & update_jacobian_grads,
6126  ExcAccessToUninitializedField("update_jacobians_grads"));
6127  AssertIndexRange(i, this->mapping_output.jacobian_grads.size());
6128  Assert(present_cell.is_initialized(), ExcNotReinited());
6129 
6130  return this->mapping_output.jacobian_grads[i];
6131 }
6132 
6133 
6134 
6135 template <int dim, int spacedim>
6136 inline const DerivativeForm<1, spacedim, dim> &
6137 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
6138 {
6139  Assert(this->update_flags & update_inverse_jacobians,
6140  ExcAccessToUninitializedField("update_inverse_jacobians"));
6141  AssertIndexRange(i, this->mapping_output.inverse_jacobians.size());
6142  Assert(present_cell.is_initialized(), ExcNotReinited());
6143 
6144  return this->mapping_output.inverse_jacobians[i];
6145 }
6146 
6147 
6148 
6149 template <int dim, int spacedim>
6150 inline const Tensor<1, spacedim> &
6151 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
6152 {
6153  Assert(this->update_flags & update_normal_vectors,
6155  "update_normal_vectors")));
6156  AssertIndexRange(i, this->mapping_output.normal_vectors.size());
6157  Assert(present_cell.is_initialized(), ExcNotReinited());
6158 
6159  return this->mapping_output.normal_vectors[i];
6160 }
6161 
6162 
6163 
6164 /*--------------------- Inline functions: FEValues --------------------------*/
6165 
6166 
6167 template <int dim, int spacedim>
6168 inline const Quadrature<dim> &
6170 {
6171  return quadrature;
6172 }
6173 
6174 
6175 
6176 template <int dim, int spacedim>
6177 inline const FEValues<dim, spacedim> &
6179 {
6180  return *this;
6181 }
6182 
6183 
6184 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
6185 
6186 
6187 template <int dim, int spacedim>
6188 inline unsigned int
6190 {
6191  return present_face_no;
6192 }
6193 
6194 
6195 template <int dim, int spacedim>
6196 inline unsigned int
6198 {
6199  return present_face_index;
6200 }
6201 
6202 
6203 /*----------------------- Inline functions: FE*FaceValues -------------------*/
6204 
6205 template <int dim, int spacedim>
6206 inline const Quadrature<dim - 1> &
6208 {
6209  return quadrature[quadrature.size() == 1 ? 0 : present_face_no];
6210 }
6211 
6212 
6213 
6214 template <int dim, int spacedim>
6215 inline const FEFaceValues<dim, spacedim> &
6217 {
6218  return *this;
6219 }
6220 
6221 
6222 
6223 template <int dim, int spacedim>
6224 inline const FESubfaceValues<dim, spacedim> &
6226 {
6227  return *this;
6228 }
6229 
6230 
6231 
6232 template <int dim, int spacedim>
6233 inline const Tensor<1, spacedim> &
6234 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
6235 {
6236  AssertIndexRange(i, this->mapping_output.boundary_forms.size());
6237  Assert(this->update_flags & update_boundary_forms,
6239  "update_boundary_forms")));
6240 
6241  return this->mapping_output.boundary_forms[i];
6242 }
6243 
6244 #endif // DOXYGEN
6245 
6247 
6248 #endif
unsigned int get_face_index() const
unsigned int present_face_no
Definition: fe_values.h:4274
const Quadrature< dim - 1 > & get_quadrature() const
unsigned int present_face_index
Definition: fe_values.h:4280
const Tensor< 1, spacedim > & boundary_form(const unsigned int i) const
unsigned int get_face_number() const
const hp::QCollection< dim - 1 > quadrature
Definition: fe_values.h:4285
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell, const unsigned int face_no)
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell, const typename Triangulation< dim, spacedim >::face_iterator &face)
const FEFaceValues< dim, spacedim > & get_present_fe_values() const
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell, const unsigned int face_no, const unsigned int subface_no)
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell, const typename Triangulation< dim, spacedim >::face_iterator &face, const typename Triangulation< dim, spacedim >::face_iterator &subface)
const FESubfaceValues< dim, spacedim > & get_present_fe_values() const
CellIteratorContainer(const TriaIterator< DoFCellAccessor< dim, spacedim, lda >> &cell)
const DoFHandler< dim, spacedim > * dof_handler
Definition: fe_values.h:3886
Triangulation< dim, spacedim >::cell_iterator cell
Definition: fe_values.h:3885
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:4004
CellIteratorContainer present_cell
Definition: fe_values.h:3895
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition: fe_values.h:4019
FEValuesBase(const FEValuesBase &)=delete
const DerivativeForm< 2, dim, spacedim > & jacobian_grad(const unsigned int quadrature_point) const
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int quadrature_point) const
boost::signals2::connection tria_listener_mesh_transform
Definition: fe_values.h:3913
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_ending_at(const unsigned int end_dof_index) const
const FEValuesViews::Vector< dim, spacedim > & operator[](const FEValuesExtractors::Vector &vector) const
internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:3955
const Tensor< 1, spacedim > & shape_grad(const unsigned int function_no, const unsigned int quadrature_point) const
const std::vector< DerivativeForm< 1, spacedim, dim > > & get_inverse_jacobians() const
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:3940
UpdateFlags get_update_flags() const
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
const DerivativeForm< 3, dim, spacedim > & jacobian_2nd_derivative(const unsigned int quadrature_point) const
const unsigned int dofs_per_cell
Definition: fe_values.h:2451
UpdateFlags update_flags
Definition: fe_values.h:3986
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
Definition: fe_values.h:3964
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int quadrature_point) const
Tensor< 1, spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
const std::vector< double > & get_JxW_values() const
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
Tensor< 3, spacedim > shape_3rd_derivative_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_starting_at(const unsigned int start_dof_index) const
FEValuesBase & operator=(const FEValuesBase &)=delete
const unsigned int n_quadrature_points
Definition: fe_values.h:2433
const Tensor< 2, spacedim > & shape_hessian(const unsigned int function_no, const unsigned int point_no) const
Tensor< 2, spacedim > shape_hessian_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
const FEValuesViews::SymmetricTensor< 2, dim, spacedim > & operator[](const FEValuesExtractors::SymmetricTensor< 2 > &tensor) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices() const
const std::vector< DerivativeForm< 2, dim, spacedim > > & get_jacobian_grads() const
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition: fe_values.h:3948
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int quadrature_point) const
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
boost::signals2::connection tria_listener_refinement
Definition: fe_values.h:3904
const std::vector< DerivativeForm< 3, dim, spacedim > > & get_jacobian_2nd_derivatives() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > quadrature_point_indices() const
const Point< spacedim > & quadrature_point(const unsigned int q) const
const Tensor< 3, spacedim > & jacobian_pushed_forward_grad(const unsigned int quadrature_point) const
double JxW(const unsigned int quadrature_point) const
const std::vector< Tensor< 4, spacedim > > & get_jacobian_pushed_forward_2nd_derivatives() const
const Tensor< 1, spacedim > & normal_vector(const unsigned int i) const
const Tensor< 5, spacedim > & jacobian_pushed_forward_3rd_derivative(const unsigned int quadrature_point) const
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3980
const Tensor< 4, spacedim > & jacobian_pushed_forward_2nd_derivative(const unsigned int quadrature_point) const
const Tensor< 3, spacedim > & shape_3rd_derivative(const unsigned int function_no, const unsigned int point_no) const
const Mapping< dim, spacedim > & get_mapping() const
const std::vector< Tensor< 5, spacedim > > & get_jacobian_pushed_forward_3rd_derivatives() const
const FiniteElement< dim, spacedim > & get_fe() const
const std::vector< DerivativeForm< 4, dim, spacedim > > & get_jacobian_3rd_derivatives() const
const std::vector< Point< spacedim > > & get_quadrature_points() const
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > fe_data
Definition: fe_values.h:3972
const FEValuesViews::Tensor< 2, dim, spacedim > & operator[](const FEValuesExtractors::Tensor< 2 > &tensor) const
const unsigned int max_n_quadrature_points
Definition: fe_values.h:2444
typename ProductType< Number, hessian_type >::type solution_hessian_type
Definition: fe_values.h:215
value_type value(const unsigned int shape_function, const unsigned int q_point) const
void get_function_values_from_local_dof_values(const InputVector &dof_values, std::vector< solution_value_type< typename InputVector::value_type >> &values) const
Definition: fe_values.cc:1578
void get_function_values(const InputVector &fe_function, std::vector< solution_value_type< typename InputVector::value_type >> &values) const
Definition: fe_values.cc:1547
const unsigned int component
Definition: fe_values.h:635
void get_function_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< solution_gradient_type< typename InputVector::value_type >> &gradients) const
Definition: fe_values.cc:1632
void get_function_third_derivatives(const InputVector &fe_function, std::vector< solution_third_derivative_type< typename InputVector::value_type >> &third_derivatives) const
Definition: fe_values.cc:1764
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:640
::Tensor< 1, spacedim > gradient_type
Definition: fe_values.h:162
typename ProductType< Number, value_type >::type solution_laplacian_type
Definition: fe_values.h:205
third_derivative_type third_derivative(const unsigned int shape_function, const unsigned int q_point) const
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
Definition: fe_values.h:225
void get_function_gradients(const InputVector &fe_function, std::vector< solution_gradient_type< typename InputVector::value_type >> &gradients) const
Definition: fe_values.cc:1602
typename ProductType< Number, value_type >::type solution_value_type
Definition: fe_values.h:185
::Tensor< 2, spacedim > hessian_type
Definition: fe_values.h:169
typename ProductType< Number, gradient_type >::type solution_gradient_type
Definition: fe_values.h:195
void get_function_hessians_from_local_dof_values(const InputVector &dof_values, std::vector< solution_hessian_type< typename InputVector::value_type >> &hessians) const
Definition: fe_values.cc:1686
Scalar & operator=(const Scalar< dim, spacedim > &)=delete
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:629
Scalar(Scalar< dim, spacedim > &&)=default
void get_function_laplacians_from_local_dof_values(const InputVector &dof_values, std::vector< solution_laplacian_type< typename InputVector::value_type >> &laplacians) const
Definition: fe_values.cc:1740
Scalar(const Scalar< dim, spacedim > &)=delete
void get_function_laplacians(const InputVector &fe_function, std::vector< solution_laplacian_type< typename InputVector::value_type >> &laplacians) const
Definition: fe_values.cc:1710
void get_function_hessians(const InputVector &fe_function, std::vector< solution_hessian_type< typename InputVector::value_type >> &hessians) const
Definition: fe_values.cc:1656
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
Scalar & operator=(Scalar< dim, spacedim > &&) noexcept=default
void get_function_third_derivatives_from_local_dof_values(const InputVector &dof_values, std::vector< solution_third_derivative_type< typename InputVector::value_type >> &third_derivatives) const
Definition: fe_values.cc:1795
::Tensor< 3, spacedim > third_derivative_type
Definition: fe_values.h:176
SymmetricTensor(const SymmetricTensor< 2, dim, spacedim > &)=delete
SymmetricTensor(SymmetricTensor< 2, dim, spacedim > &&)=default
SymmetricTensor & operator=(const SymmetricTensor< 2, dim, spacedim > &)=delete
SymmetricTensor & operator=(SymmetricTensor< 2, dim, spacedim > &&) noexcept=default
typename ProductType< Number, value_type >::type solution_value_type
Definition: fe_values.h:1506
typename ProductType< Number, divergence_type >::type solution_divergence_type
Definition: fe_values.h:1516
typename ProductType< Number, value_type >::type solution_value_type
Definition: fe_values.h:1842
typename ProductType< Number, divergence_type >::type solution_divergence_type
Definition: fe_values.h:1852
typename ProductType< Number, gradient_type >::type solution_gradient_type
Definition: fe_values.h:1862
divergence_type divergence(const unsigned int shape_function, const unsigned int q_point) const
Tensor(const Tensor< 2, dim, spacedim > &)=delete
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:2190
value_type value(const unsigned int shape_function, const unsigned int q_point) const
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
Tensor(Tensor< 2, dim, spacedim > &&)=default
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:2201
Tensor & operator=(const Tensor< 2, dim, spacedim > &)=delete
Tensor & operator=(Tensor< 2, dim, spacedim > &&)=default
Vector & operator=(Vector< dim, spacedim > &&)=default
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
Definition: fe_values.h:812
typename ProductType< Number, divergence_type >::type solution_divergence_type
Definition: fe_values.h:773
Vector(const Vector< dim, spacedim > &)=delete
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
typename ProductType< Number, hessian_type >::type solution_hessian_type
Definition: fe_values.h:802
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1436
Vector & operator=(const Vector< dim, spacedim > &)=delete
typename ProductType< Number, symmetric_gradient_type >::type solution_symmetric_gradient_type
Definition: fe_values.h:763
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
divergence_type divergence(const unsigned int shape_function, const unsigned int q_point) const
third_derivative_type third_derivative(const unsigned int shape_function, const unsigned int q_point) const
typename ProductType< Number, gradient_type >::type solution_gradient_type
Definition: fe_values.h:753
typename ProductType< Number, value_type >::type solution_value_type
Definition: fe_values.h:743
symmetric_gradient_type symmetric_gradient(const unsigned int shape_function, const unsigned int q_point) const
Vector(Vector< dim, spacedim > &&)=default
typename ::internal::CurlType< spacedim >::type curl_type
Definition: fe_values.h:720
const unsigned int first_vector_component
Definition: fe_values.h:1442
typename ProductType< Number, curl_type >::type solution_curl_type
Definition: fe_values.h:792
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1447
value_type value(const unsigned int shape_function, const unsigned int q_point) const
curl_type curl(const unsigned int shape_function, const unsigned int q_point) const
typename ProductType< Number, value_type >::type solution_laplacian_type
Definition: fe_values.h:783
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell)
const Quadrature< dim > quadrature
Definition: fe_values.h:4155
const FEValues< dim, spacedim > & get_present_fe_values() const
const Quadrature< dim > & get_quadrature() const
Definition: tensor.h:503
static constexpr TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
typename Tensor< rank_ - 1, dim, Number >::tensor_type value_type
Definition: tensor.h:536
Definition: vector.h:109
#define DEAL_II_DEPRECATED
Definition: config.h:162
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:461
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:462
#define DeclException0(Exception0)
Definition: exceptions.h:464
static ::ExceptionBase & ExcAccessToUninitializedField()
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
Definition: exceptions.h:1509
static ::ExceptionBase & ExcNotImplemented()
#define AssertIndexRange(index, range)
Definition: exceptions.h:1768
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
static ::ExceptionBase & ExcFENotPrimitive()
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:509
static ::ExceptionBase & ExcMessage(std::string arg1)
TriaIterator< CellAccessor< dim, spacedim > > cell_iterator
Definition: tria.h:1368
UpdateFlags
@ update_jacobian_pushed_forward_2nd_derivatives
@ update_jacobian_pushed_forward_grads
@ update_hessians
Second derivatives of shape functions.
@ update_jacobian_3rd_derivatives
@ update_values
Shape function values.
@ update_jacobian_grads
Gradient of volume element.
@ update_normal_vectors
Normal vectors.
@ update_3rd_derivatives
Third derivatives of shape functions.
@ update_JxW_values
Transformed quadrature weights.
@ update_jacobians
Volume element.
@ update_inverse_jacobians
Volume element.
@ update_gradients
Shape function gradients.
@ update_quadrature_points
Transformed quadrature points.
@ update_jacobian_pushed_forward_3rd_derivatives
@ update_boundary_forms
Outer normal vector, not normalized.
@ update_jacobian_2nd_derivatives
typename ::internal::FEValuesViews::ViewType< dim, spacedim, Extractor >::type View
Definition: fe_values.h:2309
static const char U
std::enable_if_t< std::is_fundamental< T >::value, std::size_t > memory_consumption(const T &t)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
unsigned int global_dof_index
Definition: types.h:81
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:258
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:242
typename ProductType< Number, typename Scalar< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:250
typename ProductType< Number, typename Scalar< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:266
typename ProductType< Number, typename Scalar< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:274
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1542
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1534
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:1896
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1880
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1888
typename ProductType< Number, typename Vector< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:885
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:861
typename ProductType< Number, typename Vector< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:877
typename ProductType< Number, typename Vector< dim, spacedim >::symmetric_gradient_type >::type symmetric_gradient_type
Definition: fe_values.h:845
typename ProductType< Number, typename Vector< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:837
typename ProductType< Number, typename Vector< dim, spacedim >::curl_type >::type curl_type
Definition: fe_values.h:869
typename ProductType< Number, typename Vector< dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:853
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:829
typename internal::ProductTypeImpl< typename std::decay< T >::type, typename std::decay< U >::type >::type type
std::vector<::FEValuesViews::Scalar< dim, spacedim > > scalars
Definition: fe_values.h:2286
std::vector<::FEValuesViews::Vector< dim, spacedim > > vectors
Definition: fe_values.h:2287
std::vector<::FEValuesViews::SymmetricTensor< 2, dim, spacedim > > symmetric_second_order_tensors
Definition: fe_values.h:2289
std::vector<::FEValuesViews::Tensor< 2, dim, spacedim > > second_order_tensors
Definition: fe_values.h:2291
typename ::FEValuesViews::SymmetricTensor< rank, dim, spacedim > type
Definition: fe_values.h:2269
typename ::FEValuesViews::Tensor< rank, dim, spacedim > type
Definition: fe_values.h:2255
constexpr SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)