Reference documentation for deal.II version 9.3.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
fe_values.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2021 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>
37 
38 #include <deal.II/grid/tria.h>
40 
42 
43 #include <algorithm>
44 #include <memory>
45 #include <type_traits>
46 
47 
48 // dummy include in order to have the
49 // definition of PetscScalar available
50 // without including other PETSc stuff
51 #ifdef DEAL_II_WITH_PETSC
52 # include <petsc.h>
53 #endif
54 
56 
57 // Forward declaration
58 #ifndef DOXYGEN
59 template <int dim, int spacedim = dim>
60 class FEValuesBase;
61 #endif
62 
63 namespace internal
64 {
69  template <int dim, class NumberType = double>
70  struct CurlType;
71 
78  template <class NumberType>
79  struct CurlType<1, NumberType>
80  {
82  };
83 
90  template <class NumberType>
91  struct CurlType<2, NumberType>
92  {
94  };
95 
102  template <class NumberType>
103  struct CurlType<3, NumberType>
104  {
106  };
107 } // namespace internal
108 
109 
110 
132 namespace FEValuesViews
133 {
145  template <int dim, int spacedim = dim>
146  class Scalar
147  {
148  public:
155 
162 
169 
176 
183  template <typename Number>
185 
192  template <typename Number>
193  using solution_gradient_type =
195 
202  template <typename Number>
205 
212  template <typename Number>
213  using solution_hessian_type =
215 
222  template <typename Number>
225 
232  template <typename Number>
234  {
239  using value_type =
240  typename ProductType<Number,
242 
247  using gradient_type = typename ProductType<
248  Number,
250 
255  using laplacian_type =
256  typename ProductType<Number,
258 
263  using hessian_type = typename ProductType<
264  Number,
266 
271  using third_derivative_type = typename ProductType<
272  Number,
274  };
275 
281  {
291 
300  unsigned int row_index;
301  };
302 
306  Scalar();
307 
313  Scalar(const FEValuesBase<dim, spacedim> &fe_values_base,
314  const unsigned int component);
315 
320  Scalar(const Scalar<dim, spacedim> &) = delete;
321 
325  // NOLINTNEXTLINE OSX does not compile with noexcept
326  Scalar(Scalar<dim, spacedim> &&) = default;
327 
331  ~Scalar() = default;
332 
337  Scalar &
338  operator=(const Scalar<dim, spacedim> &) = delete;
339 
343  Scalar &
344  operator=(Scalar<dim, spacedim> &&) noexcept = default;
345 
359  value_type
360  value(const unsigned int shape_function, const unsigned int q_point) const;
361 
373  gradient(const unsigned int shape_function,
374  const unsigned int q_point) const;
375 
387  hessian(const unsigned int shape_function,
388  const unsigned int q_point) const;
389 
401  third_derivative(const unsigned int shape_function,
402  const unsigned int q_point) const;
403 
421  template <class InputVector>
422  void
423  get_function_values(
424  const InputVector &fe_function,
425  std::vector<solution_value_type<typename InputVector::value_type>>
426  &values) const;
427 
462  template <class InputVector>
463  void
464  get_function_values_from_local_dof_values(
465  const InputVector &dof_values,
466  std::vector<solution_value_type<typename InputVector::value_type>>
467  &values) const;
468 
486  template <class InputVector>
487  void
488  get_function_gradients(
489  const InputVector &fe_function,
490  std::vector<solution_gradient_type<typename InputVector::value_type>>
491  &gradients) const;
492 
499  template <class InputVector>
500  void
501  get_function_gradients_from_local_dof_values(
502  const InputVector &dof_values,
503  std::vector<solution_gradient_type<typename InputVector::value_type>>
504  &gradients) const;
505 
523  template <class InputVector>
524  void
525  get_function_hessians(
526  const InputVector &fe_function,
527  std::vector<solution_hessian_type<typename InputVector::value_type>>
528  &hessians) const;
529 
536  template <class InputVector>
537  void
538  get_function_hessians_from_local_dof_values(
539  const InputVector &dof_values,
540  std::vector<solution_hessian_type<typename InputVector::value_type>>
541  &hessians) const;
542 
543 
562  template <class InputVector>
563  void
564  get_function_laplacians(
565  const InputVector &fe_function,
566  std::vector<solution_laplacian_type<typename InputVector::value_type>>
567  &laplacians) const;
568 
575  template <class InputVector>
576  void
577  get_function_laplacians_from_local_dof_values(
578  const InputVector &dof_values,
579  std::vector<solution_laplacian_type<typename InputVector::value_type>>
580  &laplacians) const;
581 
582 
601  template <class InputVector>
602  void
603  get_function_third_derivatives(
604  const InputVector &fe_function,
605  std::vector<
606  solution_third_derivative_type<typename InputVector::value_type>>
607  &third_derivatives) const;
608 
615  template <class InputVector>
616  void
617  get_function_third_derivatives_from_local_dof_values(
618  const InputVector &dof_values,
619  std::vector<
620  solution_third_derivative_type<typename InputVector::value_type>>
621  &third_derivatives) const;
622 
623 
624  private:
628  const SmartPointer<const FEValuesBase<dim, spacedim>> fe_values;
629 
634  const unsigned int component;
635 
639  std::vector<ShapeFunctionData> shape_function_data;
640  };
641 
642 
643 
673  template <int dim, int spacedim = dim>
674  class Vector
675  {
676  public:
683 
693 
705 
712 
719  using curl_type = typename ::internal::CurlType<spacedim>::type;
720 
727 
734 
741  template <typename Number>
743 
750  template <typename Number>
751  using solution_gradient_type =
753 
760  template <typename Number>
763 
770  template <typename Number>
773 
780  template <typename Number>
783 
790  template <typename Number>
792 
799  template <typename Number>
800  using solution_hessian_type =
802 
809  template <typename Number>
812 
819  template <typename Number>
821  {
826  using value_type =
827  typename ProductType<Number,
829 
834  using gradient_type = typename ProductType<
835  Number,
837 
842  using symmetric_gradient_type = typename ProductType<
843  Number,
845 
850  using divergence_type = typename ProductType<
851  Number,
853 
858  using laplacian_type =
859  typename ProductType<Number,
861 
866  using curl_type =
867  typename ProductType<Number,
869 
874  using hessian_type = typename ProductType<
875  Number,
877 
882  using third_derivative_type = typename ProductType<
883  Number,
885  };
886 
892  {
901  bool is_nonzero_shape_function_component[spacedim];
902 
912  unsigned int row_index[spacedim];
913 
924  };
925 
929  Vector();
930 
939  Vector(const FEValuesBase<dim, spacedim> &fe_values_base,
940  const unsigned int first_vector_component);
941 
946  Vector(const Vector<dim, spacedim> &) = delete;
947 
951  // NOLINTNEXTLINE OSX does not compile with noexcept
952  Vector(Vector<dim, spacedim> &&) = default;
953 
957  ~Vector() = default;
958 
963  Vector &
964  operator=(const Vector<dim, spacedim> &) = delete;
965 
969  // NOLINTNEXTLINE OSX does not compile with noexcept
970  Vector &
971  operator=(Vector<dim, spacedim> &&) = default; // NOLINT
972 
989  value_type
990  value(const unsigned int shape_function, const unsigned int q_point) const;
991 
1006  gradient(const unsigned int shape_function,
1007  const unsigned int q_point) const;
1008 
1025  symmetric_gradient(const unsigned int shape_function,
1026  const unsigned int q_point) const;
1027 
1039  divergence(const unsigned int shape_function,
1040  const unsigned int q_point) const;
1041 
1062  curl_type
1063  curl(const unsigned int shape_function, const unsigned int q_point) const;
1064 
1075  hessian_type
1076  hessian(const unsigned int shape_function,
1077  const unsigned int q_point) const;
1078 
1090  third_derivative(const unsigned int shape_function,
1091  const unsigned int q_point) const;
1092 
1110  template <class InputVector>
1111  void
1112  get_function_values(
1113  const InputVector &fe_function,
1115  &values) const;
1116 
1151  template <class InputVector>
1152  void
1153  get_function_values_from_local_dof_values(
1154  const InputVector &dof_values,
1156  &values) const;
1157 
1175  template <class InputVector>
1176  void
1177  get_function_gradients(
1178  const InputVector &fe_function,
1180  &gradients) const;
1181 
1188  template <class InputVector>
1189  void
1190  get_function_gradients_from_local_dof_values(
1191  const InputVector &dof_values,
1193  &gradients) const;
1194 
1218  template <class InputVector>
1219  void
1220  get_function_symmetric_gradients(
1221  const InputVector &fe_function,
1222  std::vector<
1224  &symmetric_gradients) const;
1225 
1232  template <class InputVector>
1233  void
1234  get_function_symmetric_gradients_from_local_dof_values(
1235  const InputVector &dof_values,
1236  std::vector<
1238  &symmetric_gradients) const;
1239 
1258  template <class InputVector>
1259  void
1260  get_function_divergences(
1261  const InputVector &fe_function,
1263  &divergences) const;
1264 
1271  template <class InputVector>
1272  void
1273  get_function_divergences_from_local_dof_values(
1274  const InputVector &dof_values,
1276  &divergences) const;
1277 
1296  template <class InputVector>
1297  void
1298  get_function_curls(
1299  const InputVector &fe_function,
1301  const;
1302 
1309  template <class InputVector>
1310  void
1311  get_function_curls_from_local_dof_values(
1312  const InputVector &dof_values,
1314  const;
1315 
1333  template <class InputVector>
1334  void
1335  get_function_hessians(
1336  const InputVector &fe_function,
1338  &hessians) const;
1339 
1346  template <class InputVector>
1347  void
1348  get_function_hessians_from_local_dof_values(
1349  const InputVector &dof_values,
1351  &hessians) const;
1352 
1371  template <class InputVector>
1372  void
1373  get_function_laplacians(
1374  const InputVector &fe_function,
1376  &laplacians) const;
1377 
1384  template <class InputVector>
1385  void
1386  get_function_laplacians_from_local_dof_values(
1387  const InputVector &dof_values,
1389  &laplacians) const;
1390 
1409  template <class InputVector>
1410  void
1411  get_function_third_derivatives(
1412  const InputVector &fe_function,
1413  std::vector<
1415  &third_derivatives) const;
1416 
1423  template <class InputVector>
1424  void
1425  get_function_third_derivatives_from_local_dof_values(
1426  const InputVector &dof_values,
1427  std::vector<
1429  &third_derivatives) const;
1430 
1431  private:
1436 
1441  const unsigned int first_vector_component;
1442 
1446  std::vector<ShapeFunctionData> shape_function_data;
1447  };
1448 
1449 
1450  template <int rank, int dim, int spacedim = dim>
1452 
1475  template <int dim, int spacedim>
1476  class SymmetricTensor<2, dim, spacedim>
1477  {
1478  public:
1486 
1497 
1504  template <typename Number>
1506 
1513  template <typename Number>
1514  using solution_divergence_type =
1516 
1517 
1524  template <typename Number>
1525  struct DEAL_II_DEPRECATED_EARLY OutputType
1526  {
1531  using value_type = typename ProductType<
1532  Number,
1534 
1539  using divergence_type = typename ProductType<
1540  Number,
1542  };
1543 
1548  struct ShapeFunctionData
1549  {
1558  bool is_nonzero_shape_function_component
1559  [value_type::n_independent_components];
1560 
1570  unsigned int row_index[value_type::n_independent_components];
1571 
1581 
1586  };
1587 
1591  SymmetricTensor();
1592 
1602  SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1603  const unsigned int first_tensor_component);
1604 
1610 
1614  // NOLINTNEXTLINE OSX does not compile with noexcept
1616 
1621  SymmetricTensor &
1622  operator=(const SymmetricTensor<2, dim, spacedim> &) = delete;
1623 
1627  SymmetricTensor &
1628  operator=(SymmetricTensor<2, dim, spacedim> &&) noexcept = default;
1629 
1647  value_type
1648  value(const unsigned int shape_function, const unsigned int q_point) const;
1649 
1664  divergence(const unsigned int shape_function,
1665  const unsigned int q_point) const;
1666 
1684  template <class InputVector>
1685  void
1686  get_function_values(
1687  const InputVector &fe_function,
1688  std::vector<solution_value_type<typename InputVector::value_type>>
1689  &values) const;
1690 
1725  template <class InputVector>
1726  void
1727  get_function_values_from_local_dof_values(
1728  const InputVector &dof_values,
1729  std::vector<solution_value_type<typename InputVector::value_type>>
1730  &values) const;
1731 
1753  template <class InputVector>
1754  void
1755  get_function_divergences(
1756  const InputVector &fe_function,
1757  std::vector<solution_divergence_type<typename InputVector::value_type>>
1758  &divergences) const;
1759 
1766  template <class InputVector>
1767  void
1768  get_function_divergences_from_local_dof_values(
1769  const InputVector &dof_values,
1770  std::vector<solution_divergence_type<typename InputVector::value_type>>
1771  &divergences) const;
1772 
1773  private:
1777  const SmartPointer<const FEValuesBase<dim, spacedim>> fe_values;
1778 
1783  const unsigned int first_tensor_component;
1784 
1788  std::vector<ShapeFunctionData> shape_function_data;
1789  };
1790 
1791 
1792  template <int rank, int dim, int spacedim = dim>
1793  class Tensor;
1794 
1813  template <int dim, int spacedim>
1814  class Tensor<2, dim, spacedim>
1815  {
1816  public:
1822 
1827 
1833 
1840  template <typename Number>
1842 
1849  template <typename Number>
1850  using solution_divergence_type =
1852 
1859  template <typename Number>
1860  using solution_gradient_type =
1862 
1863 
1870  template <typename Number>
1871  struct DEAL_II_DEPRECATED_EARLY OutputType
1872  {
1877  using value_type = typename ProductType<
1878  Number,
1880 
1885  using divergence_type = typename ProductType<
1886  Number,
1888 
1893  using gradient_type = typename ProductType<
1894  Number,
1896  };
1897 
1902  struct ShapeFunctionData
1903  {
1912  bool is_nonzero_shape_function_component
1913  [value_type::n_independent_components];
1914 
1924  unsigned int row_index[value_type::n_independent_components];
1925 
1935 
1940  };
1941 
1945  Tensor();
1946 
1951  Tensor(const Tensor<2, dim, spacedim> &) = delete;
1952 
1956  // NOLINTNEXTLINE OSX does not compile with noexcept
1957  Tensor(Tensor<2, dim, spacedim> &&) = default;
1958 
1962  ~Tensor() = default;
1963 
1973  Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1974  const unsigned int first_tensor_component);
1975 
1976 
1981  Tensor &
1982  operator=(const Tensor<2, dim, spacedim> &) = delete;
1983 
1987  // NOLINTNEXTLINE
1988  Tensor &operator=(Tensor<2, dim, spacedim> &&) = default;
1989 
2006  value_type
2007  value(const unsigned int shape_function, const unsigned int q_point) const;
2008 
2023  divergence(const unsigned int shape_function,
2024  const unsigned int q_point) const;
2025 
2040  gradient(const unsigned int shape_function,
2041  const unsigned int q_point) const;
2042 
2060  template <class InputVector>
2061  void
2062  get_function_values(
2063  const InputVector &fe_function,
2065  &values) const;
2066 
2101  template <class InputVector>
2102  void
2103  get_function_values_from_local_dof_values(
2104  const InputVector &dof_values,
2106  &values) const;
2107 
2129  template <class InputVector>
2130  void
2131  get_function_divergences(
2132  const InputVector &fe_function,
2134  &divergences) const;
2135 
2142  template <class InputVector>
2143  void
2144  get_function_divergences_from_local_dof_values(
2145  const InputVector &dof_values,
2147  &divergences) const;
2148 
2165  template <class InputVector>
2166  void
2167  get_function_gradients(
2168  const InputVector &fe_function,
2170  &gradients) const;
2171 
2178  template <class InputVector>
2179  void
2180  get_function_gradients_from_local_dof_values(
2181  const InputVector &dof_values,
2183  &gradients) const;
2184 
2185  private:
2190 
2195  const unsigned int first_tensor_component;
2196 
2200  std::vector<ShapeFunctionData> shape_function_data;
2201  };
2202 
2203 } // namespace FEValuesViews
2204 
2205 
2206 namespace internal
2207 {
2208  namespace FEValuesViews
2209  {
2214  template <int dim, int spacedim, typename Extractor>
2215  struct ViewType
2216  {};
2217 
2225  template <int dim, int spacedim>
2226  struct ViewType<dim, spacedim, FEValuesExtractors::Scalar>
2227  {
2228  using type = typename ::FEValuesViews::Scalar<dim, spacedim>;
2229  };
2230 
2238  template <int dim, int spacedim>
2239  struct ViewType<dim, spacedim, FEValuesExtractors::Vector>
2240  {
2241  using type = typename ::FEValuesViews::Vector<dim, spacedim>;
2242  };
2243 
2251  template <int dim, int spacedim, int rank>
2252  struct ViewType<dim, spacedim, FEValuesExtractors::Tensor<rank>>
2253  {
2254  using type = typename ::FEValuesViews::Tensor<rank, dim, spacedim>;
2255  };
2256 
2264  template <int dim, int spacedim, int rank>
2265  struct ViewType<dim, spacedim, FEValuesExtractors::SymmetricTensor<rank>>
2266  {
2267  using type =
2268  typename ::FEValuesViews::SymmetricTensor<rank, dim, spacedim>;
2269  };
2270 
2278  template <int dim, int spacedim>
2279  struct Cache
2280  {
2285  std::vector<::FEValuesViews::Scalar<dim, spacedim>> scalars;
2286  std::vector<::FEValuesViews::Vector<dim, spacedim>> vectors;
2287  std::vector<::FEValuesViews::SymmetricTensor<2, dim, spacedim>>
2289  std::vector<::FEValuesViews::Tensor<2, dim, spacedim>>
2291 
2295  Cache(const FEValuesBase<dim, spacedim> &fe_values);
2296  };
2297  } // namespace FEValuesViews
2298 } // namespace internal
2299 
2300 namespace FEValuesViews
2301 {
2306  template <int dim, int spacedim, typename Extractor>
2307  using View = typename ::internal::FEValuesViews::
2308  ViewType<dim, spacedim, Extractor>::type;
2309 } // namespace FEValuesViews
2310 
2311 
2411 template <int dim, int spacedim>
2412 class FEValuesBase : public Subscriptor
2413 {
2414 public:
2418  static const unsigned int dimension = dim;
2419 
2423  static const unsigned int space_dimension = spacedim;
2424 
2432  const unsigned int n_quadrature_points;
2433 
2443  const unsigned int max_n_quadrature_points;
2444 
2450  const unsigned int dofs_per_cell;
2451 
2452 
2460  FEValuesBase(const unsigned int n_q_points,
2461  const unsigned int dofs_per_cell,
2462  const UpdateFlags update_flags,
2463  const Mapping<dim, spacedim> & mapping,
2464  const FiniteElement<dim, spacedim> &fe);
2465 
2470  FEValuesBase &
2471  operator=(const FEValuesBase &) = delete;
2472 
2477  FEValuesBase(const FEValuesBase &) = delete;
2478 
2482  virtual ~FEValuesBase() override;
2483 
2484 
2488 
2489 
2510  const double &
2511  shape_value(const unsigned int function_no,
2512  const unsigned int point_no) const;
2513 
2534  double
2535  shape_value_component(const unsigned int function_no,
2536  const unsigned int point_no,
2537  const unsigned int component) const;
2538 
2564  const Tensor<1, spacedim> &
2565  shape_grad(const unsigned int function_no,
2566  const unsigned int quadrature_point) const;
2567 
2585  shape_grad_component(const unsigned int function_no,
2586  const unsigned int point_no,
2587  const unsigned int component) const;
2588 
2608  const Tensor<2, spacedim> &
2609  shape_hessian(const unsigned int function_no,
2610  const unsigned int point_no) const;
2611 
2629  shape_hessian_component(const unsigned int function_no,
2630  const unsigned int point_no,
2631  const unsigned int component) const;
2632 
2652  const Tensor<3, spacedim> &
2653  shape_3rd_derivative(const unsigned int function_no,
2654  const unsigned int point_no) const;
2655 
2673  shape_3rd_derivative_component(const unsigned int function_no,
2674  const unsigned int point_no,
2675  const unsigned int component) const;
2676 
2678 
2680 
2717  template <class InputVector>
2718  void
2719  get_function_values(
2720  const InputVector & fe_function,
2721  std::vector<typename InputVector::value_type> &values) const;
2722 
2736  template <class InputVector>
2737  void
2738  get_function_values(
2739  const InputVector & fe_function,
2740  std::vector<Vector<typename InputVector::value_type>> &values) const;
2741 
2798  template <class InputVector>
2799  void
2800  get_function_values(
2801  const InputVector & fe_function,
2803  std::vector<typename InputVector::value_type> & values) const;
2804 
2813  template <class InputVector>
2814  void
2815  get_function_values(
2816  const InputVector & fe_function,
2818  std::vector<Vector<typename InputVector::value_type>> &values) const;
2819 
2820 
2842  template <class InputVector>
2843  void
2844  get_function_values(
2845  const InputVector & fe_function,
2847  ArrayView<std::vector<typename InputVector::value_type>> values,
2848  const bool quadrature_points_fastest) const;
2849 
2851 
2853 
2890  template <class InputVector>
2891  void
2892  get_function_gradients(
2893  const InputVector &fe_function,
2895  &gradients) const;
2896 
2913  template <class InputVector>
2914  void
2915  get_function_gradients(
2916  const InputVector &fe_function,
2917  std::vector<
2919  &gradients) const;
2920 
2929  template <class InputVector>
2930  void
2931  get_function_gradients(
2932  const InputVector & fe_function,
2935  &gradients) const;
2936 
2945  template <class InputVector>
2946  void
2947  get_function_gradients(
2948  const InputVector & fe_function,
2950  ArrayView<
2952  gradients,
2953  const bool quadrature_points_fastest = false) const;
2954 
2956 
2960 
2998  template <class InputVector>
2999  void
3000  get_function_hessians(
3001  const InputVector &fe_function,
3003  &hessians) const;
3004 
3022  template <class InputVector>
3023  void
3024  get_function_hessians(
3025  const InputVector &fe_function,
3026  std::vector<
3028  & hessians,
3029  const bool quadrature_points_fastest = false) const;
3030 
3039  template <class InputVector>
3040  void
3041  get_function_hessians(
3042  const InputVector & fe_function,
3045  &hessians) const;
3046 
3055  template <class InputVector>
3056  void
3057  get_function_hessians(
3058  const InputVector & fe_function,
3060  ArrayView<
3062  hessians,
3063  const bool quadrature_points_fastest = false) const;
3064 
3105  template <class InputVector>
3106  void
3107  get_function_laplacians(
3108  const InputVector & fe_function,
3109  std::vector<typename InputVector::value_type> &laplacians) const;
3110 
3130  template <class InputVector>
3131  void
3132  get_function_laplacians(
3133  const InputVector & fe_function,
3134  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
3135 
3144  template <class InputVector>
3145  void
3146  get_function_laplacians(
3147  const InputVector & fe_function,
3149  std::vector<typename InputVector::value_type> & laplacians) const;
3150 
3159  template <class InputVector>
3160  void
3161  get_function_laplacians(
3162  const InputVector & fe_function,
3164  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
3165 
3174  template <class InputVector>
3175  void
3176  get_function_laplacians(
3177  const InputVector & fe_function,
3179  std::vector<std::vector<typename InputVector::value_type>> &laplacians,
3180  const bool quadrature_points_fastest = false) const;
3181 
3183 
3185 
3224  template <class InputVector>
3225  void
3226  get_function_third_derivatives(
3227  const InputVector &fe_function,
3229  &third_derivatives) const;
3230 
3249  template <class InputVector>
3250  void
3251  get_function_third_derivatives(
3252  const InputVector &fe_function,
3253  std::vector<
3255  & third_derivatives,
3256  const bool quadrature_points_fastest = false) const;
3257 
3266  template <class InputVector>
3267  void
3268  get_function_third_derivatives(
3269  const InputVector & fe_function,
3272  &third_derivatives) const;
3273 
3282  template <class InputVector>
3283  void
3284  get_function_third_derivatives(
3285  const InputVector & fe_function,
3287  ArrayView<
3289  third_derivatives,
3290  const bool quadrature_points_fastest = false) const;
3292 
3294 
3295 
3320  dof_indices() const;
3321 
3355  dof_indices_starting_at(const unsigned int start_dof_index) const;
3356 
3388  dof_indices_ending_at(const unsigned int end_dof_index) const;
3389 
3391 
3393 
3394 
3416  quadrature_point_indices() const;
3417 
3423  const Point<spacedim> &
3424  quadrature_point(const unsigned int q) const;
3425 
3431  const std::vector<Point<spacedim>> &
3432  get_quadrature_points() const;
3433 
3449  double
3450  JxW(const unsigned int quadrature_point) const;
3451 
3455  const std::vector<double> &
3456  get_JxW_values() const;
3457 
3465  jacobian(const unsigned int quadrature_point) const;
3466 
3473  const std::vector<DerivativeForm<1, dim, spacedim>> &
3474  get_jacobians() const;
3475 
3484  jacobian_grad(const unsigned int quadrature_point) const;
3485 
3492  const std::vector<DerivativeForm<2, dim, spacedim>> &
3493  get_jacobian_grads() const;
3494 
3503  const Tensor<3, spacedim> &
3504  jacobian_pushed_forward_grad(const unsigned int quadrature_point) const;
3505 
3512  const std::vector<Tensor<3, spacedim>> &
3513  get_jacobian_pushed_forward_grads() const;
3514 
3523  jacobian_2nd_derivative(const unsigned int quadrature_point) const;
3524 
3531  const std::vector<DerivativeForm<3, dim, spacedim>> &
3532  get_jacobian_2nd_derivatives() const;
3533 
3543  const Tensor<4, spacedim> &
3544  jacobian_pushed_forward_2nd_derivative(
3545  const unsigned int quadrature_point) const;
3546 
3553  const std::vector<Tensor<4, spacedim>> &
3554  get_jacobian_pushed_forward_2nd_derivatives() const;
3555 
3565  jacobian_3rd_derivative(const unsigned int quadrature_point) const;
3566 
3573  const std::vector<DerivativeForm<4, dim, spacedim>> &
3574  get_jacobian_3rd_derivatives() const;
3575 
3585  const Tensor<5, spacedim> &
3586  jacobian_pushed_forward_3rd_derivative(
3587  const unsigned int quadrature_point) const;
3588 
3595  const std::vector<Tensor<5, spacedim>> &
3596  get_jacobian_pushed_forward_3rd_derivatives() const;
3597 
3605  inverse_jacobian(const unsigned int quadrature_point) const;
3606 
3613  const std::vector<DerivativeForm<1, spacedim, dim>> &
3614  get_inverse_jacobians() const;
3615 
3635  const Tensor<1, spacedim> &
3636  normal_vector(const unsigned int i) const;
3637 
3645  const std::vector<Tensor<1, spacedim>> &
3646  get_normal_vectors() const;
3647 
3649 
3651 
3652 
3662  operator[](const FEValuesExtractors::Scalar &scalar) const;
3663 
3673  operator[](const FEValuesExtractors::Vector &vector) const;
3674 
3685  operator[](const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
3686 
3687 
3697  operator[](const FEValuesExtractors::Tensor<2> &tensor) const;
3698 
3700 
3702 
3703 
3707  const Mapping<dim, spacedim> &
3708  get_mapping() const;
3709 
3714  get_fe() const;
3715 
3719  UpdateFlags
3720  get_update_flags() const;
3721 
3726  get_cell() const;
3727 
3734  get_cell_similarity() const;
3735 
3740  std::size_t
3741  memory_consumption() const;
3743 
3744 
3753  std::string,
3754  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3755  << "object for which this kind of information has not been computed. What "
3756  << "information these objects compute is determined by the update_* flags you "
3757  << "pass to the constructor. Here, the operation you are attempting requires "
3758  << "the <" << arg1
3759  << "> flag to be set, but it was apparently not specified "
3760  << "upon construction.");
3761 
3769  ExcFEDontMatch,
3770  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3771  "to the DoFHandler that provided the cell iterator do not match.");
3777  DeclException1(ExcShapeFunctionNotPrimitive,
3778  int,
3779  << "The shape function with index " << arg1
3780  << " is not primitive, i.e. it is vector-valued and "
3781  << "has more than one non-zero vector component. This "
3782  << "function cannot be called for these shape functions. "
3783  << "Maybe you want to use the same function with the "
3784  << "_component suffix?");
3785 
3794  "The given FiniteElement is not a primitive element but the requested operation "
3795  "only works for those. See FiniteElement::is_primitive() for more information.");
3796 
3797 protected:
3828  class CellIteratorBase;
3829 
3834  template <typename CI>
3837 
3843  std::unique_ptr<const CellIteratorBase> present_cell;
3844 
3852  boost::signals2::connection tria_listener_refinement;
3853 
3861  boost::signals2::connection tria_listener_mesh_transform;
3862 
3868  void
3869  invalidate_present_cell();
3870 
3880  void
3881  maybe_invalidate_previous_present_cell(
3882  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3883 
3889 
3895  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3897 
3904 
3905 
3913 
3919  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3921 
3927  spacedim>
3929 
3930 
3935 
3944  UpdateFlags
3945  compute_update_flags(const UpdateFlags update_flags) const;
3946 
3953 
3959  void
3960  check_cell_similarity(
3961  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3962 
3963 private:
3968 
3969  // Make the view classes friends of this class, since they access internal
3970  // data.
3971  template <int, int>
3973  template <int, int>
3975  template <int, int, int>
3977  template <int, int, int>
3979 };
3980 
3981 
3982 
3992 template <int dim, int spacedim = dim>
3993 class FEValues : public FEValuesBase<dim, spacedim>
3994 {
3995 public:
4000  static const unsigned int integral_dimension = dim;
4001 
4006  FEValues(const Mapping<dim, spacedim> & mapping,
4007  const FiniteElement<dim, spacedim> &fe,
4008  const Quadrature<dim> & quadrature,
4009  const UpdateFlags update_flags);
4010 
4017  FEValues(const Mapping<dim, spacedim> & mapping,
4018  const FiniteElement<dim, spacedim> &fe,
4019  const hp::QCollection<dim> & quadrature,
4020  const UpdateFlags update_flags);
4021 
4028  const Quadrature<dim> & quadrature,
4029  const UpdateFlags update_flags);
4030 
4038  const hp::QCollection<dim> & quadrature,
4039  const UpdateFlags update_flags);
4040 
4047  template <bool level_dof_access>
4048  void
4049  reinit(
4051 
4065  void
4067 
4072  const Quadrature<dim> &
4073  get_quadrature() const;
4074 
4079  std::size_t
4080  memory_consumption() const;
4081 
4096  const FEValues<dim, spacedim> &
4097  get_present_fe_values() const;
4098 
4099 private:
4104 
4108  void
4109  initialize(const UpdateFlags update_flags);
4110 
4117  void
4118  do_reinit();
4119 };
4120 
4121 
4131 template <int dim, int spacedim = dim>
4132 class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
4133 {
4134 public:
4139  static const unsigned int integral_dimension = dim - 1;
4140 
4152  FEFaceValuesBase(const unsigned int dofs_per_cell,
4153  const UpdateFlags update_flags,
4154  const Mapping<dim, spacedim> & mapping,
4155  const FiniteElement<dim, spacedim> &fe,
4156  const Quadrature<dim - 1> & quadrature);
4157 
4164  FEFaceValuesBase(const unsigned int dofs_per_cell,
4165  const UpdateFlags update_flags,
4166  const Mapping<dim, spacedim> & mapping,
4167  const FiniteElement<dim, spacedim> &fe,
4168  const hp::QCollection<dim - 1> & quadrature);
4169 
4177  const Tensor<1, spacedim> &
4178  boundary_form(const unsigned int i) const;
4179 
4186  const std::vector<Tensor<1, spacedim>> &
4187  get_boundary_forms() const;
4188 
4193  unsigned int
4194  get_face_index() const;
4195 
4200  const Quadrature<dim - 1> &
4201  get_quadrature() const;
4202 
4207  std::size_t
4208  memory_consumption() const;
4209 
4210 protected:
4215  unsigned int present_face_no;
4216 
4221  unsigned int present_face_index;
4222 
4226  const hp::QCollection<dim - 1> quadrature;
4227 };
4228 
4229 
4230 
4244 template <int dim, int spacedim = dim>
4245 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
4246 {
4247 public:
4252  static const unsigned int dimension = dim;
4253 
4254  static const unsigned int space_dimension = spacedim;
4255 
4260  static const unsigned int integral_dimension = dim - 1;
4261 
4266  FEFaceValues(const Mapping<dim, spacedim> & mapping,
4267  const FiniteElement<dim, spacedim> &fe,
4268  const Quadrature<dim - 1> & quadrature,
4269  const UpdateFlags update_flags);
4270 
4277  FEFaceValues(const Mapping<dim, spacedim> & mapping,
4278  const FiniteElement<dim, spacedim> &fe,
4279  const hp::QCollection<dim - 1> & quadrature,
4280  const UpdateFlags update_flags);
4281 
4288  const Quadrature<dim - 1> & quadrature,
4289  const UpdateFlags update_flags);
4290 
4298  const hp::QCollection<dim - 1> & quadrature,
4299  const UpdateFlags update_flags);
4300 
4305  template <bool level_dof_access>
4306  void
4307  reinit(
4309  const unsigned int face_no);
4310 
4317  template <bool level_dof_access>
4318  void
4319  reinit(
4321  const typename Triangulation<dim, spacedim>::face_iterator & face);
4322 
4336  void
4338  const unsigned int face_no);
4339 
4340  /*
4341  * Reinitialize the gradients, Jacobi determinants, etc for the given face
4342  * on a given cell of type "iterator into a Triangulation object", and the
4343  * given finite element. Since iterators into a triangulation alone only
4344  * convey information about the geometry of a cell, but not about degrees of
4345  * freedom possibly associated with this cell, you will not be able to call
4346  * some functions of this class if they need information about degrees of
4347  * freedom. These functions are, above all, the
4348  * <tt>get_function_value/gradients/hessians/third_derivatives</tt>
4349  * functions. If you want to call these functions, you have to call the @p
4350  * reinit variants that take iterators into DoFHandler or other DoF handler
4351  * type objects.
4352  *
4353  * @note @p face must be one of @p cell's face iterators.
4354  */
4355  void
4357  const typename Triangulation<dim, spacedim>::face_iterator &face);
4358 
4374  get_present_fe_values() const;
4375 
4376 private:
4380  void
4381  initialize(const UpdateFlags update_flags);
4382 
4389  void
4390  do_reinit(const unsigned int face_no);
4391 };
4392 
4393 
4410 template <int dim, int spacedim = dim>
4411 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
4412 {
4413 public:
4417  static const unsigned int dimension = dim;
4418 
4422  static const unsigned int space_dimension = spacedim;
4423 
4428  static const unsigned int integral_dimension = dim - 1;
4429 
4434  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
4435  const FiniteElement<dim, spacedim> &fe,
4436  const Quadrature<dim - 1> & face_quadrature,
4437  const UpdateFlags update_flags);
4438 
4445  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
4446  const FiniteElement<dim, spacedim> &fe,
4447  const hp::QCollection<dim - 1> & face_quadrature,
4448  const UpdateFlags update_flags);
4449 
4456  const Quadrature<dim - 1> & face_quadrature,
4457  const UpdateFlags update_flags);
4458 
4466  const hp::QCollection<dim - 1> & face_quadrature,
4467  const UpdateFlags update_flags);
4468 
4475  template <bool level_dof_access>
4476  void
4477  reinit(
4479  const unsigned int face_no,
4480  const unsigned int subface_no);
4481 
4486  template <bool level_dof_access>
4487  void
4488  reinit(
4490  const typename Triangulation<dim, spacedim>::face_iterator & face,
4491  const typename Triangulation<dim, spacedim>::face_iterator &subface);
4492 
4506  void
4508  const unsigned int face_no,
4509  const unsigned int subface_no);
4510 
4530  void
4532  const typename Triangulation<dim, spacedim>::face_iterator &face,
4533  const typename Triangulation<dim, spacedim>::face_iterator &subface);
4534 
4550  get_present_fe_values() const;
4551 
4557  DeclException0(ExcReinitCalledWithBoundaryFace);
4558 
4564  DeclException0(ExcFaceHasNoSubfaces);
4565 
4566 private:
4570  void
4571  initialize(const UpdateFlags update_flags);
4572 
4579  void
4580  do_reinit(const unsigned int face_no, const unsigned int subface_no);
4581 };
4582 
4583 
4584 #ifndef DOXYGEN
4585 
4586 
4587 /*------------------------ Inline functions: namespace FEValuesViews --------*/
4588 
4589 namespace FEValuesViews
4590 {
4591  template <int dim, int spacedim>
4592  inline typename Scalar<dim, spacedim>::value_type
4593  Scalar<dim, spacedim>::value(const unsigned int shape_function,
4594  const unsigned int q_point) const
4595  {
4596  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4597  Assert(
4598  fe_values->update_flags & update_values,
4600  "update_values"))));
4601 
4602  // an adaptation of the FEValuesBase::shape_value_component function
4603  // except that here we know the component as fixed and we have
4604  // pre-computed and cached a bunch of information. See the comments there.
4605  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4606  return fe_values->finite_element_output.shape_values(
4607  shape_function_data[shape_function].row_index, q_point);
4608  else
4609  return 0;
4610  }
4611 
4612 
4613 
4614  template <int dim, int spacedim>
4615  inline typename Scalar<dim, spacedim>::gradient_type
4616  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
4617  const unsigned int q_point) const
4618  {
4619  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4620  Assert(fe_values->update_flags & update_gradients,
4622  "update_gradients")));
4623 
4624  // an adaptation of the FEValuesBase::shape_grad_component
4625  // function except that here we know the component as fixed and we have
4626  // pre-computed and cached a bunch of information. See the comments there.
4627  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4628  return fe_values->finite_element_output
4629  .shape_gradients[shape_function_data[shape_function].row_index]
4630  [q_point];
4631  else
4632  return gradient_type();
4633  }
4634 
4635 
4636 
4637  template <int dim, int spacedim>
4638  inline typename Scalar<dim, spacedim>::hessian_type
4639  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
4640  const unsigned int q_point) const
4641  {
4642  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4643  Assert(fe_values->update_flags & update_hessians,
4645  "update_hessians")));
4646 
4647  // an adaptation of the FEValuesBase::shape_hessian_component
4648  // function except that here we know the component as fixed and we have
4649  // pre-computed and cached a bunch of information. See the comments there.
4650  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4651  return fe_values->finite_element_output
4652  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
4653  else
4654  return hessian_type();
4655  }
4656 
4657 
4658 
4659  template <int dim, int spacedim>
4660  inline typename Scalar<dim, spacedim>::third_derivative_type
4661  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
4662  const unsigned int q_point) const
4663  {
4664  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4665  Assert(fe_values->update_flags & update_3rd_derivatives,
4667  "update_3rd_derivatives")));
4668 
4669  // an adaptation of the FEValuesBase::shape_3rdderivative_component
4670  // function except that here we know the component as fixed and we have
4671  // pre-computed and cached a bunch of information. See the comments there.
4672  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4673  return fe_values->finite_element_output
4674  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
4675  [q_point];
4676  else
4677  return third_derivative_type();
4678  }
4679 
4680 
4681 
4682  template <int dim, int spacedim>
4683  inline typename Vector<dim, spacedim>::value_type
4684  Vector<dim, spacedim>::value(const unsigned int shape_function,
4685  const unsigned int q_point) const
4686  {
4687  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4688  Assert(fe_values->update_flags & update_values,
4690  "update_values")));
4691 
4692  // same as for the scalar case except that we have one more index
4693  const int snc =
4694  shape_function_data[shape_function].single_nonzero_component;
4695  if (snc == -2)
4696  return value_type();
4697  else if (snc != -1)
4698  {
4699  value_type return_value;
4700  return_value[shape_function_data[shape_function]
4701  .single_nonzero_component_index] =
4702  fe_values->finite_element_output.shape_values(snc, q_point);
4703  return return_value;
4704  }
4705  else
4706  {
4707  value_type return_value;
4708  for (unsigned int d = 0; d < dim; ++d)
4709  if (shape_function_data[shape_function]
4710  .is_nonzero_shape_function_component[d])
4711  return_value[d] = fe_values->finite_element_output.shape_values(
4712  shape_function_data[shape_function].row_index[d], q_point);
4713 
4714  return return_value;
4715  }
4716  }
4717 
4718 
4719 
4720  template <int dim, int spacedim>
4721  inline typename Vector<dim, spacedim>::gradient_type
4722  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4723  const unsigned int q_point) const
4724  {
4725  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4726  Assert(fe_values->update_flags & update_gradients,
4728  "update_gradients")));
4729 
4730  // same as for the scalar case except that we have one more index
4731  const int snc =
4732  shape_function_data[shape_function].single_nonzero_component;
4733  if (snc == -2)
4734  return gradient_type();
4735  else if (snc != -1)
4736  {
4737  gradient_type return_value;
4738  return_value[shape_function_data[shape_function]
4739  .single_nonzero_component_index] =
4740  fe_values->finite_element_output.shape_gradients[snc][q_point];
4741  return return_value;
4742  }
4743  else
4744  {
4745  gradient_type return_value;
4746  for (unsigned int d = 0; d < dim; ++d)
4747  if (shape_function_data[shape_function]
4748  .is_nonzero_shape_function_component[d])
4749  return_value[d] =
4750  fe_values->finite_element_output.shape_gradients
4751  [shape_function_data[shape_function].row_index[d]][q_point];
4752 
4753  return return_value;
4754  }
4755  }
4756 
4757 
4758 
4759  template <int dim, int spacedim>
4760  inline typename Vector<dim, spacedim>::divergence_type
4761  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4762  const unsigned int q_point) const
4763  {
4764  // this function works like in the case above
4765  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4766  Assert(fe_values->update_flags & update_gradients,
4768  "update_gradients")));
4769 
4770  // same as for the scalar case except that we have one more index
4771  const int snc =
4772  shape_function_data[shape_function].single_nonzero_component;
4773  if (snc == -2)
4774  return divergence_type();
4775  else if (snc != -1)
4776  return fe_values->finite_element_output
4777  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4778  .single_nonzero_component_index];
4779  else
4780  {
4781  divergence_type return_value = 0;
4782  for (unsigned int d = 0; d < dim; ++d)
4783  if (shape_function_data[shape_function]
4784  .is_nonzero_shape_function_component[d])
4785  return_value +=
4786  fe_values->finite_element_output.shape_gradients
4787  [shape_function_data[shape_function].row_index[d]][q_point][d];
4788 
4789  return return_value;
4790  }
4791  }
4792 
4793 
4794 
4795  template <int dim, int spacedim>
4796  inline typename Vector<dim, spacedim>::curl_type
4797  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4798  const unsigned int q_point) const
4799  {
4800  // this function works like in the case above
4801 
4802  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4803  Assert(fe_values->update_flags & update_gradients,
4805  "update_gradients")));
4806  // same as for the scalar case except that we have one more index
4807  const int snc =
4808  shape_function_data[shape_function].single_nonzero_component;
4809 
4810  if (snc == -2)
4811  return curl_type();
4812 
4813  else
4814  switch (dim)
4815  {
4816  case 1:
4817  {
4818  Assert(false,
4819  ExcMessage(
4820  "Computing the curl in 1d is not a useful operation"));
4821  return curl_type();
4822  }
4823 
4824  case 2:
4825  {
4826  if (snc != -1)
4827  {
4828  curl_type return_value;
4829 
4830  // the single nonzero component can only be zero or one in 2d
4831  if (shape_function_data[shape_function]
4832  .single_nonzero_component_index == 0)
4833  return_value[0] =
4834  -1.0 * fe_values->finite_element_output
4835  .shape_gradients[snc][q_point][1];
4836  else
4837  return_value[0] = fe_values->finite_element_output
4838  .shape_gradients[snc][q_point][0];
4839 
4840  return return_value;
4841  }
4842 
4843  else
4844  {
4845  curl_type return_value;
4846 
4847  return_value[0] = 0.0;
4848 
4849  if (shape_function_data[shape_function]
4850  .is_nonzero_shape_function_component[0])
4851  return_value[0] -=
4852  fe_values->finite_element_output
4853  .shape_gradients[shape_function_data[shape_function]
4854  .row_index[0]][q_point][1];
4855 
4856  if (shape_function_data[shape_function]
4857  .is_nonzero_shape_function_component[1])
4858  return_value[0] +=
4859  fe_values->finite_element_output
4860  .shape_gradients[shape_function_data[shape_function]
4861  .row_index[1]][q_point][0];
4862 
4863  return return_value;
4864  }
4865  }
4866 
4867  case 3:
4868  {
4869  if (snc != -1)
4870  {
4871  curl_type return_value;
4872 
4873  switch (shape_function_data[shape_function]
4874  .single_nonzero_component_index)
4875  {
4876  case 0:
4877  {
4878  return_value[0] = 0;
4879  return_value[1] = fe_values->finite_element_output
4880  .shape_gradients[snc][q_point][2];
4881  return_value[2] =
4882  -1.0 * fe_values->finite_element_output
4883  .shape_gradients[snc][q_point][1];
4884  return return_value;
4885  }
4886 
4887  case 1:
4888  {
4889  return_value[0] =
4890  -1.0 * fe_values->finite_element_output
4891  .shape_gradients[snc][q_point][2];
4892  return_value[1] = 0;
4893  return_value[2] = fe_values->finite_element_output
4894  .shape_gradients[snc][q_point][0];
4895  return return_value;
4896  }
4897 
4898  default:
4899  {
4900  return_value[0] = fe_values->finite_element_output
4901  .shape_gradients[snc][q_point][1];
4902  return_value[1] =
4903  -1.0 * fe_values->finite_element_output
4904  .shape_gradients[snc][q_point][0];
4905  return_value[2] = 0;
4906  return return_value;
4907  }
4908  }
4909  }
4910 
4911  else
4912  {
4913  curl_type return_value;
4914 
4915  for (unsigned int i = 0; i < dim; ++i)
4916  return_value[i] = 0.0;
4917 
4918  if (shape_function_data[shape_function]
4919  .is_nonzero_shape_function_component[0])
4920  {
4921  return_value[1] +=
4922  fe_values->finite_element_output
4923  .shape_gradients[shape_function_data[shape_function]
4924  .row_index[0]][q_point][2];
4925  return_value[2] -=
4926  fe_values->finite_element_output
4927  .shape_gradients[shape_function_data[shape_function]
4928  .row_index[0]][q_point][1];
4929  }
4930 
4931  if (shape_function_data[shape_function]
4932  .is_nonzero_shape_function_component[1])
4933  {
4934  return_value[0] -=
4935  fe_values->finite_element_output
4936  .shape_gradients[shape_function_data[shape_function]
4937  .row_index[1]][q_point][2];
4938  return_value[2] +=
4939  fe_values->finite_element_output
4940  .shape_gradients[shape_function_data[shape_function]
4941  .row_index[1]][q_point][0];
4942  }
4943 
4944  if (shape_function_data[shape_function]
4945  .is_nonzero_shape_function_component[2])
4946  {
4947  return_value[0] +=
4948  fe_values->finite_element_output
4949  .shape_gradients[shape_function_data[shape_function]
4950  .row_index[2]][q_point][1];
4951  return_value[1] -=
4952  fe_values->finite_element_output
4953  .shape_gradients[shape_function_data[shape_function]
4954  .row_index[2]][q_point][0];
4955  }
4956 
4957  return return_value;
4958  }
4959  }
4960  }
4961  // should not end up here
4962  Assert(false, ExcInternalError());
4963  return curl_type();
4964  }
4965 
4966 
4967 
4968  template <int dim, int spacedim>
4969  inline typename Vector<dim, spacedim>::hessian_type
4970  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
4971  const unsigned int q_point) const
4972  {
4973  // this function works like in the case above
4974  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4975  Assert(fe_values->update_flags & update_hessians,
4977  "update_hessians")));
4978 
4979  // same as for the scalar case except that we have one more index
4980  const int snc =
4981  shape_function_data[shape_function].single_nonzero_component;
4982  if (snc == -2)
4983  return hessian_type();
4984  else if (snc != -1)
4985  {
4986  hessian_type return_value;
4987  return_value[shape_function_data[shape_function]
4988  .single_nonzero_component_index] =
4989  fe_values->finite_element_output.shape_hessians[snc][q_point];
4990  return return_value;
4991  }
4992  else
4993  {
4994  hessian_type return_value;
4995  for (unsigned int d = 0; d < dim; ++d)
4996  if (shape_function_data[shape_function]
4997  .is_nonzero_shape_function_component[d])
4998  return_value[d] =
4999  fe_values->finite_element_output.shape_hessians
5000  [shape_function_data[shape_function].row_index[d]][q_point];
5001 
5002  return return_value;
5003  }
5004  }
5005 
5006 
5007 
5008  template <int dim, int spacedim>
5009  inline typename Vector<dim, spacedim>::third_derivative_type
5010  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
5011  const unsigned int q_point) const
5012  {
5013  // this function works like in the case above
5014  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5015  Assert(fe_values->update_flags & update_3rd_derivatives,
5017  "update_3rd_derivatives")));
5018 
5019  // same as for the scalar case except that we have one more index
5020  const int snc =
5021  shape_function_data[shape_function].single_nonzero_component;
5022  if (snc == -2)
5023  return third_derivative_type();
5024  else if (snc != -1)
5025  {
5026  third_derivative_type return_value;
5027  return_value[shape_function_data[shape_function]
5028  .single_nonzero_component_index] =
5029  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
5030  return return_value;
5031  }
5032  else
5033  {
5034  third_derivative_type return_value;
5035  for (unsigned int d = 0; d < dim; ++d)
5036  if (shape_function_data[shape_function]
5037  .is_nonzero_shape_function_component[d])
5038  return_value[d] =
5039  fe_values->finite_element_output.shape_3rd_derivatives
5040  [shape_function_data[shape_function].row_index[d]][q_point];
5041 
5042  return return_value;
5043  }
5044  }
5045 
5046 
5047 
5048  namespace internal
5049  {
5054  inline ::SymmetricTensor<2, 1>
5055  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
5056  {
5057  AssertIndexRange(n, 1);
5058  (void)n;
5059 
5060  return {{t[0]}};
5061  }
5062 
5063 
5064 
5065  inline ::SymmetricTensor<2, 2>
5066  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
5067  {
5068  switch (n)
5069  {
5070  case 0:
5071  {
5072  return {{t[0], 0, t[1] / 2}};
5073  }
5074  case 1:
5075  {
5076  return {{0, t[1], t[0] / 2}};
5077  }
5078  default:
5079  {
5080  AssertIndexRange(n, 2);
5081  return {};
5082  }
5083  }
5084  }
5085 
5086 
5087 
5088  inline ::SymmetricTensor<2, 3>
5089  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
5090  {
5091  switch (n)
5092  {
5093  case 0:
5094  {
5095  return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
5096  }
5097  case 1:
5098  {
5099  return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
5100  }
5101  case 2:
5102  {
5103  return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
5104  }
5105  default:
5106  {
5107  AssertIndexRange(n, 3);
5108  return {};
5109  }
5110  }
5111  }
5112  } // namespace internal
5113 
5114 
5115 
5116  template <int dim, int spacedim>
5117  inline typename Vector<dim, spacedim>::symmetric_gradient_type
5118  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
5119  const unsigned int q_point) const
5120  {
5121  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5122  Assert(fe_values->update_flags & update_gradients,
5124  "update_gradients")));
5125 
5126  // same as for the scalar case except that we have one more index
5127  const int snc =
5128  shape_function_data[shape_function].single_nonzero_component;
5129  if (snc == -2)
5130  return symmetric_gradient_type();
5131  else if (snc != -1)
5132  return internal::symmetrize_single_row(
5133  shape_function_data[shape_function].single_nonzero_component_index,
5134  fe_values->finite_element_output.shape_gradients[snc][q_point]);
5135  else
5136  {
5137  gradient_type return_value;
5138  for (unsigned int d = 0; d < dim; ++d)
5139  if (shape_function_data[shape_function]
5140  .is_nonzero_shape_function_component[d])
5141  return_value[d] =
5142  fe_values->finite_element_output.shape_gradients
5143  [shape_function_data[shape_function].row_index[d]][q_point];
5144 
5145  return symmetrize(return_value);
5146  }
5147  }
5148 
5149 
5150 
5151  template <int dim, int spacedim>
5153  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
5154  const unsigned int q_point) const
5155  {
5156  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5157  Assert(fe_values->update_flags & update_values,
5159  "update_values")));
5160 
5161  // similar to the vector case where we have more then one index and we need
5162  // to convert between unrolled and component indexing for tensors
5163  const int snc =
5164  shape_function_data[shape_function].single_nonzero_component;
5165 
5166  if (snc == -2)
5167  {
5168  // shape function is zero for the selected components
5169  return value_type();
5170  }
5171  else if (snc != -1)
5172  {
5173  value_type return_value;
5174  const unsigned int comp =
5175  shape_function_data[shape_function].single_nonzero_component_index;
5176  return_value[value_type::unrolled_to_component_indices(comp)] =
5177  fe_values->finite_element_output.shape_values(snc, q_point);
5178  return return_value;
5179  }
5180  else
5181  {
5182  value_type return_value;
5183  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
5184  if (shape_function_data[shape_function]
5185  .is_nonzero_shape_function_component[d])
5186  return_value[value_type::unrolled_to_component_indices(d)] =
5187  fe_values->finite_element_output.shape_values(
5188  shape_function_data[shape_function].row_index[d], q_point);
5189  return return_value;
5190  }
5191  }
5192 
5193 
5194 
5195  template <int dim, int spacedim>
5198  const unsigned int shape_function,
5199  const unsigned int q_point) const
5200  {
5201  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5202  Assert(fe_values->update_flags & update_gradients,
5204  "update_gradients")));
5205 
5206  const int snc =
5207  shape_function_data[shape_function].single_nonzero_component;
5208 
5209  if (snc == -2)
5210  {
5211  // shape function is zero for the selected components
5212  return divergence_type();
5213  }
5214  else if (snc != -1)
5215  {
5216  // we have a single non-zero component when the symmetric tensor is
5217  // represented in unrolled form. this implies we potentially have
5218  // two non-zero components when represented in component form! we
5219  // will only have one non-zero entry if the non-zero component lies on
5220  // the diagonal of the tensor.
5221  //
5222  // the divergence of a second-order tensor is a first order tensor.
5223  //
5224  // assume the second-order tensor is A with components A_{ij}. then
5225  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
5226  // entries in the tensorial representation. define the
5227  // divergence as:
5228  // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
5229  // (which is incidentally also
5230  // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
5231  // In both cases, a sum is implied.
5232  //
5233  // Now, we know the nonzero component in unrolled form: it is indicated
5234  // by 'snc'. we can figure out which tensor components belong to this:
5235  const unsigned int comp =
5236  shape_function_data[shape_function].single_nonzero_component_index;
5237  const unsigned int ii =
5238  value_type::unrolled_to_component_indices(comp)[0];
5239  const unsigned int jj =
5240  value_type::unrolled_to_component_indices(comp)[1];
5241 
5242  // given the form of the divergence above, if ii=jj there is only a
5243  // single nonzero component of the full tensor and the gradient
5244  // equals
5245  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
5246  // all other entries of 'b' are zero
5247  //
5248  // on the other hand, if ii!=jj, then there are two nonzero entries in
5249  // the full tensor and
5250  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
5251  // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
5252  // again, all other entries of 'b' are zero
5253  const ::Tensor<1, spacedim> &phi_grad =
5254  fe_values->finite_element_output.shape_gradients[snc][q_point];
5255 
5256  divergence_type return_value;
5257  return_value[ii] = phi_grad[jj];
5258 
5259  if (ii != jj)
5260  return_value[jj] = phi_grad[ii];
5261 
5262  return return_value;
5263  }
5264  else
5265  {
5266  Assert(false, ExcNotImplemented());
5267  divergence_type return_value;
5268  return return_value;
5269  }
5270  }
5271 
5272 
5273 
5274  template <int dim, int spacedim>
5275  inline typename Tensor<2, dim, spacedim>::value_type
5276  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
5277  const unsigned int q_point) const
5278  {
5279  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5280  Assert(fe_values->update_flags & update_values,
5282  "update_values")));
5283 
5284  // similar to the vector case where we have more then one index and we need
5285  // to convert between unrolled and component indexing for tensors
5286  const int snc =
5287  shape_function_data[shape_function].single_nonzero_component;
5288 
5289  if (snc == -2)
5290  {
5291  // shape function is zero for the selected components
5292  return value_type();
5293  }
5294  else if (snc != -1)
5295  {
5296  value_type return_value;
5297  const unsigned int comp =
5298  shape_function_data[shape_function].single_nonzero_component_index;
5299  const TableIndices<2> indices =
5301  return_value[indices] =
5302  fe_values->finite_element_output.shape_values(snc, q_point);
5303  return return_value;
5304  }
5305  else
5306  {
5307  value_type return_value;
5308  for (unsigned int d = 0; d < dim * dim; ++d)
5309  if (shape_function_data[shape_function]
5310  .is_nonzero_shape_function_component[d])
5311  {
5312  const TableIndices<2> indices =
5314  return_value[indices] =
5315  fe_values->finite_element_output.shape_values(
5316  shape_function_data[shape_function].row_index[d], q_point);
5317  }
5318  return return_value;
5319  }
5320  }
5321 
5322 
5323 
5324  template <int dim, int spacedim>
5326  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
5327  const unsigned int q_point) const
5328  {
5329  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5330  Assert(fe_values->update_flags & update_gradients,
5332  "update_gradients")));
5333 
5334  const int snc =
5335  shape_function_data[shape_function].single_nonzero_component;
5336 
5337  if (snc == -2)
5338  {
5339  // shape function is zero for the selected components
5340  return divergence_type();
5341  }
5342  else if (snc != -1)
5343  {
5344  // we have a single non-zero component when the tensor is
5345  // represented in unrolled form.
5346  //
5347  // the divergence of a second-order tensor is a first order tensor.
5348  //
5349  // assume the second-order tensor is A with components A_{ij},
5350  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
5351  //
5352  // Now, we know the nonzero component in unrolled form: it is indicated
5353  // by 'snc'. we can figure out which tensor components belong to this:
5354  const unsigned int comp =
5355  shape_function_data[shape_function].single_nonzero_component_index;
5356  const TableIndices<2> indices =
5358  const unsigned int ii = indices[0];
5359  const unsigned int jj = indices[1];
5360 
5361  const ::Tensor<1, spacedim> &phi_grad =
5362  fe_values->finite_element_output.shape_gradients[snc][q_point];
5363 
5364  divergence_type return_value;
5365  // note that we contract \nabla from the right
5366  return_value[ii] = phi_grad[jj];
5367 
5368  return return_value;
5369  }
5370  else
5371  {
5372  Assert(false, ExcNotImplemented());
5373  divergence_type return_value;
5374  return return_value;
5375  }
5376  }
5377 
5378 
5379 
5380  template <int dim, int spacedim>
5382  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
5383  const unsigned int q_point) const
5384  {
5385  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5386  Assert(fe_values->update_flags & update_gradients,
5388  "update_gradients")));
5389 
5390  const int snc =
5391  shape_function_data[shape_function].single_nonzero_component;
5392 
5393  if (snc == -2)
5394  {
5395  // shape function is zero for the selected components
5396  return gradient_type();
5397  }
5398  else if (snc != -1)
5399  {
5400  // we have a single non-zero component when the tensor is
5401  // represented in unrolled form.
5402  //
5403  // the gradient of a second-order tensor is a third order tensor.
5404  //
5405  // assume the second-order tensor is A with components A_{ij},
5406  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
5407  //
5408  // Now, we know the nonzero component in unrolled form: it is indicated
5409  // by 'snc'. we can figure out which tensor components belong to this:
5410  const unsigned int comp =
5411  shape_function_data[shape_function].single_nonzero_component_index;
5412  const TableIndices<2> indices =
5414  const unsigned int ii = indices[0];
5415  const unsigned int jj = indices[1];
5416 
5417  const ::Tensor<1, spacedim> &phi_grad =
5418  fe_values->finite_element_output.shape_gradients[snc][q_point];
5419 
5420  gradient_type return_value;
5421  return_value[ii][jj] = phi_grad;
5422 
5423  return return_value;
5424  }
5425  else
5426  {
5427  Assert(false, ExcNotImplemented());
5428  gradient_type return_value;
5429  return return_value;
5430  }
5431  }
5432 
5433 } // namespace FEValuesViews
5434 
5435 
5436 
5437 /*---------------------- Inline functions: FEValuesBase ---------------------*/
5438 
5439 
5440 
5441 template <int dim, int spacedim>
5443  operator[](const FEValuesExtractors::Scalar &scalar) const
5444 {
5445  AssertIndexRange(scalar.component, fe_values_views_cache.scalars.size());
5446 
5447  return fe_values_views_cache.scalars[scalar.component];
5448 }
5449 
5450 
5451 
5452 template <int dim, int spacedim>
5454  operator[](const FEValuesExtractors::Vector &vector) const
5455 {
5457  fe_values_views_cache.vectors.size());
5458 
5459  return fe_values_views_cache.vectors[vector.first_vector_component];
5460 }
5461 
5462 
5463 
5464 template <int dim, int spacedim>
5468 {
5469  Assert(
5470  tensor.first_tensor_component <
5471  fe_values_views_cache.symmetric_second_order_tensors.size(),
5473  0,
5474  fe_values_views_cache.symmetric_second_order_tensors.size()));
5475 
5476  return fe_values_views_cache
5477  .symmetric_second_order_tensors[tensor.first_tensor_component];
5478 }
5479 
5480 
5481 
5482 template <int dim, int spacedim>
5485  operator[](const FEValuesExtractors::Tensor<2> &tensor) const
5486 {
5488  fe_values_views_cache.second_order_tensors.size());
5489 
5490  return fe_values_views_cache
5491  .second_order_tensors[tensor.first_tensor_component];
5492 }
5493 
5494 
5495 
5496 template <int dim, int spacedim>
5497 inline const double &
5498 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
5499  const unsigned int j) const
5500 {
5501  AssertIndexRange(i, fe->n_dofs_per_cell());
5502  Assert(this->update_flags & update_values,
5503  ExcAccessToUninitializedField("update_values"));
5504  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5505  Assert(present_cell.get() != nullptr,
5506  ExcMessage("FEValues object is not reinit'ed to any cell"));
5507  // if the entire FE is primitive,
5508  // then we can take a short-cut:
5509  if (fe->is_primitive())
5510  return this->finite_element_output.shape_values(i, j);
5511  else
5512  {
5513  // otherwise, use the mapping
5514  // between shape function
5515  // numbers and rows. note that
5516  // by the assertions above, we
5517  // know that this particular
5518  // shape function is primitive,
5519  // so we can call
5520  // system_to_component_index
5521  const unsigned int row =
5522  this->finite_element_output
5523  .shape_function_to_row_table[i * fe->n_components() +
5524  fe->system_to_component_index(i).first];
5525  return this->finite_element_output.shape_values(row, j);
5526  }
5527 }
5528 
5529 
5530 
5531 template <int dim, int spacedim>
5532 inline double
5534  const unsigned int i,
5535  const unsigned int j,
5536  const unsigned int component) const
5537 {
5538  AssertIndexRange(i, fe->n_dofs_per_cell());
5539  Assert(this->update_flags & update_values,
5540  ExcAccessToUninitializedField("update_values"));
5541  AssertIndexRange(component, fe->n_components());
5542  Assert(present_cell.get() != nullptr,
5543  ExcMessage("FEValues object is not reinit'ed to any cell"));
5544 
5545  // check whether the shape function
5546  // is non-zero at all within
5547  // this component:
5548  if (fe->get_nonzero_components(i)[component] == false)
5549  return 0;
5550 
5551  // look up the right row in the
5552  // table and take the data from
5553  // there
5554  const unsigned int row =
5555  this->finite_element_output
5556  .shape_function_to_row_table[i * fe->n_components() + component];
5557  return this->finite_element_output.shape_values(row, j);
5558 }
5559 
5560 
5561 
5562 template <int dim, int spacedim>
5563 inline const Tensor<1, spacedim> &
5564 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
5565  const unsigned int j) const
5566 {
5567  AssertIndexRange(i, fe->n_dofs_per_cell());
5568  Assert(this->update_flags & update_gradients,
5569  ExcAccessToUninitializedField("update_gradients"));
5570  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5571  Assert(present_cell.get() != nullptr,
5572  ExcMessage("FEValues object is not reinit'ed to any cell"));
5573  // if the entire FE is primitive,
5574  // then we can take a short-cut:
5575  if (fe->is_primitive())
5576  return this->finite_element_output.shape_gradients[i][j];
5577  else
5578  {
5579  // otherwise, use the mapping
5580  // between shape function
5581  // numbers and rows. note that
5582  // by the assertions above, we
5583  // know that this particular
5584  // shape function is primitive,
5585  // so we can call
5586  // system_to_component_index
5587  const unsigned int row =
5588  this->finite_element_output
5589  .shape_function_to_row_table[i * fe->n_components() +
5590  fe->system_to_component_index(i).first];
5591  return this->finite_element_output.shape_gradients[row][j];
5592  }
5593 }
5594 
5595 
5596 
5597 template <int dim, int spacedim>
5598 inline Tensor<1, spacedim>
5600  const unsigned int i,
5601  const unsigned int j,
5602  const unsigned int component) const
5603 {
5604  AssertIndexRange(i, fe->n_dofs_per_cell());
5605  Assert(this->update_flags & update_gradients,
5606  ExcAccessToUninitializedField("update_gradients"));
5607  AssertIndexRange(component, fe->n_components());
5608  Assert(present_cell.get() != nullptr,
5609  ExcMessage("FEValues object is not reinit'ed to any cell"));
5610  // check whether the shape function
5611  // is non-zero at all within
5612  // this component:
5613  if (fe->get_nonzero_components(i)[component] == false)
5614  return Tensor<1, spacedim>();
5615 
5616  // look up the right row in the
5617  // table and take the data from
5618  // there
5619  const unsigned int row =
5620  this->finite_element_output
5621  .shape_function_to_row_table[i * fe->n_components() + component];
5622  return this->finite_element_output.shape_gradients[row][j];
5623 }
5624 
5625 
5626 
5627 template <int dim, int spacedim>
5628 inline const Tensor<2, spacedim> &
5629 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
5630  const unsigned int j) const
5631 {
5632  AssertIndexRange(i, fe->n_dofs_per_cell());
5633  Assert(this->update_flags & update_hessians,
5634  ExcAccessToUninitializedField("update_hessians"));
5635  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5636  Assert(present_cell.get() != nullptr,
5637  ExcMessage("FEValues object is not reinit'ed to any cell"));
5638  // if the entire FE is primitive,
5639  // then we can take a short-cut:
5640  if (fe->is_primitive())
5641  return this->finite_element_output.shape_hessians[i][j];
5642  else
5643  {
5644  // otherwise, use the mapping
5645  // between shape function
5646  // numbers and rows. note that
5647  // by the assertions above, we
5648  // know that this particular
5649  // shape function is primitive,
5650  // so we can call
5651  // system_to_component_index
5652  const unsigned int row =
5653  this->finite_element_output
5654  .shape_function_to_row_table[i * fe->n_components() +
5655  fe->system_to_component_index(i).first];
5656  return this->finite_element_output.shape_hessians[row][j];
5657  }
5658 }
5659 
5660 
5661 
5662 template <int dim, int spacedim>
5663 inline Tensor<2, spacedim>
5665  const unsigned int i,
5666  const unsigned int j,
5667  const unsigned int component) const
5668 {
5669  AssertIndexRange(i, fe->n_dofs_per_cell());
5670  Assert(this->update_flags & update_hessians,
5671  ExcAccessToUninitializedField("update_hessians"));
5672  AssertIndexRange(component, fe->n_components());
5673  Assert(present_cell.get() != nullptr,
5674  ExcMessage("FEValues object is not reinit'ed to any cell"));
5675  // check whether the shape function
5676  // is non-zero at all within
5677  // this component:
5678  if (fe->get_nonzero_components(i)[component] == false)
5679  return Tensor<2, spacedim>();
5680 
5681  // look up the right row in the
5682  // table and take the data from
5683  // there
5684  const unsigned int row =
5685  this->finite_element_output
5686  .shape_function_to_row_table[i * fe->n_components() + component];
5687  return this->finite_element_output.shape_hessians[row][j];
5688 }
5689 
5690 
5691 
5692 template <int dim, int spacedim>
5693 inline const Tensor<3, spacedim> &
5695  const unsigned int j) const
5696 {
5697  AssertIndexRange(i, fe->n_dofs_per_cell());
5698  Assert(this->update_flags & update_3rd_derivatives,
5699  ExcAccessToUninitializedField("update_3rd_derivatives"));
5700  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5701  Assert(present_cell.get() != nullptr,
5702  ExcMessage("FEValues object is not reinit'ed to any cell"));
5703  // if the entire FE is primitive,
5704  // then we can take a short-cut:
5705  if (fe->is_primitive())
5706  return this->finite_element_output.shape_3rd_derivatives[i][j];
5707  else
5708  {
5709  // otherwise, use the mapping
5710  // between shape function
5711  // numbers and rows. note that
5712  // by the assertions above, we
5713  // know that this particular
5714  // shape function is primitive,
5715  // so we can call
5716  // system_to_component_index
5717  const unsigned int row =
5718  this->finite_element_output
5719  .shape_function_to_row_table[i * fe->n_components() +
5720  fe->system_to_component_index(i).first];
5721  return this->finite_element_output.shape_3rd_derivatives[row][j];
5722  }
5723 }
5724 
5725 
5726 
5727 template <int dim, int spacedim>
5728 inline Tensor<3, spacedim>
5730  const unsigned int i,
5731  const unsigned int j,
5732  const unsigned int component) const
5733 {
5734  AssertIndexRange(i, fe->n_dofs_per_cell());
5735  Assert(this->update_flags & update_3rd_derivatives,
5736  ExcAccessToUninitializedField("update_3rd_derivatives"));
5737  AssertIndexRange(component, fe->n_components());
5738  Assert(present_cell.get() != nullptr,
5739  ExcMessage("FEValues object is not reinit'ed to any cell"));
5740  // check whether the shape function
5741  // is non-zero at all within
5742  // this component:
5743  if (fe->get_nonzero_components(i)[component] == false)
5744  return Tensor<3, spacedim>();
5745 
5746  // look up the right row in the
5747  // table and take the data from
5748  // there
5749  const unsigned int row =
5750  this->finite_element_output
5751  .shape_function_to_row_table[i * fe->n_components() + component];
5752  return this->finite_element_output.shape_3rd_derivatives[row][j];
5753 }
5754 
5755 
5756 
5757 template <int dim, int spacedim>
5758 inline const FiniteElement<dim, spacedim> &
5760 {
5761  return *fe;
5762 }
5763 
5764 
5765 
5766 template <int dim, int spacedim>
5767 inline const Mapping<dim, spacedim> &
5769 {
5770  return *mapping;
5771 }
5772 
5773 
5774 
5775 template <int dim, int spacedim>
5776 inline UpdateFlags
5778 {
5779  return this->update_flags;
5780 }
5781 
5782 
5783 
5784 template <int dim, int spacedim>
5785 inline const std::vector<Point<spacedim>> &
5787 {
5788  Assert(this->update_flags & update_quadrature_points,
5789  ExcAccessToUninitializedField("update_quadrature_points"));
5790  Assert(present_cell.get() != nullptr,
5791  ExcMessage("FEValues object is not reinit'ed to any cell"));
5792  return this->mapping_output.quadrature_points;
5793 }
5794 
5795 
5796 
5797 template <int dim, int spacedim>
5798 inline const std::vector<double> &
5800 {
5801  Assert(this->update_flags & update_JxW_values,
5802  ExcAccessToUninitializedField("update_JxW_values"));
5803  Assert(present_cell.get() != nullptr,
5804  ExcMessage("FEValues object is not reinit'ed to any cell"));
5805  return this->mapping_output.JxW_values;
5806 }
5807 
5808 
5809 
5810 template <int dim, int spacedim>
5811 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5813 {
5814  Assert(this->update_flags & update_jacobians,
5815  ExcAccessToUninitializedField("update_jacobians"));
5816  Assert(present_cell.get() != nullptr,
5817  ExcMessage("FEValues object is not reinit'ed to any cell"));
5818  return this->mapping_output.jacobians;
5819 }
5820 
5821 
5822 
5823 template <int dim, int spacedim>
5824 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5826 {
5827  Assert(this->update_flags & update_jacobian_grads,
5828  ExcAccessToUninitializedField("update_jacobians_grads"));
5829  Assert(present_cell.get() != nullptr,
5830  ExcMessage("FEValues object is not reinit'ed to any cell"));
5831  return this->mapping_output.jacobian_grads;
5832 }
5833 
5834 
5835 
5836 template <int dim, int spacedim>
5837 inline const Tensor<3, spacedim> &
5839  const unsigned int i) const
5840 {
5841  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5842  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5843  Assert(present_cell.get() != nullptr,
5844  ExcMessage("FEValues object is not reinit'ed to any cell"));
5845  return this->mapping_output.jacobian_pushed_forward_grads[i];
5846 }
5847 
5848 
5849 
5850 template <int dim, int spacedim>
5851 inline const std::vector<Tensor<3, spacedim>> &
5853 {
5854  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5855  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5856  Assert(present_cell.get() != nullptr,
5857  ExcMessage("FEValues object is not reinit'ed to any cell"));
5858  return this->mapping_output.jacobian_pushed_forward_grads;
5859 }
5860 
5861 
5862 
5863 template <int dim, int spacedim>
5864 inline const DerivativeForm<3, dim, spacedim> &
5865 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5866 {
5867  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5868  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5869  Assert(present_cell.get() != nullptr,
5870  ExcMessage("FEValues object is not reinit'ed to any cell"));
5871  return this->mapping_output.jacobian_2nd_derivatives[i];
5872 }
5873 
5874 
5875 
5876 template <int dim, int spacedim>
5877 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5879 {
5880  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5881  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5882  Assert(present_cell.get() != nullptr,
5883  ExcMessage("FEValues object is not reinit'ed to any cell"));
5884  return this->mapping_output.jacobian_2nd_derivatives;
5885 }
5886 
5887 
5888 
5889 template <int dim, int spacedim>
5890 inline const Tensor<4, spacedim> &
5892  const unsigned int i) const
5893 {
5896  "update_jacobian_pushed_forward_2nd_derivatives"));
5897  Assert(present_cell.get() != nullptr,
5898  ExcMessage("FEValues object is not reinit'ed to any cell"));
5899  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5900 }
5901 
5902 
5903 
5904 template <int dim, int spacedim>
5905 inline const std::vector<Tensor<4, spacedim>> &
5907 {
5908  Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
5910  "update_jacobian_pushed_forward_2nd_derivatives"));
5911  Assert(present_cell.get() != nullptr,
5912  ExcMessage("FEValues object is not reinit'ed to any cell"));
5913  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5914 }
5915 
5916 
5917 
5918 template <int dim, int spacedim>
5919 inline const DerivativeForm<4, dim, spacedim> &
5920 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5921 {
5922  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5923  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5924  Assert(present_cell.get() != nullptr,
5925  ExcMessage("FEValues object is not reinit'ed to any cell"));
5926  return this->mapping_output.jacobian_3rd_derivatives[i];
5927 }
5928 
5929 
5930 
5931 template <int dim, int spacedim>
5932 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5934 {
5935  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5936  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5937  Assert(present_cell.get() != nullptr,
5938  ExcMessage("FEValues object is not reinit'ed to any cell"));
5939  return this->mapping_output.jacobian_3rd_derivatives;
5940 }
5941 
5942 
5943 
5944 template <int dim, int spacedim>
5945 inline const Tensor<5, spacedim> &
5947  const unsigned int i) const
5948 {
5951  "update_jacobian_pushed_forward_3rd_derivatives"));
5952  Assert(present_cell.get() != nullptr,
5953  ExcMessage("FEValues object is not reinit'ed to any cell"));
5954  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
5955 }
5956 
5957 
5958 
5959 template <int dim, int spacedim>
5960 inline const std::vector<Tensor<5, spacedim>> &
5962 {
5963  Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
5965  "update_jacobian_pushed_forward_3rd_derivatives"));
5966  Assert(present_cell.get() != nullptr,
5967  ExcMessage("FEValues object is not reinit'ed to any cell"));
5968  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
5969 }
5970 
5971 
5972 
5973 template <int dim, int spacedim>
5974 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
5976 {
5977  Assert(this->update_flags & update_inverse_jacobians,
5978  ExcAccessToUninitializedField("update_inverse_jacobians"));
5979  Assert(present_cell.get() != nullptr,
5980  ExcMessage("FEValues object is not reinit'ed to any cell"));
5981  return this->mapping_output.inverse_jacobians;
5982 }
5983 
5984 
5985 
5986 template <int dim, int spacedim>
5989 {
5990  return {0U, dofs_per_cell};
5991 }
5992 
5993 
5994 
5995 template <int dim, int spacedim>
5998  const unsigned int start_dof_index) const
5999 {
6000  Assert(start_dof_index <= dofs_per_cell,
6001  ExcIndexRange(start_dof_index, 0, dofs_per_cell + 1));
6002  return {start_dof_index, dofs_per_cell};
6003 }
6004 
6005 
6006 
6007 template <int dim, int spacedim>
6010  const unsigned int end_dof_index) const
6011 {
6012  Assert(end_dof_index < dofs_per_cell,
6013  ExcIndexRange(end_dof_index, 0, dofs_per_cell));
6014  return {0U, end_dof_index + 1};
6015 }
6016 
6017 
6018 
6019 template <int dim, int spacedim>
6022 {
6023  return {0U, n_quadrature_points};
6024 }
6025 
6026 
6027 
6028 template <int dim, int spacedim>
6029 inline const Point<spacedim> &
6030 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
6031 {
6032  Assert(this->update_flags & update_quadrature_points,
6033  ExcAccessToUninitializedField("update_quadrature_points"));
6034  AssertIndexRange(i, this->mapping_output.quadrature_points.size());
6035  Assert(present_cell.get() != nullptr,
6036  ExcMessage("FEValues object is not reinit'ed to any cell"));
6037 
6038  return this->mapping_output.quadrature_points[i];
6039 }
6040 
6041 
6042 
6043 template <int dim, int spacedim>
6044 inline double
6045 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
6046 {
6047  Assert(this->update_flags & update_JxW_values,
6048  ExcAccessToUninitializedField("update_JxW_values"));
6049  AssertIndexRange(i, this->mapping_output.JxW_values.size());
6050  Assert(present_cell.get() != nullptr,
6051  ExcMessage("FEValues object is not reinit'ed to any cell"));
6052 
6053  return this->mapping_output.JxW_values[i];
6054 }
6055 
6056 
6057 
6058 template <int dim, int spacedim>
6059 inline const DerivativeForm<1, dim, spacedim> &
6060 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
6061 {
6062  Assert(this->update_flags & update_jacobians,
6063  ExcAccessToUninitializedField("update_jacobians"));
6064  AssertIndexRange(i, this->mapping_output.jacobians.size());
6065  Assert(present_cell.get() != nullptr,
6066  ExcMessage("FEValues object is not reinit'ed to any cell"));
6067 
6068  return this->mapping_output.jacobians[i];
6069 }
6070 
6071 
6072 
6073 template <int dim, int spacedim>
6074 inline const DerivativeForm<2, dim, spacedim> &
6075 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
6076 {
6077  Assert(this->update_flags & update_jacobian_grads,
6078  ExcAccessToUninitializedField("update_jacobians_grads"));
6079  AssertIndexRange(i, this->mapping_output.jacobian_grads.size());
6080  Assert(present_cell.get() != nullptr,
6081  ExcMessage("FEValues object is not reinit'ed to any cell"));
6082 
6083  return this->mapping_output.jacobian_grads[i];
6084 }
6085 
6086 
6087 
6088 template <int dim, int spacedim>
6089 inline const DerivativeForm<1, spacedim, dim> &
6090 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
6091 {
6092  Assert(this->update_flags & update_inverse_jacobians,
6093  ExcAccessToUninitializedField("update_inverse_jacobians"));
6094  AssertIndexRange(i, this->mapping_output.inverse_jacobians.size());
6095  Assert(present_cell.get() != nullptr,
6096  ExcMessage("FEValues object is not reinit'ed to any cell"));
6097 
6098  return this->mapping_output.inverse_jacobians[i];
6099 }
6100 
6101 
6102 
6103 template <int dim, int spacedim>
6104 inline const Tensor<1, spacedim> &
6105 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
6106 {
6107  Assert(this->update_flags & update_normal_vectors,
6109  "update_normal_vectors")));
6110  AssertIndexRange(i, this->mapping_output.normal_vectors.size());
6111  Assert(present_cell.get() != nullptr,
6112  ExcMessage("FEValues object is not reinit'ed to any cell"));
6113 
6114  return this->mapping_output.normal_vectors[i];
6115 }
6116 
6117 
6118 
6119 /*--------------------- Inline functions: FEValues --------------------------*/
6120 
6121 
6122 template <int dim, int spacedim>
6123 inline const Quadrature<dim> &
6125 {
6126  return quadrature;
6127 }
6128 
6129 
6130 
6131 template <int dim, int spacedim>
6132 inline const FEValues<dim, spacedim> &
6134 {
6135  return *this;
6136 }
6137 
6138 
6139 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
6140 
6141 
6142 template <int dim, int spacedim>
6143 inline unsigned int
6145 {
6146  return present_face_index;
6147 }
6148 
6149 
6150 /*----------------------- Inline functions: FE*FaceValues -------------------*/
6151 
6152 template <int dim, int spacedim>
6153 inline const Quadrature<dim - 1> &
6155 {
6156  return quadrature[quadrature.size() == 1 ? 0 : present_face_no];
6157 }
6158 
6159 
6160 
6161 template <int dim, int spacedim>
6162 inline const FEFaceValues<dim, spacedim> &
6164 {
6165  return *this;
6166 }
6167 
6168 
6169 
6170 template <int dim, int spacedim>
6171 inline const FESubfaceValues<dim, spacedim> &
6173 {
6174  return *this;
6175 }
6176 
6177 
6178 
6179 template <int dim, int spacedim>
6180 inline const Tensor<1, spacedim> &
6181 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
6182 {
6183  AssertIndexRange(i, this->mapping_output.boundary_forms.size());
6184  Assert(this->update_flags & update_boundary_forms,
6186  "update_boundary_forms")));
6187 
6188  return this->mapping_output.boundary_forms[i];
6189 }
6190 
6191 #endif // DOXYGEN
6192 
6194 
6195 #endif
Transformed quadrature weights.
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
Definition: fe_values.h:811
Shape function values.
typename ProductType< Number, typename Vector< dim, spacedim >::curl_type >::type curl_type
Definition: fe_values.h:868
typename ProductType< Number, divergence_type >::type solution_divergence_type
Definition: fe_values.h:1515
const FEFaceValues< dim, spacedim > & get_present_fe_values() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_starting_at(const unsigned int start_dof_index) const
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
const Tensor< 3, spacedim > & jacobian_pushed_forward_grad(const unsigned int quadrature_point) const
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > fe_data
Definition: fe_values.h:3920
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:3952
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:1895
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1533
unsigned int present_face_no
Definition: fe_values.h:4215
unsigned int present_face_index
Definition: fe_values.h:4221
std::vector<::FEValuesViews::Vector< dim, spacedim > > vectors
Definition: fe_values.h:2286
typename ::internal::FEValuesViews::ViewType< dim, spacedim, Extractor >::type View
Definition: fe_values.h:2308
#define DEAL_II_DEPRECATED_EARLY
Definition: config.h:165
constexpr SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
static ::ExceptionBase & ExcAccessToUninitializedField()
typename internal::ProductTypeImpl< typename std::decay< T >::type, typename std::decay< U >::type >::type type
const unsigned int dofs_per_cell
Definition: fe_values.h:2450
typename ::internal::CurlType< spacedim >::type curl_type
Definition: fe_values.h:719
typename ProductType< Number, value_type >::type solution_value_type
Definition: fe_values.h:1505
Volume element.
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
Definition: fe_values.h:224
typename ::FEValuesViews::SymmetricTensor< rank, dim, spacedim > type
Definition: fe_values.h:2268
#define AssertIndexRange(index, range)
Definition: exceptions.h:1690
const Mapping< dim, spacedim > & get_mapping() const
Outer normal vector, not normalized.
typename ProductType< Number, typename Scalar< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:265
typename ProductType< Number, divergence_type >::type solution_divergence_type
Definition: fe_values.h:1851
const DerivativeForm< 2, dim, spacedim > & jacobian_grad(const unsigned int quadrature_point) const
const FiniteElement< dim, spacedim > & get_fe() const
std::unique_ptr< const CellIteratorBase > present_cell
Definition: fe_values.h:3836
UpdateFlags get_update_flags() const
STL namespace.
static const char U
Transformed quadrature points.
std::vector< double > get_quadrature_points(const unsigned int n)
typename ProductType< Number, value_type >::type solution_value_type
Definition: fe_values.h:1841
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:241
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:3888
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int quadrature_point) const
const DerivativeForm< 3, dim, spacedim > & jacobian_2nd_derivative(const unsigned int quadrature_point) const
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition: fe_values.h:3967
std_cxx20::ranges::iota_view< unsigned int, unsigned int > quadrature_point_indices() const
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
static ::ExceptionBase & ExcFENotPrimitive()
typename ProductType< Number, value_type >::type solution_value_type
Definition: fe_values.h:742
Tensor< 2, spacedim > shape_hessian_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
std::vector<::FEValuesViews::SymmetricTensor< 2, dim, spacedim > > symmetric_second_order_tensors
Definition: fe_values.h:2288
typename ProductType< Number, typename Vector< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:836
const Quadrature< dim - 1 > & get_quadrature() const
const std::vector< Point< spacedim > > & get_quadrature_points() const
const Point< spacedim > & quadrature_point(const unsigned int q) const
const std::vector< Tensor< 4, spacedim > > & get_jacobian_pushed_forward_2nd_derivatives() const
typename ProductType< Number, typename Scalar< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:249
const hp::QCollection< dim - 1 > quadrature
Definition: fe_values.h:4226
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1879
static ::ExceptionBase & ExcMessage(std::string arg1)
typename ProductType< Number, curl_type >::type solution_curl_type
Definition: fe_values.h:791
Number value_type
Definition: vector.h:122
typename ProductType< Number, typename Vector< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:884
static constexpr TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:515
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_ending_at(const unsigned int end_dof_index) const
Tensor< 3, spacedim > shape_3rd_derivative_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
typename ProductType< Number, hessian_type >::type solution_hessian_type
Definition: fe_values.h:801
Third derivatives of shape functions.
const FEValues< dim, spacedim > & get_present_fe_values() const
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:2200
unsigned int get_face_index() const
#define Assert(cond, exc)
Definition: exceptions.h:1465
UpdateFlags
const Tensor< 2, spacedim > & shape_hessian(const unsigned int function_no, const unsigned int point_no) const
Abstract base class for mapping classes.
Definition: mapping.h:303
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1446
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
const Quadrature< dim > quadrature
Definition: fe_values.h:4103
const unsigned int first_vector_component
Definition: fe_values.h:1441
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:395
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition: fe_values.h:3896
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int quadrature_point) const
typename ProductType< Number, value_type >::type solution_laplacian_type
Definition: fe_values.h:204
typename Tensor< rank_ - 1, dim, Number >::tensor_type value_type
Definition: tensor.h:496
typename ProductType< Number, gradient_type >::type solution_gradient_type
Definition: fe_values.h:194
const std::vector< DerivativeForm< 4, dim, spacedim > > & get_jacobian_3rd_derivatives() const
Second derivatives of shape functions.
Gradient of volume element.
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1541
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
std::vector<::FEValuesViews::Scalar< dim, spacedim > > scalars
Definition: fe_values.h:2285
const Tensor< 1, spacedim > & boundary_form(const unsigned int i) const
const Quadrature< dim > & get_quadrature() const
const unsigned int n_quadrature_points
Definition: fe_values.h:2432
typename ProductType< Number, hessian_type >::type solution_hessian_type
Definition: fe_values.h:214
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices() const
Tensor< 1, spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
typename ProductType< Number, value_type >::type solution_value_type
Definition: fe_values.h:184
typename ProductType< Number, symmetric_gradient_type >::type solution_symmetric_gradient_type
Definition: fe_values.h:762
boost::signals2::connection tria_listener_mesh_transform
Definition: fe_values.h:3861
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:828
typename ::FEValuesViews::Tensor< rank, dim, spacedim > type
Definition: fe_values.h:2254
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
Definition: tensor.h:462
const std::vector< double > & get_JxW_values() const
double JxW(const unsigned int quadrature_point) const
typename ProductType< Number, typename Vector< dim, spacedim >::symmetric_gradient_type >::type symmetric_gradient_type
Definition: fe_values.h:844
typename ProductType< Number, typename Vector< dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:852
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:394
typename ProductType< Number, typename Vector< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:876
Shape function gradients.
Normal vectors.
const std::vector< DerivativeForm< 2, dim, spacedim > > & get_jacobian_grads() const
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1887
const Tensor< 4, spacedim > & jacobian_pushed_forward_2nd_derivative(const unsigned int quadrature_point) const
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:2189
Definition: fe.h:38
typename ProductType< Number, value_type >::type solution_laplacian_type
Definition: fe_values.h:782
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int quadrature_point) const
static ::ExceptionBase & ExcNotImplemented()
const std::vector< DerivativeForm< 3, dim, spacedim > > & get_jacobian_2nd_derivatives() const
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1435
const Tensor< 1, spacedim > & shape_grad(const unsigned int function_no, const unsigned int quadrature_point) const
boost::signals2::connection tria_listener_refinement
Definition: fe_values.h:3852
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3928
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:860
typename ProductType< Number, gradient_type >::type solution_gradient_type
Definition: fe_values.h:752
::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:3903
const std::vector< Tensor< 5, spacedim > > & get_jacobian_pushed_forward_3rd_derivatives() const
const FESubfaceValues< dim, spacedim > & get_present_fe_values() const
const unsigned int max_n_quadrature_points
Definition: fe_values.h:2443
std::vector<::FEValuesViews::Tensor< 2, dim, spacedim > > second_order_tensors
Definition: fe_values.h:2290
const std::vector< DerivativeForm< 1, spacedim, dim > > & get_inverse_jacobians() const
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:257
typename ProductType< Number, typename Scalar< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:273
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
UpdateFlags update_flags
Definition: fe_values.h:3934
const Tensor< 3, spacedim > & shape_3rd_derivative(const unsigned int function_no, const unsigned int point_no) const
static ::ExceptionBase & ExcInternalError()
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
Definition: fe_values.h:3912
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
typename ProductType< Number, divergence_type >::type solution_divergence_type
Definition: fe_values.h:772
typename ProductType< Number, gradient_type >::type solution_gradient_type
Definition: fe_values.h:1861