Reference documentation for deal.II version Git 3f1f337db3 2021-10-23 13:19:02 -0600
\(\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 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 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  Tensor &
1988  operator=(Tensor<2, dim, spacedim> &&) = default; // NOLINT
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_number() const;
4195 
4200  unsigned int
4201  get_face_index() const;
4202 
4207  const Quadrature<dim - 1> &
4208  get_quadrature() const;
4209 
4214  std::size_t
4215  memory_consumption() const;
4216 
4217 protected:
4222  unsigned int present_face_no;
4223 
4228  unsigned int present_face_index;
4229 
4233  const hp::QCollection<dim - 1> quadrature;
4234 };
4235 
4236 
4237 
4251 template <int dim, int spacedim = dim>
4252 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
4253 {
4254 public:
4259  static const unsigned int dimension = dim;
4260 
4261  static const unsigned int space_dimension = spacedim;
4262 
4267  static const unsigned int integral_dimension = dim - 1;
4268 
4273  FEFaceValues(const Mapping<dim, spacedim> & mapping,
4274  const FiniteElement<dim, spacedim> &fe,
4275  const Quadrature<dim - 1> & quadrature,
4276  const UpdateFlags update_flags);
4277 
4284  FEFaceValues(const Mapping<dim, spacedim> & mapping,
4285  const FiniteElement<dim, spacedim> &fe,
4286  const hp::QCollection<dim - 1> & quadrature,
4287  const UpdateFlags update_flags);
4288 
4295  const Quadrature<dim - 1> & quadrature,
4296  const UpdateFlags update_flags);
4297 
4305  const hp::QCollection<dim - 1> & quadrature,
4306  const UpdateFlags update_flags);
4307 
4312  template <bool level_dof_access>
4313  void
4314  reinit(
4316  const unsigned int face_no);
4317 
4324  template <bool level_dof_access>
4325  void
4326  reinit(
4328  const typename Triangulation<dim, spacedim>::face_iterator & face);
4329 
4343  void
4345  const unsigned int face_no);
4346 
4347  /*
4348  * Reinitialize the gradients, Jacobi determinants, etc for the given face
4349  * on a given cell of type "iterator into a Triangulation object", and the
4350  * given finite element. Since iterators into a triangulation alone only
4351  * convey information about the geometry of a cell, but not about degrees of
4352  * freedom possibly associated with this cell, you will not be able to call
4353  * some functions of this class if they need information about degrees of
4354  * freedom. These functions are, above all, the
4355  * <tt>get_function_value/gradients/hessians/third_derivatives</tt>
4356  * functions. If you want to call these functions, you have to call the @p
4357  * reinit variants that take iterators into DoFHandler or other DoF handler
4358  * type objects.
4359  *
4360  * @note @p face must be one of @p cell's face iterators.
4361  */
4362  void
4364  const typename Triangulation<dim, spacedim>::face_iterator &face);
4365 
4381  get_present_fe_values() const;
4382 
4383 private:
4387  void
4388  initialize(const UpdateFlags update_flags);
4389 
4396  void
4397  do_reinit(const unsigned int face_no);
4398 };
4399 
4400 
4417 template <int dim, int spacedim = dim>
4418 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
4419 {
4420 public:
4424  static const unsigned int dimension = dim;
4425 
4429  static const unsigned int space_dimension = spacedim;
4430 
4435  static const unsigned int integral_dimension = dim - 1;
4436 
4441  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
4442  const FiniteElement<dim, spacedim> &fe,
4443  const Quadrature<dim - 1> & face_quadrature,
4444  const UpdateFlags update_flags);
4445 
4452  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
4453  const FiniteElement<dim, spacedim> &fe,
4454  const hp::QCollection<dim - 1> & face_quadrature,
4455  const UpdateFlags update_flags);
4456 
4463  const Quadrature<dim - 1> & face_quadrature,
4464  const UpdateFlags update_flags);
4465 
4473  const hp::QCollection<dim - 1> & face_quadrature,
4474  const UpdateFlags update_flags);
4475 
4482  template <bool level_dof_access>
4483  void
4484  reinit(
4486  const unsigned int face_no,
4487  const unsigned int subface_no);
4488 
4493  template <bool level_dof_access>
4494  void
4495  reinit(
4497  const typename Triangulation<dim, spacedim>::face_iterator & face,
4498  const typename Triangulation<dim, spacedim>::face_iterator &subface);
4499 
4513  void
4515  const unsigned int face_no,
4516  const unsigned int subface_no);
4517 
4537  void
4539  const typename Triangulation<dim, spacedim>::face_iterator &face,
4540  const typename Triangulation<dim, spacedim>::face_iterator &subface);
4541 
4557  get_present_fe_values() const;
4558 
4564  DeclException0(ExcReinitCalledWithBoundaryFace);
4565 
4571  DeclException0(ExcFaceHasNoSubfaces);
4572 
4573 private:
4577  void
4578  initialize(const UpdateFlags update_flags);
4579 
4586  void
4587  do_reinit(const unsigned int face_no, const unsigned int subface_no);
4588 };
4589 
4590 
4591 #ifndef DOXYGEN
4592 
4593 
4594 /*------------------------ Inline functions: namespace FEValuesViews --------*/
4595 
4596 namespace FEValuesViews
4597 {
4598  template <int dim, int spacedim>
4599  inline typename Scalar<dim, spacedim>::value_type
4600  Scalar<dim, spacedim>::value(const unsigned int shape_function,
4601  const unsigned int q_point) const
4602  {
4603  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4604  Assert(
4605  fe_values->update_flags & update_values,
4607  "update_values"))));
4608 
4609  // an adaptation of the FEValuesBase::shape_value_component function
4610  // except that here we know the component as fixed and we have
4611  // pre-computed and cached a bunch of information. See the comments there.
4612  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4613  return fe_values->finite_element_output.shape_values(
4614  shape_function_data[shape_function].row_index, q_point);
4615  else
4616  return 0;
4617  }
4618 
4619 
4620 
4621  template <int dim, int spacedim>
4622  inline typename Scalar<dim, spacedim>::gradient_type
4623  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
4624  const unsigned int q_point) const
4625  {
4626  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4627  Assert(fe_values->update_flags & update_gradients,
4629  "update_gradients")));
4630 
4631  // an adaptation of the FEValuesBase::shape_grad_component
4632  // function except that here we know the component as fixed and we have
4633  // pre-computed and cached a bunch of information. See the comments there.
4634  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4635  return fe_values->finite_element_output
4636  .shape_gradients[shape_function_data[shape_function].row_index]
4637  [q_point];
4638  else
4639  return gradient_type();
4640  }
4641 
4642 
4643 
4644  template <int dim, int spacedim>
4645  inline typename Scalar<dim, spacedim>::hessian_type
4646  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
4647  const unsigned int q_point) const
4648  {
4649  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4650  Assert(fe_values->update_flags & update_hessians,
4652  "update_hessians")));
4653 
4654  // an adaptation of the FEValuesBase::shape_hessian_component
4655  // function except that here we know the component as fixed and we have
4656  // pre-computed and cached a bunch of information. See the comments there.
4657  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4658  return fe_values->finite_element_output
4659  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
4660  else
4661  return hessian_type();
4662  }
4663 
4664 
4665 
4666  template <int dim, int spacedim>
4667  inline typename Scalar<dim, spacedim>::third_derivative_type
4668  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
4669  const unsigned int q_point) const
4670  {
4671  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4672  Assert(fe_values->update_flags & update_3rd_derivatives,
4674  "update_3rd_derivatives")));
4675 
4676  // an adaptation of the FEValuesBase::shape_3rdderivative_component
4677  // function except that here we know the component as fixed and we have
4678  // pre-computed and cached a bunch of information. See the comments there.
4679  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4680  return fe_values->finite_element_output
4681  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
4682  [q_point];
4683  else
4684  return third_derivative_type();
4685  }
4686 
4687 
4688 
4689  template <int dim, int spacedim>
4690  inline typename Vector<dim, spacedim>::value_type
4691  Vector<dim, spacedim>::value(const unsigned int shape_function,
4692  const unsigned int q_point) const
4693  {
4694  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4695  Assert(fe_values->update_flags & update_values,
4697  "update_values")));
4698 
4699  // same as for the scalar case except that we have one more index
4700  const int snc =
4701  shape_function_data[shape_function].single_nonzero_component;
4702  if (snc == -2)
4703  return value_type();
4704  else if (snc != -1)
4705  {
4706  value_type return_value;
4707  return_value[shape_function_data[shape_function]
4708  .single_nonzero_component_index] =
4709  fe_values->finite_element_output.shape_values(snc, q_point);
4710  return return_value;
4711  }
4712  else
4713  {
4714  value_type return_value;
4715  for (unsigned int d = 0; d < dim; ++d)
4716  if (shape_function_data[shape_function]
4717  .is_nonzero_shape_function_component[d])
4718  return_value[d] = fe_values->finite_element_output.shape_values(
4719  shape_function_data[shape_function].row_index[d], q_point);
4720 
4721  return return_value;
4722  }
4723  }
4724 
4725 
4726 
4727  template <int dim, int spacedim>
4728  inline typename Vector<dim, spacedim>::gradient_type
4729  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4730  const unsigned int q_point) const
4731  {
4732  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4733  Assert(fe_values->update_flags & update_gradients,
4735  "update_gradients")));
4736 
4737  // same as for the scalar case except that we have one more index
4738  const int snc =
4739  shape_function_data[shape_function].single_nonzero_component;
4740  if (snc == -2)
4741  return gradient_type();
4742  else if (snc != -1)
4743  {
4744  gradient_type return_value;
4745  return_value[shape_function_data[shape_function]
4746  .single_nonzero_component_index] =
4747  fe_values->finite_element_output.shape_gradients[snc][q_point];
4748  return return_value;
4749  }
4750  else
4751  {
4752  gradient_type return_value;
4753  for (unsigned int d = 0; d < dim; ++d)
4754  if (shape_function_data[shape_function]
4755  .is_nonzero_shape_function_component[d])
4756  return_value[d] =
4757  fe_values->finite_element_output.shape_gradients
4758  [shape_function_data[shape_function].row_index[d]][q_point];
4759 
4760  return return_value;
4761  }
4762  }
4763 
4764 
4765 
4766  template <int dim, int spacedim>
4767  inline typename Vector<dim, spacedim>::divergence_type
4768  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4769  const unsigned int q_point) const
4770  {
4771  // this function works like in the case above
4772  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4773  Assert(fe_values->update_flags & update_gradients,
4775  "update_gradients")));
4776 
4777  // same as for the scalar case except that we have one more index
4778  const int snc =
4779  shape_function_data[shape_function].single_nonzero_component;
4780  if (snc == -2)
4781  return divergence_type();
4782  else if (snc != -1)
4783  return fe_values->finite_element_output
4784  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4785  .single_nonzero_component_index];
4786  else
4787  {
4788  divergence_type return_value = 0;
4789  for (unsigned int d = 0; d < dim; ++d)
4790  if (shape_function_data[shape_function]
4791  .is_nonzero_shape_function_component[d])
4792  return_value +=
4793  fe_values->finite_element_output.shape_gradients
4794  [shape_function_data[shape_function].row_index[d]][q_point][d];
4795 
4796  return return_value;
4797  }
4798  }
4799 
4800 
4801 
4802  template <int dim, int spacedim>
4803  inline typename Vector<dim, spacedim>::curl_type
4804  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4805  const unsigned int q_point) const
4806  {
4807  // this function works like in the case above
4808 
4809  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4810  Assert(fe_values->update_flags & update_gradients,
4812  "update_gradients")));
4813  // same as for the scalar case except that we have one more index
4814  const int snc =
4815  shape_function_data[shape_function].single_nonzero_component;
4816 
4817  if (snc == -2)
4818  return curl_type();
4819 
4820  else
4821  switch (dim)
4822  {
4823  case 1:
4824  {
4825  Assert(false,
4826  ExcMessage(
4827  "Computing the curl in 1d is not a useful operation"));
4828  return curl_type();
4829  }
4830 
4831  case 2:
4832  {
4833  if (snc != -1)
4834  {
4835  curl_type return_value;
4836 
4837  // the single nonzero component can only be zero or one in 2d
4838  if (shape_function_data[shape_function]
4839  .single_nonzero_component_index == 0)
4840  return_value[0] =
4841  -1.0 * fe_values->finite_element_output
4842  .shape_gradients[snc][q_point][1];
4843  else
4844  return_value[0] = fe_values->finite_element_output
4845  .shape_gradients[snc][q_point][0];
4846 
4847  return return_value;
4848  }
4849 
4850  else
4851  {
4852  curl_type return_value;
4853 
4854  return_value[0] = 0.0;
4855 
4856  if (shape_function_data[shape_function]
4857  .is_nonzero_shape_function_component[0])
4858  return_value[0] -=
4859  fe_values->finite_element_output
4860  .shape_gradients[shape_function_data[shape_function]
4861  .row_index[0]][q_point][1];
4862 
4863  if (shape_function_data[shape_function]
4864  .is_nonzero_shape_function_component[1])
4865  return_value[0] +=
4866  fe_values->finite_element_output
4867  .shape_gradients[shape_function_data[shape_function]
4868  .row_index[1]][q_point][0];
4869 
4870  return return_value;
4871  }
4872  }
4873 
4874  case 3:
4875  {
4876  if (snc != -1)
4877  {
4878  curl_type return_value;
4879 
4880  switch (shape_function_data[shape_function]
4881  .single_nonzero_component_index)
4882  {
4883  case 0:
4884  {
4885  return_value[0] = 0;
4886  return_value[1] = fe_values->finite_element_output
4887  .shape_gradients[snc][q_point][2];
4888  return_value[2] =
4889  -1.0 * fe_values->finite_element_output
4890  .shape_gradients[snc][q_point][1];
4891  return return_value;
4892  }
4893 
4894  case 1:
4895  {
4896  return_value[0] =
4897  -1.0 * fe_values->finite_element_output
4898  .shape_gradients[snc][q_point][2];
4899  return_value[1] = 0;
4900  return_value[2] = fe_values->finite_element_output
4901  .shape_gradients[snc][q_point][0];
4902  return return_value;
4903  }
4904 
4905  default:
4906  {
4907  return_value[0] = fe_values->finite_element_output
4908  .shape_gradients[snc][q_point][1];
4909  return_value[1] =
4910  -1.0 * fe_values->finite_element_output
4911  .shape_gradients[snc][q_point][0];
4912  return_value[2] = 0;
4913  return return_value;
4914  }
4915  }
4916  }
4917 
4918  else
4919  {
4920  curl_type return_value;
4921 
4922  for (unsigned int i = 0; i < dim; ++i)
4923  return_value[i] = 0.0;
4924 
4925  if (shape_function_data[shape_function]
4926  .is_nonzero_shape_function_component[0])
4927  {
4928  return_value[1] +=
4929  fe_values->finite_element_output
4930  .shape_gradients[shape_function_data[shape_function]
4931  .row_index[0]][q_point][2];
4932  return_value[2] -=
4933  fe_values->finite_element_output
4934  .shape_gradients[shape_function_data[shape_function]
4935  .row_index[0]][q_point][1];
4936  }
4937 
4938  if (shape_function_data[shape_function]
4939  .is_nonzero_shape_function_component[1])
4940  {
4941  return_value[0] -=
4942  fe_values->finite_element_output
4943  .shape_gradients[shape_function_data[shape_function]
4944  .row_index[1]][q_point][2];
4945  return_value[2] +=
4946  fe_values->finite_element_output
4947  .shape_gradients[shape_function_data[shape_function]
4948  .row_index[1]][q_point][0];
4949  }
4950 
4951  if (shape_function_data[shape_function]
4952  .is_nonzero_shape_function_component[2])
4953  {
4954  return_value[0] +=
4955  fe_values->finite_element_output
4956  .shape_gradients[shape_function_data[shape_function]
4957  .row_index[2]][q_point][1];
4958  return_value[1] -=
4959  fe_values->finite_element_output
4960  .shape_gradients[shape_function_data[shape_function]
4961  .row_index[2]][q_point][0];
4962  }
4963 
4964  return return_value;
4965  }
4966  }
4967  }
4968  // should not end up here
4969  Assert(false, ExcInternalError());
4970  return curl_type();
4971  }
4972 
4973 
4974 
4975  template <int dim, int spacedim>
4976  inline typename Vector<dim, spacedim>::hessian_type
4977  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
4978  const unsigned int q_point) const
4979  {
4980  // this function works like in the case above
4981  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4982  Assert(fe_values->update_flags & update_hessians,
4984  "update_hessians")));
4985 
4986  // same as for the scalar case except that we have one more index
4987  const int snc =
4988  shape_function_data[shape_function].single_nonzero_component;
4989  if (snc == -2)
4990  return hessian_type();
4991  else if (snc != -1)
4992  {
4993  hessian_type return_value;
4994  return_value[shape_function_data[shape_function]
4995  .single_nonzero_component_index] =
4996  fe_values->finite_element_output.shape_hessians[snc][q_point];
4997  return return_value;
4998  }
4999  else
5000  {
5001  hessian_type return_value;
5002  for (unsigned int d = 0; d < dim; ++d)
5003  if (shape_function_data[shape_function]
5004  .is_nonzero_shape_function_component[d])
5005  return_value[d] =
5006  fe_values->finite_element_output.shape_hessians
5007  [shape_function_data[shape_function].row_index[d]][q_point];
5008 
5009  return return_value;
5010  }
5011  }
5012 
5013 
5014 
5015  template <int dim, int spacedim>
5016  inline typename Vector<dim, spacedim>::third_derivative_type
5017  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
5018  const unsigned int q_point) const
5019  {
5020  // this function works like in the case above
5021  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5022  Assert(fe_values->update_flags & update_3rd_derivatives,
5024  "update_3rd_derivatives")));
5025 
5026  // same as for the scalar case except that we have one more index
5027  const int snc =
5028  shape_function_data[shape_function].single_nonzero_component;
5029  if (snc == -2)
5030  return third_derivative_type();
5031  else if (snc != -1)
5032  {
5033  third_derivative_type return_value;
5034  return_value[shape_function_data[shape_function]
5035  .single_nonzero_component_index] =
5036  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
5037  return return_value;
5038  }
5039  else
5040  {
5041  third_derivative_type return_value;
5042  for (unsigned int d = 0; d < dim; ++d)
5043  if (shape_function_data[shape_function]
5044  .is_nonzero_shape_function_component[d])
5045  return_value[d] =
5046  fe_values->finite_element_output.shape_3rd_derivatives
5047  [shape_function_data[shape_function].row_index[d]][q_point];
5048 
5049  return return_value;
5050  }
5051  }
5052 
5053 
5054 
5055  namespace internal
5056  {
5061  inline ::SymmetricTensor<2, 1>
5062  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
5063  {
5064  AssertIndexRange(n, 1);
5065  (void)n;
5066 
5067  return {{t[0]}};
5068  }
5069 
5070 
5071 
5072  inline ::SymmetricTensor<2, 2>
5073  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
5074  {
5075  switch (n)
5076  {
5077  case 0:
5078  {
5079  return {{t[0], 0, t[1] / 2}};
5080  }
5081  case 1:
5082  {
5083  return {{0, t[1], t[0] / 2}};
5084  }
5085  default:
5086  {
5087  AssertIndexRange(n, 2);
5088  return {};
5089  }
5090  }
5091  }
5092 
5093 
5094 
5095  inline ::SymmetricTensor<2, 3>
5096  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
5097  {
5098  switch (n)
5099  {
5100  case 0:
5101  {
5102  return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
5103  }
5104  case 1:
5105  {
5106  return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
5107  }
5108  case 2:
5109  {
5110  return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
5111  }
5112  default:
5113  {
5114  AssertIndexRange(n, 3);
5115  return {};
5116  }
5117  }
5118  }
5119  } // namespace internal
5120 
5121 
5122 
5123  template <int dim, int spacedim>
5124  inline typename Vector<dim, spacedim>::symmetric_gradient_type
5125  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
5126  const unsigned int q_point) const
5127  {
5128  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5129  Assert(fe_values->update_flags & update_gradients,
5131  "update_gradients")));
5132 
5133  // same as for the scalar case except that we have one more index
5134  const int snc =
5135  shape_function_data[shape_function].single_nonzero_component;
5136  if (snc == -2)
5137  return symmetric_gradient_type();
5138  else if (snc != -1)
5139  return internal::symmetrize_single_row(
5140  shape_function_data[shape_function].single_nonzero_component_index,
5141  fe_values->finite_element_output.shape_gradients[snc][q_point]);
5142  else
5143  {
5144  gradient_type return_value;
5145  for (unsigned int d = 0; d < dim; ++d)
5146  if (shape_function_data[shape_function]
5147  .is_nonzero_shape_function_component[d])
5148  return_value[d] =
5149  fe_values->finite_element_output.shape_gradients
5150  [shape_function_data[shape_function].row_index[d]][q_point];
5151 
5152  return symmetrize(return_value);
5153  }
5154  }
5155 
5156 
5157 
5158  template <int dim, int spacedim>
5160  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
5161  const unsigned int q_point) const
5162  {
5163  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5164  Assert(fe_values->update_flags & update_values,
5166  "update_values")));
5167 
5168  // similar to the vector case where we have more then one index and we need
5169  // to convert between unrolled and component indexing for tensors
5170  const int snc =
5171  shape_function_data[shape_function].single_nonzero_component;
5172 
5173  if (snc == -2)
5174  {
5175  // shape function is zero for the selected components
5176  return value_type();
5177  }
5178  else if (snc != -1)
5179  {
5180  value_type return_value;
5181  const unsigned int comp =
5182  shape_function_data[shape_function].single_nonzero_component_index;
5183  return_value[value_type::unrolled_to_component_indices(comp)] =
5184  fe_values->finite_element_output.shape_values(snc, q_point);
5185  return return_value;
5186  }
5187  else
5188  {
5189  value_type return_value;
5190  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
5191  if (shape_function_data[shape_function]
5192  .is_nonzero_shape_function_component[d])
5193  return_value[value_type::unrolled_to_component_indices(d)] =
5194  fe_values->finite_element_output.shape_values(
5195  shape_function_data[shape_function].row_index[d], q_point);
5196  return return_value;
5197  }
5198  }
5199 
5200 
5201 
5202  template <int dim, int spacedim>
5205  const unsigned int shape_function,
5206  const unsigned int q_point) const
5207  {
5208  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5209  Assert(fe_values->update_flags & update_gradients,
5211  "update_gradients")));
5212 
5213  const int snc =
5214  shape_function_data[shape_function].single_nonzero_component;
5215 
5216  if (snc == -2)
5217  {
5218  // shape function is zero for the selected components
5219  return divergence_type();
5220  }
5221  else if (snc != -1)
5222  {
5223  // we have a single non-zero component when the symmetric tensor is
5224  // represented in unrolled form. this implies we potentially have
5225  // two non-zero components when represented in component form! we
5226  // will only have one non-zero entry if the non-zero component lies on
5227  // the diagonal of the tensor.
5228  //
5229  // the divergence of a second-order tensor is a first order tensor.
5230  //
5231  // assume the second-order tensor is A with components A_{ij}. then
5232  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
5233  // entries in the tensorial representation. define the
5234  // divergence as:
5235  // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
5236  // (which is incidentally also
5237  // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
5238  // In both cases, a sum is implied.
5239  //
5240  // Now, we know the nonzero component in unrolled form: it is indicated
5241  // by 'snc'. we can figure out which tensor components belong to this:
5242  const unsigned int comp =
5243  shape_function_data[shape_function].single_nonzero_component_index;
5244  const unsigned int ii =
5245  value_type::unrolled_to_component_indices(comp)[0];
5246  const unsigned int jj =
5247  value_type::unrolled_to_component_indices(comp)[1];
5248 
5249  // given the form of the divergence above, if ii=jj there is only a
5250  // single nonzero component of the full tensor and the gradient
5251  // equals
5252  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
5253  // all other entries of 'b' are zero
5254  //
5255  // on the other hand, if ii!=jj, then there are two nonzero entries in
5256  // the full tensor and
5257  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
5258  // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
5259  // again, all other entries of 'b' are zero
5260  const ::Tensor<1, spacedim> &phi_grad =
5261  fe_values->finite_element_output.shape_gradients[snc][q_point];
5262 
5263  divergence_type return_value;
5264  return_value[ii] = phi_grad[jj];
5265 
5266  if (ii != jj)
5267  return_value[jj] = phi_grad[ii];
5268 
5269  return return_value;
5270  }
5271  else
5272  {
5273  Assert(false, ExcNotImplemented());
5274  divergence_type return_value;
5275  return return_value;
5276  }
5277  }
5278 
5279 
5280 
5281  template <int dim, int spacedim>
5282  inline typename Tensor<2, dim, spacedim>::value_type
5283  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
5284  const unsigned int q_point) const
5285  {
5286  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5287  Assert(fe_values->update_flags & update_values,
5289  "update_values")));
5290 
5291  // similar to the vector case where we have more then one index and we need
5292  // to convert between unrolled and component indexing for tensors
5293  const int snc =
5294  shape_function_data[shape_function].single_nonzero_component;
5295 
5296  if (snc == -2)
5297  {
5298  // shape function is zero for the selected components
5299  return value_type();
5300  }
5301  else if (snc != -1)
5302  {
5303  value_type return_value;
5304  const unsigned int comp =
5305  shape_function_data[shape_function].single_nonzero_component_index;
5306  const TableIndices<2> indices =
5308  return_value[indices] =
5309  fe_values->finite_element_output.shape_values(snc, q_point);
5310  return return_value;
5311  }
5312  else
5313  {
5314  value_type return_value;
5315  for (unsigned int d = 0; d < dim * dim; ++d)
5316  if (shape_function_data[shape_function]
5317  .is_nonzero_shape_function_component[d])
5318  {
5319  const TableIndices<2> indices =
5321  return_value[indices] =
5322  fe_values->finite_element_output.shape_values(
5323  shape_function_data[shape_function].row_index[d], q_point);
5324  }
5325  return return_value;
5326  }
5327  }
5328 
5329 
5330 
5331  template <int dim, int spacedim>
5333  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
5334  const unsigned int q_point) const
5335  {
5336  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5337  Assert(fe_values->update_flags & update_gradients,
5339  "update_gradients")));
5340 
5341  const int snc =
5342  shape_function_data[shape_function].single_nonzero_component;
5343 
5344  if (snc == -2)
5345  {
5346  // shape function is zero for the selected components
5347  return divergence_type();
5348  }
5349  else if (snc != -1)
5350  {
5351  // we have a single non-zero component when the tensor is
5352  // represented in unrolled form.
5353  //
5354  // the divergence of a second-order tensor is a first order tensor.
5355  //
5356  // assume the second-order tensor is A with components A_{ij},
5357  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
5358  //
5359  // Now, we know the nonzero component in unrolled form: it is indicated
5360  // by 'snc'. we can figure out which tensor components belong to this:
5361  const unsigned int comp =
5362  shape_function_data[shape_function].single_nonzero_component_index;
5363  const TableIndices<2> indices =
5365  const unsigned int ii = indices[0];
5366  const unsigned int jj = indices[1];
5367 
5368  const ::Tensor<1, spacedim> &phi_grad =
5369  fe_values->finite_element_output.shape_gradients[snc][q_point];
5370 
5371  divergence_type return_value;
5372  // note that we contract \nabla from the right
5373  return_value[ii] = phi_grad[jj];
5374 
5375  return return_value;
5376  }
5377  else
5378  {
5379  Assert(false, ExcNotImplemented());
5380  divergence_type return_value;
5381  return return_value;
5382  }
5383  }
5384 
5385 
5386 
5387  template <int dim, int spacedim>
5389  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
5390  const unsigned int q_point) const
5391  {
5392  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5393  Assert(fe_values->update_flags & update_gradients,
5395  "update_gradients")));
5396 
5397  const int snc =
5398  shape_function_data[shape_function].single_nonzero_component;
5399 
5400  if (snc == -2)
5401  {
5402  // shape function is zero for the selected components
5403  return gradient_type();
5404  }
5405  else if (snc != -1)
5406  {
5407  // we have a single non-zero component when the tensor is
5408  // represented in unrolled form.
5409  //
5410  // the gradient of a second-order tensor is a third order tensor.
5411  //
5412  // assume the second-order tensor is A with components A_{ij},
5413  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
5414  //
5415  // Now, we know the nonzero component in unrolled form: it is indicated
5416  // by 'snc'. we can figure out which tensor components belong to this:
5417  const unsigned int comp =
5418  shape_function_data[shape_function].single_nonzero_component_index;
5419  const TableIndices<2> indices =
5421  const unsigned int ii = indices[0];
5422  const unsigned int jj = indices[1];
5423 
5424  const ::Tensor<1, spacedim> &phi_grad =
5425  fe_values->finite_element_output.shape_gradients[snc][q_point];
5426 
5427  gradient_type return_value;
5428  return_value[ii][jj] = phi_grad;
5429 
5430  return return_value;
5431  }
5432  else
5433  {
5434  Assert(false, ExcNotImplemented());
5435  gradient_type return_value;
5436  return return_value;
5437  }
5438  }
5439 
5440 } // namespace FEValuesViews
5441 
5442 
5443 
5444 /*---------------------- Inline functions: FEValuesBase ---------------------*/
5445 
5446 
5447 
5448 template <int dim, int spacedim>
5451  const FEValuesExtractors::Scalar &scalar) const
5452 {
5453  AssertIndexRange(scalar.component, fe_values_views_cache.scalars.size());
5454 
5455  return fe_values_views_cache.scalars[scalar.component];
5456 }
5457 
5458 
5459 
5460 template <int dim, int spacedim>
5463  const FEValuesExtractors::Vector &vector) const
5464 {
5466  fe_values_views_cache.vectors.size());
5467 
5468  return fe_values_views_cache.vectors[vector.first_vector_component];
5469 }
5470 
5471 
5472 
5473 template <int dim, int spacedim>
5476  const FEValuesExtractors::SymmetricTensor<2> &tensor) const
5477 {
5478  Assert(
5479  tensor.first_tensor_component <
5480  fe_values_views_cache.symmetric_second_order_tensors.size(),
5482  0,
5483  fe_values_views_cache.symmetric_second_order_tensors.size()));
5484 
5485  return fe_values_views_cache
5486  .symmetric_second_order_tensors[tensor.first_tensor_component];
5487 }
5488 
5489 
5490 
5491 template <int dim, int spacedim>
5494  const FEValuesExtractors::Tensor<2> &tensor) const
5495 {
5497  fe_values_views_cache.second_order_tensors.size());
5498 
5499  return fe_values_views_cache
5500  .second_order_tensors[tensor.first_tensor_component];
5501 }
5502 
5503 
5504 
5505 template <int dim, int spacedim>
5506 inline const double &
5507 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
5508  const unsigned int j) const
5509 {
5510  AssertIndexRange(i, fe->n_dofs_per_cell());
5511  Assert(this->update_flags & update_values,
5512  ExcAccessToUninitializedField("update_values"));
5513  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5514  Assert(present_cell.get() != nullptr,
5515  ExcMessage("FEValues object is not reinit'ed to any cell"));
5516  // if the entire FE is primitive,
5517  // then we can take a short-cut:
5518  if (fe->is_primitive())
5519  return this->finite_element_output.shape_values(i, j);
5520  else
5521  {
5522  // otherwise, use the mapping
5523  // between shape function
5524  // numbers and rows. note that
5525  // by the assertions above, we
5526  // know that this particular
5527  // shape function is primitive,
5528  // so we can call
5529  // system_to_component_index
5530  const unsigned int row =
5531  this->finite_element_output
5532  .shape_function_to_row_table[i * fe->n_components() +
5533  fe->system_to_component_index(i).first];
5534  return this->finite_element_output.shape_values(row, j);
5535  }
5536 }
5537 
5538 
5539 
5540 template <int dim, int spacedim>
5541 inline double
5543  const unsigned int i,
5544  const unsigned int j,
5545  const unsigned int component) const
5546 {
5547  AssertIndexRange(i, fe->n_dofs_per_cell());
5548  Assert(this->update_flags & update_values,
5549  ExcAccessToUninitializedField("update_values"));
5550  AssertIndexRange(component, fe->n_components());
5551  Assert(present_cell.get() != nullptr,
5552  ExcMessage("FEValues object is not reinit'ed to any cell"));
5553 
5554  // check whether the shape function
5555  // is non-zero at all within
5556  // this component:
5557  if (fe->get_nonzero_components(i)[component] == false)
5558  return 0;
5559 
5560  // look up the right row in the
5561  // table and take the data from
5562  // there
5563  const unsigned int row =
5564  this->finite_element_output
5565  .shape_function_to_row_table[i * fe->n_components() + component];
5566  return this->finite_element_output.shape_values(row, j);
5567 }
5568 
5569 
5570 
5571 template <int dim, int spacedim>
5572 inline const Tensor<1, spacedim> &
5573 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
5574  const unsigned int j) const
5575 {
5576  AssertIndexRange(i, fe->n_dofs_per_cell());
5577  Assert(this->update_flags & update_gradients,
5578  ExcAccessToUninitializedField("update_gradients"));
5579  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5580  Assert(present_cell.get() != nullptr,
5581  ExcMessage("FEValues object is not reinit'ed to any cell"));
5582  // if the entire FE is primitive,
5583  // then we can take a short-cut:
5584  if (fe->is_primitive())
5585  return this->finite_element_output.shape_gradients[i][j];
5586  else
5587  {
5588  // otherwise, use the mapping
5589  // between shape function
5590  // numbers and rows. note that
5591  // by the assertions above, we
5592  // know that this particular
5593  // shape function is primitive,
5594  // so we can call
5595  // system_to_component_index
5596  const unsigned int row =
5597  this->finite_element_output
5598  .shape_function_to_row_table[i * fe->n_components() +
5599  fe->system_to_component_index(i).first];
5600  return this->finite_element_output.shape_gradients[row][j];
5601  }
5602 }
5603 
5604 
5605 
5606 template <int dim, int spacedim>
5607 inline Tensor<1, spacedim>
5609  const unsigned int i,
5610  const unsigned int j,
5611  const unsigned int component) const
5612 {
5613  AssertIndexRange(i, fe->n_dofs_per_cell());
5614  Assert(this->update_flags & update_gradients,
5615  ExcAccessToUninitializedField("update_gradients"));
5616  AssertIndexRange(component, fe->n_components());
5617  Assert(present_cell.get() != nullptr,
5618  ExcMessage("FEValues object is not reinit'ed to any cell"));
5619  // check whether the shape function
5620  // is non-zero at all within
5621  // this component:
5622  if (fe->get_nonzero_components(i)[component] == false)
5623  return Tensor<1, spacedim>();
5624 
5625  // look up the right row in the
5626  // table and take the data from
5627  // there
5628  const unsigned int row =
5629  this->finite_element_output
5630  .shape_function_to_row_table[i * fe->n_components() + component];
5631  return this->finite_element_output.shape_gradients[row][j];
5632 }
5633 
5634 
5635 
5636 template <int dim, int spacedim>
5637 inline const Tensor<2, spacedim> &
5638 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
5639  const unsigned int j) const
5640 {
5641  AssertIndexRange(i, fe->n_dofs_per_cell());
5642  Assert(this->update_flags & update_hessians,
5643  ExcAccessToUninitializedField("update_hessians"));
5644  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5645  Assert(present_cell.get() != nullptr,
5646  ExcMessage("FEValues object is not reinit'ed to any cell"));
5647  // if the entire FE is primitive,
5648  // then we can take a short-cut:
5649  if (fe->is_primitive())
5650  return this->finite_element_output.shape_hessians[i][j];
5651  else
5652  {
5653  // otherwise, use the mapping
5654  // between shape function
5655  // numbers and rows. note that
5656  // by the assertions above, we
5657  // know that this particular
5658  // shape function is primitive,
5659  // so we can call
5660  // system_to_component_index
5661  const unsigned int row =
5662  this->finite_element_output
5663  .shape_function_to_row_table[i * fe->n_components() +
5664  fe->system_to_component_index(i).first];
5665  return this->finite_element_output.shape_hessians[row][j];
5666  }
5667 }
5668 
5669 
5670 
5671 template <int dim, int spacedim>
5672 inline Tensor<2, spacedim>
5674  const unsigned int i,
5675  const unsigned int j,
5676  const unsigned int component) const
5677 {
5678  AssertIndexRange(i, fe->n_dofs_per_cell());
5679  Assert(this->update_flags & update_hessians,
5680  ExcAccessToUninitializedField("update_hessians"));
5681  AssertIndexRange(component, fe->n_components());
5682  Assert(present_cell.get() != nullptr,
5683  ExcMessage("FEValues object is not reinit'ed to any cell"));
5684  // check whether the shape function
5685  // is non-zero at all within
5686  // this component:
5687  if (fe->get_nonzero_components(i)[component] == false)
5688  return Tensor<2, spacedim>();
5689 
5690  // look up the right row in the
5691  // table and take the data from
5692  // there
5693  const unsigned int row =
5694  this->finite_element_output
5695  .shape_function_to_row_table[i * fe->n_components() + component];
5696  return this->finite_element_output.shape_hessians[row][j];
5697 }
5698 
5699 
5700 
5701 template <int dim, int spacedim>
5702 inline const Tensor<3, spacedim> &
5704  const unsigned int j) const
5705 {
5706  AssertIndexRange(i, fe->n_dofs_per_cell());
5707  Assert(this->update_flags & update_3rd_derivatives,
5708  ExcAccessToUninitializedField("update_3rd_derivatives"));
5709  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5710  Assert(present_cell.get() != nullptr,
5711  ExcMessage("FEValues object is not reinit'ed to any cell"));
5712  // if the entire FE is primitive,
5713  // then we can take a short-cut:
5714  if (fe->is_primitive())
5715  return this->finite_element_output.shape_3rd_derivatives[i][j];
5716  else
5717  {
5718  // otherwise, use the mapping
5719  // between shape function
5720  // numbers and rows. note that
5721  // by the assertions above, we
5722  // know that this particular
5723  // shape function is primitive,
5724  // so we can call
5725  // system_to_component_index
5726  const unsigned int row =
5727  this->finite_element_output
5728  .shape_function_to_row_table[i * fe->n_components() +
5729  fe->system_to_component_index(i).first];
5730  return this->finite_element_output.shape_3rd_derivatives[row][j];
5731  }
5732 }
5733 
5734 
5735 
5736 template <int dim, int spacedim>
5737 inline Tensor<3, spacedim>
5739  const unsigned int i,
5740  const unsigned int j,
5741  const unsigned int component) const
5742 {
5743  AssertIndexRange(i, fe->n_dofs_per_cell());
5744  Assert(this->update_flags & update_3rd_derivatives,
5745  ExcAccessToUninitializedField("update_3rd_derivatives"));
5746  AssertIndexRange(component, fe->n_components());
5747  Assert(present_cell.get() != nullptr,
5748  ExcMessage("FEValues object is not reinit'ed to any cell"));
5749  // check whether the shape function
5750  // is non-zero at all within
5751  // this component:
5752  if (fe->get_nonzero_components(i)[component] == false)
5753  return Tensor<3, spacedim>();
5754 
5755  // look up the right row in the
5756  // table and take the data from
5757  // there
5758  const unsigned int row =
5759  this->finite_element_output
5760  .shape_function_to_row_table[i * fe->n_components() + component];
5761  return this->finite_element_output.shape_3rd_derivatives[row][j];
5762 }
5763 
5764 
5765 
5766 template <int dim, int spacedim>
5767 inline const FiniteElement<dim, spacedim> &
5769 {
5770  return *fe;
5771 }
5772 
5773 
5774 
5775 template <int dim, int spacedim>
5776 inline const Mapping<dim, spacedim> &
5778 {
5779  return *mapping;
5780 }
5781 
5782 
5783 
5784 template <int dim, int spacedim>
5785 inline UpdateFlags
5787 {
5788  return this->update_flags;
5789 }
5790 
5791 
5792 
5793 template <int dim, int spacedim>
5794 inline const std::vector<Point<spacedim>> &
5796 {
5797  Assert(this->update_flags & update_quadrature_points,
5798  ExcAccessToUninitializedField("update_quadrature_points"));
5799  Assert(present_cell.get() != nullptr,
5800  ExcMessage("FEValues object is not reinit'ed to any cell"));
5801  return this->mapping_output.quadrature_points;
5802 }
5803 
5804 
5805 
5806 template <int dim, int spacedim>
5807 inline const std::vector<double> &
5809 {
5810  Assert(this->update_flags & update_JxW_values,
5811  ExcAccessToUninitializedField("update_JxW_values"));
5812  Assert(present_cell.get() != nullptr,
5813  ExcMessage("FEValues object is not reinit'ed to any cell"));
5814  return this->mapping_output.JxW_values;
5815 }
5816 
5817 
5818 
5819 template <int dim, int spacedim>
5820 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5822 {
5823  Assert(this->update_flags & update_jacobians,
5824  ExcAccessToUninitializedField("update_jacobians"));
5825  Assert(present_cell.get() != nullptr,
5826  ExcMessage("FEValues object is not reinit'ed to any cell"));
5827  return this->mapping_output.jacobians;
5828 }
5829 
5830 
5831 
5832 template <int dim, int spacedim>
5833 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5835 {
5836  Assert(this->update_flags & update_jacobian_grads,
5837  ExcAccessToUninitializedField("update_jacobians_grads"));
5838  Assert(present_cell.get() != nullptr,
5839  ExcMessage("FEValues object is not reinit'ed to any cell"));
5840  return this->mapping_output.jacobian_grads;
5841 }
5842 
5843 
5844 
5845 template <int dim, int spacedim>
5846 inline const Tensor<3, spacedim> &
5848  const unsigned int i) const
5849 {
5850  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5851  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5852  Assert(present_cell.get() != nullptr,
5853  ExcMessage("FEValues object is not reinit'ed to any cell"));
5854  return this->mapping_output.jacobian_pushed_forward_grads[i];
5855 }
5856 
5857 
5858 
5859 template <int dim, int spacedim>
5860 inline const std::vector<Tensor<3, spacedim>> &
5862 {
5863  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5864  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5865  Assert(present_cell.get() != nullptr,
5866  ExcMessage("FEValues object is not reinit'ed to any cell"));
5867  return this->mapping_output.jacobian_pushed_forward_grads;
5868 }
5869 
5870 
5871 
5872 template <int dim, int spacedim>
5873 inline const DerivativeForm<3, dim, spacedim> &
5874 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5875 {
5876  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5877  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5878  Assert(present_cell.get() != nullptr,
5879  ExcMessage("FEValues object is not reinit'ed to any cell"));
5880  return this->mapping_output.jacobian_2nd_derivatives[i];
5881 }
5882 
5883 
5884 
5885 template <int dim, int spacedim>
5886 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5888 {
5889  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5890  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5891  Assert(present_cell.get() != nullptr,
5892  ExcMessage("FEValues object is not reinit'ed to any cell"));
5893  return this->mapping_output.jacobian_2nd_derivatives;
5894 }
5895 
5896 
5897 
5898 template <int dim, int spacedim>
5899 inline const Tensor<4, spacedim> &
5901  const unsigned int i) const
5902 {
5905  "update_jacobian_pushed_forward_2nd_derivatives"));
5906  Assert(present_cell.get() != nullptr,
5907  ExcMessage("FEValues object is not reinit'ed to any cell"));
5908  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5909 }
5910 
5911 
5912 
5913 template <int dim, int spacedim>
5914 inline const std::vector<Tensor<4, spacedim>> &
5916 {
5917  Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
5919  "update_jacobian_pushed_forward_2nd_derivatives"));
5920  Assert(present_cell.get() != nullptr,
5921  ExcMessage("FEValues object is not reinit'ed to any cell"));
5922  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5923 }
5924 
5925 
5926 
5927 template <int dim, int spacedim>
5928 inline const DerivativeForm<4, dim, spacedim> &
5929 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5930 {
5931  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5932  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5933  Assert(present_cell.get() != nullptr,
5934  ExcMessage("FEValues object is not reinit'ed to any cell"));
5935  return this->mapping_output.jacobian_3rd_derivatives[i];
5936 }
5937 
5938 
5939 
5940 template <int dim, int spacedim>
5941 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5943 {
5944  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5945  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5946  Assert(present_cell.get() != nullptr,
5947  ExcMessage("FEValues object is not reinit'ed to any cell"));
5948  return this->mapping_output.jacobian_3rd_derivatives;
5949 }
5950 
5951 
5952 
5953 template <int dim, int spacedim>
5954 inline const Tensor<5, spacedim> &
5956  const unsigned int i) const
5957 {
5960  "update_jacobian_pushed_forward_3rd_derivatives"));
5961  Assert(present_cell.get() != nullptr,
5962  ExcMessage("FEValues object is not reinit'ed to any cell"));
5963  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
5964 }
5965 
5966 
5967 
5968 template <int dim, int spacedim>
5969 inline const std::vector<Tensor<5, spacedim>> &
5971 {
5972  Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
5974  "update_jacobian_pushed_forward_3rd_derivatives"));
5975  Assert(present_cell.get() != nullptr,
5976  ExcMessage("FEValues object is not reinit'ed to any cell"));
5977  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
5978 }
5979 
5980 
5981 
5982 template <int dim, int spacedim>
5983 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
5985 {
5986  Assert(this->update_flags & update_inverse_jacobians,
5987  ExcAccessToUninitializedField("update_inverse_jacobians"));
5988  Assert(present_cell.get() != nullptr,
5989  ExcMessage("FEValues object is not reinit'ed to any cell"));
5990  return this->mapping_output.inverse_jacobians;
5991 }
5992 
5993 
5994 
5995 template <int dim, int spacedim>
5998 {
5999  return {0U, dofs_per_cell};
6000 }
6001 
6002 
6003 
6004 template <int dim, int spacedim>
6007  const unsigned int start_dof_index) const
6008 {
6009  Assert(start_dof_index <= dofs_per_cell,
6010  ExcIndexRange(start_dof_index, 0, dofs_per_cell + 1));
6011  return {start_dof_index, dofs_per_cell};
6012 }
6013 
6014 
6015 
6016 template <int dim, int spacedim>
6019  const unsigned int end_dof_index) const
6020 {
6021  Assert(end_dof_index < dofs_per_cell,
6022  ExcIndexRange(end_dof_index, 0, dofs_per_cell));
6023  return {0U, end_dof_index + 1};
6024 }
6025 
6026 
6027 
6028 template <int dim, int spacedim>
6031 {
6032  return {0U, n_quadrature_points};
6033 }
6034 
6035 
6036 
6037 template <int dim, int spacedim>
6038 inline const Point<spacedim> &
6039 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
6040 {
6041  Assert(this->update_flags & update_quadrature_points,
6042  ExcAccessToUninitializedField("update_quadrature_points"));
6043  AssertIndexRange(i, this->mapping_output.quadrature_points.size());
6044  Assert(present_cell.get() != nullptr,
6045  ExcMessage("FEValues object is not reinit'ed to any cell"));
6046 
6047  return this->mapping_output.quadrature_points[i];
6048 }
6049 
6050 
6051 
6052 template <int dim, int spacedim>
6053 inline double
6054 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
6055 {
6056  Assert(this->update_flags & update_JxW_values,
6057  ExcAccessToUninitializedField("update_JxW_values"));
6058  AssertIndexRange(i, this->mapping_output.JxW_values.size());
6059  Assert(present_cell.get() != nullptr,
6060  ExcMessage("FEValues object is not reinit'ed to any cell"));
6061 
6062  return this->mapping_output.JxW_values[i];
6063 }
6064 
6065 
6066 
6067 template <int dim, int spacedim>
6068 inline const DerivativeForm<1, dim, spacedim> &
6069 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
6070 {
6071  Assert(this->update_flags & update_jacobians,
6072  ExcAccessToUninitializedField("update_jacobians"));
6073  AssertIndexRange(i, this->mapping_output.jacobians.size());
6074  Assert(present_cell.get() != nullptr,
6075  ExcMessage("FEValues object is not reinit'ed to any cell"));
6076 
6077  return this->mapping_output.jacobians[i];
6078 }
6079 
6080 
6081 
6082 template <int dim, int spacedim>
6083 inline const DerivativeForm<2, dim, spacedim> &
6084 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
6085 {
6086  Assert(this->update_flags & update_jacobian_grads,
6087  ExcAccessToUninitializedField("update_jacobians_grads"));
6088  AssertIndexRange(i, this->mapping_output.jacobian_grads.size());
6089  Assert(present_cell.get() != nullptr,
6090  ExcMessage("FEValues object is not reinit'ed to any cell"));
6091 
6092  return this->mapping_output.jacobian_grads[i];
6093 }
6094 
6095 
6096 
6097 template <int dim, int spacedim>
6098 inline const DerivativeForm<1, spacedim, dim> &
6099 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
6100 {
6101  Assert(this->update_flags & update_inverse_jacobians,
6102  ExcAccessToUninitializedField("update_inverse_jacobians"));
6103  AssertIndexRange(i, this->mapping_output.inverse_jacobians.size());
6104  Assert(present_cell.get() != nullptr,
6105  ExcMessage("FEValues object is not reinit'ed to any cell"));
6106 
6107  return this->mapping_output.inverse_jacobians[i];
6108 }
6109 
6110 
6111 
6112 template <int dim, int spacedim>
6113 inline const Tensor<1, spacedim> &
6114 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
6115 {
6116  Assert(this->update_flags & update_normal_vectors,
6118  "update_normal_vectors")));
6119  AssertIndexRange(i, this->mapping_output.normal_vectors.size());
6120  Assert(present_cell.get() != nullptr,
6121  ExcMessage("FEValues object is not reinit'ed to any cell"));
6122 
6123  return this->mapping_output.normal_vectors[i];
6124 }
6125 
6126 
6127 
6128 /*--------------------- Inline functions: FEValues --------------------------*/
6129 
6130 
6131 template <int dim, int spacedim>
6132 inline const Quadrature<dim> &
6134 {
6135  return quadrature;
6136 }
6137 
6138 
6139 
6140 template <int dim, int spacedim>
6141 inline const FEValues<dim, spacedim> &
6143 {
6144  return *this;
6145 }
6146 
6147 
6148 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
6149 
6150 
6151 template <int dim, int spacedim>
6152 inline unsigned int
6154 {
6155  return present_face_no;
6156 }
6157 
6158 
6159 template <int dim, int spacedim>
6160 inline unsigned int
6162 {
6163  return present_face_index;
6164 }
6165 
6166 
6167 /*----------------------- Inline functions: FE*FaceValues -------------------*/
6168 
6169 template <int dim, int spacedim>
6170 inline const Quadrature<dim - 1> &
6172 {
6173  return quadrature[quadrature.size() == 1 ? 0 : present_face_no];
6174 }
6175 
6176 
6177 
6178 template <int dim, int spacedim>
6179 inline const FEFaceValues<dim, spacedim> &
6181 {
6182  return *this;
6183 }
6184 
6185 
6186 
6187 template <int dim, int spacedim>
6188 inline const FESubfaceValues<dim, spacedim> &
6190 {
6191  return *this;
6192 }
6193 
6194 
6195 
6196 template <int dim, int spacedim>
6197 inline const Tensor<1, spacedim> &
6198 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
6199 {
6200  AssertIndexRange(i, this->mapping_output.boundary_forms.size());
6201  Assert(this->update_flags & update_boundary_forms,
6203  "update_boundary_forms")));
6204 
6205  return this->mapping_output.boundary_forms[i];
6206 }
6207 
6208 #endif // DOXYGEN
6209 
6211 
6212 #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:4222
unsigned int present_face_index
Definition: fe_values.h:4228
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
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:1720
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
::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:4233
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:509
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:1461
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:487
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:464
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
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:540
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
unsigned int get_face_number() const
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
Definition: tensor.h:506
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:400
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
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
#define DEAL_II_DEPRECATED
Definition: config.h:160
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