Reference documentation for deal.II version Git 71165045b9 2021-01-23 20:21:06 -0500
\(\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 - 2020 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.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 
181  template <typename Number>
182  struct OutputType
183  {
188  using value_type =
189  typename ProductType<Number,
191 
196  using gradient_type = typename ProductType<
197  Number,
199 
204  using laplacian_type =
205  typename ProductType<Number,
207 
212  using hessian_type = typename ProductType<
213  Number,
215 
220  using third_derivative_type = typename ProductType<
221  Number,
223  };
224 
230  {
240 
249  unsigned int row_index;
250  };
251 
255  Scalar();
256 
262  Scalar(const FEValuesBase<dim, spacedim> &fe_values_base,
263  const unsigned int component);
264 
269  Scalar &
270  operator=(const Scalar<dim, spacedim> &) = delete;
271 
285  value_type
286  value(const unsigned int shape_function, const unsigned int q_point) const;
287 
299  gradient(const unsigned int shape_function,
300  const unsigned int q_point) const;
301 
313  hessian(const unsigned int shape_function,
314  const unsigned int q_point) const;
315 
327  third_derivative(const unsigned int shape_function,
328  const unsigned int q_point) const;
329 
347  template <class InputVector>
348  void
349  get_function_values(
350  const InputVector &fe_function,
351  std::vector<typename ProductType<value_type,
352  typename InputVector::value_type>::type>
353  &values) const;
354 
389  template <class InputVector>
390  void
391  get_function_values_from_local_dof_values(
392  const InputVector &dof_values,
393  std::vector<
395  &values) const;
396 
414  template <class InputVector>
415  void
416  get_function_gradients(
417  const InputVector &fe_function,
418  std::vector<typename ProductType<gradient_type,
419  typename InputVector::value_type>::type>
420  &gradients) const;
421 
428  template <class InputVector>
429  void
430  get_function_gradients_from_local_dof_values(
431  const InputVector &dof_values,
432  std::vector<
434  &gradients) const;
435 
453  template <class InputVector>
454  void
455  get_function_hessians(
456  const InputVector &fe_function,
457  std::vector<typename ProductType<hessian_type,
458  typename InputVector::value_type>::type>
459  &hessians) const;
460 
467  template <class InputVector>
468  void
469  get_function_hessians_from_local_dof_values(
470  const InputVector &dof_values,
471  std::vector<
473  &hessians) const;
474 
475 
494  template <class InputVector>
495  void
496  get_function_laplacians(
497  const InputVector &fe_function,
498  std::vector<typename ProductType<value_type,
499  typename InputVector::value_type>::type>
500  &laplacians) const;
501 
508  template <class InputVector>
509  void
510  get_function_laplacians_from_local_dof_values(
511  const InputVector &dof_values,
512  std::vector<
514  &laplacians) const;
515 
516 
535  template <class InputVector>
536  void
537  get_function_third_derivatives(
538  const InputVector &fe_function,
539  std::vector<typename ProductType<third_derivative_type,
540  typename InputVector::value_type>::type>
541  &third_derivatives) const;
542 
549  template <class InputVector>
550  void
551  get_function_third_derivatives_from_local_dof_values(
552  const InputVector & dof_values,
553  std::vector<typename OutputType<typename InputVector::value_type>::
554  third_derivative_type> &third_derivatives) const;
555 
556 
557  private:
562 
567  const unsigned int component;
568 
572  std::vector<ShapeFunctionData> shape_function_data;
573  };
574 
575 
576 
606  template <int dim, int spacedim = dim>
607  class Vector
608  {
609  public:
616 
626 
638 
645 
652  using curl_type = typename ::internal::CurlType<spacedim>::type;
653 
660 
667 
672  template <typename Number>
673  struct OutputType
674  {
679  using value_type =
680  typename ProductType<Number,
682 
687  using gradient_type = typename ProductType<
688  Number,
690 
695  using symmetric_gradient_type = typename ProductType<
696  Number,
698 
703  using divergence_type = typename ProductType<
704  Number,
706 
711  using laplacian_type =
712  typename ProductType<Number,
714 
719  using curl_type =
720  typename ProductType<Number,
722 
727  using hessian_type = typename ProductType<
728  Number,
730 
735  using third_derivative_type = typename ProductType<
736  Number,
738  };
739 
745  {
754  bool is_nonzero_shape_function_component[spacedim];
755 
765  unsigned int row_index[spacedim];
766 
777  };
778 
782  Vector();
783 
792  Vector(const FEValuesBase<dim, spacedim> &fe_values_base,
793  const unsigned int first_vector_component);
794 
799  Vector &
800  operator=(const Vector<dim, spacedim> &) = delete;
801 
818  value_type
819  value(const unsigned int shape_function, const unsigned int q_point) const;
820 
835  gradient(const unsigned int shape_function,
836  const unsigned int q_point) const;
837 
854  symmetric_gradient(const unsigned int shape_function,
855  const unsigned int q_point) const;
856 
868  divergence(const unsigned int shape_function,
869  const unsigned int q_point) const;
870 
891  curl_type
892  curl(const unsigned int shape_function, const unsigned int q_point) const;
893 
905  hessian(const unsigned int shape_function,
906  const unsigned int q_point) const;
907 
919  third_derivative(const unsigned int shape_function,
920  const unsigned int q_point) const;
921 
939  template <class InputVector>
940  void
941  get_function_values(
942  const InputVector &fe_function,
943  std::vector<typename ProductType<value_type,
944  typename InputVector::value_type>::type>
945  &values) const;
946 
981  template <class InputVector>
982  void
983  get_function_values_from_local_dof_values(
984  const InputVector &dof_values,
985  std::vector<
987  &values) const;
988 
1006  template <class InputVector>
1007  void
1008  get_function_gradients(
1009  const InputVector &fe_function,
1010  std::vector<typename ProductType<gradient_type,
1011  typename InputVector::value_type>::type>
1012  &gradients) const;
1013 
1020  template <class InputVector>
1021  void
1022  get_function_gradients_from_local_dof_values(
1023  const InputVector &dof_values,
1024  std::vector<
1026  &gradients) const;
1027 
1051  template <class InputVector>
1052  void
1053  get_function_symmetric_gradients(
1054  const InputVector &fe_function,
1055  std::vector<typename ProductType<symmetric_gradient_type,
1056  typename InputVector::value_type>::type>
1057  &symmetric_gradients) const;
1058 
1065  template <class InputVector>
1066  void
1067  get_function_symmetric_gradients_from_local_dof_values(
1068  const InputVector & dof_values,
1069  std::vector<typename OutputType<typename InputVector::value_type>::
1070  symmetric_gradient_type> &symmetric_gradients) const;
1071 
1090  template <class InputVector>
1091  void
1092  get_function_divergences(
1093  const InputVector &fe_function,
1094  std::vector<typename ProductType<divergence_type,
1095  typename InputVector::value_type>::type>
1096  &divergences) const;
1097 
1104  template <class InputVector>
1105  void
1106  get_function_divergences_from_local_dof_values(
1107  const InputVector &dof_values,
1108  std::vector<
1110  &divergences) const;
1111 
1130  template <class InputVector>
1131  void
1132  get_function_curls(
1133  const InputVector &fe_function,
1134  std::vector<
1136  &curls) const;
1137 
1144  template <class InputVector>
1145  void
1146  get_function_curls_from_local_dof_values(
1147  const InputVector &dof_values,
1148  std::vector<
1150  &curls) const;
1151 
1169  template <class InputVector>
1170  void
1171  get_function_hessians(
1172  const InputVector &fe_function,
1173  std::vector<typename ProductType<hessian_type,
1174  typename InputVector::value_type>::type>
1175  &hessians) const;
1176 
1183  template <class InputVector>
1184  void
1185  get_function_hessians_from_local_dof_values(
1186  const InputVector &dof_values,
1187  std::vector<
1189  &hessians) const;
1190 
1209  template <class InputVector>
1210  void
1211  get_function_laplacians(
1212  const InputVector &fe_function,
1213  std::vector<typename ProductType<value_type,
1214  typename InputVector::value_type>::type>
1215  &laplacians) const;
1216 
1223  template <class InputVector>
1224  void
1225  get_function_laplacians_from_local_dof_values(
1226  const InputVector &dof_values,
1227  std::vector<
1229  &laplacians) const;
1230 
1249  template <class InputVector>
1250  void
1251  get_function_third_derivatives(
1252  const InputVector &fe_function,
1253  std::vector<typename ProductType<third_derivative_type,
1254  typename InputVector::value_type>::type>
1255  &third_derivatives) const;
1256 
1263  template <class InputVector>
1264  void
1265  get_function_third_derivatives_from_local_dof_values(
1266  const InputVector & dof_values,
1267  std::vector<typename OutputType<typename InputVector::value_type>::
1268  third_derivative_type> &third_derivatives) const;
1269 
1270  private:
1275 
1280  const unsigned int first_vector_component;
1281 
1285  std::vector<ShapeFunctionData> shape_function_data;
1286  };
1287 
1288 
1289  template <int rank, int dim, int spacedim = dim>
1291 
1314  template <int dim, int spacedim>
1315  class SymmetricTensor<2, dim, spacedim>
1316  {
1317  public:
1325 
1336 
1341  template <typename Number>
1342  struct OutputType
1343  {
1348  using value_type = typename ProductType<
1349  Number,
1351 
1356  using divergence_type = typename ProductType<
1357  Number,
1359  };
1360 
1365  struct ShapeFunctionData
1366  {
1375  bool is_nonzero_shape_function_component
1376  [value_type::n_independent_components];
1377 
1387  unsigned int row_index[value_type::n_independent_components];
1388 
1398 
1403  };
1404 
1408  SymmetricTensor();
1409 
1419  SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1420  const unsigned int first_tensor_component);
1421 
1426  SymmetricTensor &
1427  operator=(const SymmetricTensor<2, dim, spacedim> &) = delete;
1428 
1446  value_type
1447  value(const unsigned int shape_function, const unsigned int q_point) const;
1448 
1463  divergence(const unsigned int shape_function,
1464  const unsigned int q_point) const;
1465 
1483  template <class InputVector>
1484  void
1485  get_function_values(
1486  const InputVector &fe_function,
1487  std::vector<typename ProductType<value_type,
1488  typename InputVector::value_type>::type>
1489  &values) const;
1490 
1525  template <class InputVector>
1526  void
1527  get_function_values_from_local_dof_values(
1528  const InputVector &dof_values,
1529  std::vector<
1531  &values) const;
1532 
1554  template <class InputVector>
1555  void
1556  get_function_divergences(
1557  const InputVector &fe_function,
1558  std::vector<typename ProductType<divergence_type,
1559  typename InputVector::value_type>::type>
1560  &divergences) const;
1561 
1568  template <class InputVector>
1569  void
1570  get_function_divergences_from_local_dof_values(
1571  const InputVector &dof_values,
1572  std::vector<
1574  &divergences) const;
1575 
1576  private:
1581 
1586  const unsigned int first_tensor_component;
1587 
1591  std::vector<ShapeFunctionData> shape_function_data;
1592  };
1593 
1594 
1595  template <int rank, int dim, int spacedim = dim>
1596  class Tensor;
1597 
1616  template <int dim, int spacedim>
1617  class Tensor<2, dim, spacedim>
1618  {
1619  public:
1625 
1630 
1636 
1641  template <typename Number>
1642  struct OutputType
1643  {
1648  using value_type = typename ProductType<
1649  Number,
1651 
1656  using divergence_type = typename ProductType<
1657  Number,
1659 
1664  using gradient_type = typename ProductType<
1665  Number,
1667  };
1668 
1673  struct ShapeFunctionData
1674  {
1683  bool is_nonzero_shape_function_component
1684  [value_type::n_independent_components];
1685 
1695  unsigned int row_index[value_type::n_independent_components];
1696 
1706 
1711  };
1712 
1716  Tensor();
1717 
1727  Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1728  const unsigned int first_tensor_component);
1729 
1730 
1735  Tensor &
1736  operator=(const Tensor<2, dim, spacedim> &) = delete;
1737 
1754  value_type
1755  value(const unsigned int shape_function, const unsigned int q_point) const;
1756 
1771  divergence(const unsigned int shape_function,
1772  const unsigned int q_point) const;
1773 
1788  gradient(const unsigned int shape_function,
1789  const unsigned int q_point) const;
1790 
1808  template <class InputVector>
1809  void
1810  get_function_values(
1811  const InputVector &fe_function,
1812  std::vector<typename ProductType<value_type,
1813  typename InputVector::value_type>::type>
1814  &values) const;
1815 
1850  template <class InputVector>
1851  void
1852  get_function_values_from_local_dof_values(
1853  const InputVector &dof_values,
1854  std::vector<
1856  &values) const;
1857 
1879  template <class InputVector>
1880  void
1881  get_function_divergences(
1882  const InputVector &fe_function,
1883  std::vector<typename ProductType<divergence_type,
1884  typename InputVector::value_type>::type>
1885  &divergences) const;
1886 
1893  template <class InputVector>
1894  void
1895  get_function_divergences_from_local_dof_values(
1896  const InputVector &dof_values,
1897  std::vector<
1899  &divergences) const;
1900 
1917  template <class InputVector>
1918  void
1919  get_function_gradients(
1920  const InputVector &fe_function,
1921  std::vector<typename ProductType<gradient_type,
1922  typename InputVector::value_type>::type>
1923  &gradients) const;
1924 
1931  template <class InputVector>
1932  void
1933  get_function_gradients_from_local_dof_values(
1934  const InputVector &dof_values,
1935  std::vector<
1937  &gradients) const;
1938 
1939  private:
1944 
1949  const unsigned int first_tensor_component;
1950 
1954  std::vector<ShapeFunctionData> shape_function_data;
1955  };
1956 
1957 } // namespace FEValuesViews
1958 
1959 
1960 namespace internal
1961 {
1962  namespace FEValuesViews
1963  {
1968  template <int dim, int spacedim, typename Extractor>
1969  struct ViewType
1970  {};
1971 
1979  template <int dim, int spacedim>
1980  struct ViewType<dim, spacedim, FEValuesExtractors::Scalar>
1981  {
1982  using type = typename ::FEValuesViews::Scalar<dim, spacedim>;
1983  };
1984 
1992  template <int dim, int spacedim>
1993  struct ViewType<dim, spacedim, FEValuesExtractors::Vector>
1994  {
1995  using type = typename ::FEValuesViews::Vector<dim, spacedim>;
1996  };
1997 
2005  template <int dim, int spacedim, int rank>
2006  struct ViewType<dim, spacedim, FEValuesExtractors::Tensor<rank>>
2007  {
2008  using type = typename ::FEValuesViews::Tensor<rank, dim, spacedim>;
2009  };
2010 
2018  template <int dim, int spacedim, int rank>
2019  struct ViewType<dim, spacedim, FEValuesExtractors::SymmetricTensor<rank>>
2020  {
2021  using type =
2022  typename ::FEValuesViews::SymmetricTensor<rank, dim, spacedim>;
2023  };
2024 
2032  template <int dim, int spacedim>
2033  struct Cache
2034  {
2039  std::vector<::FEValuesViews::Scalar<dim, spacedim>> scalars;
2040  std::vector<::FEValuesViews::Vector<dim, spacedim>> vectors;
2041  std::vector<::FEValuesViews::SymmetricTensor<2, dim, spacedim>>
2043  std::vector<::FEValuesViews::Tensor<2, dim, spacedim>>
2045 
2049  Cache(const FEValuesBase<dim, spacedim> &fe_values);
2050  };
2051  } // namespace FEValuesViews
2052 } // namespace internal
2053 
2054 namespace FEValuesViews
2055 {
2060  template <int dim, int spacedim, typename Extractor>
2061  using View = typename ::internal::FEValuesViews::
2062  ViewType<dim, spacedim, Extractor>::type;
2063 } // namespace FEValuesViews
2064 
2065 
2165 template <int dim, int spacedim>
2166 class FEValuesBase : public Subscriptor
2167 {
2168 public:
2172  static const unsigned int dimension = dim;
2173 
2177  static const unsigned int space_dimension = spacedim;
2178 
2186  const unsigned int n_quadrature_points;
2187 
2197  const unsigned int max_n_quadrature_points;
2198 
2204  const unsigned int dofs_per_cell;
2205 
2206 
2214  FEValuesBase(const unsigned int n_q_points,
2215  const unsigned int dofs_per_cell,
2216  const UpdateFlags update_flags,
2217  const Mapping<dim, spacedim> & mapping,
2218  const FiniteElement<dim, spacedim> &fe);
2219 
2224  FEValuesBase &
2225  operator=(const FEValuesBase &) = delete;
2226 
2231  FEValuesBase(const FEValuesBase &) = delete;
2232 
2236  virtual ~FEValuesBase() override;
2237 
2238 
2242 
2243 
2264  const double &
2265  shape_value(const unsigned int function_no,
2266  const unsigned int point_no) const;
2267 
2288  double
2289  shape_value_component(const unsigned int function_no,
2290  const unsigned int point_no,
2291  const unsigned int component) const;
2292 
2318  const Tensor<1, spacedim> &
2319  shape_grad(const unsigned int function_no,
2320  const unsigned int quadrature_point) const;
2321 
2339  shape_grad_component(const unsigned int function_no,
2340  const unsigned int point_no,
2341  const unsigned int component) const;
2342 
2362  const Tensor<2, spacedim> &
2363  shape_hessian(const unsigned int function_no,
2364  const unsigned int point_no) const;
2365 
2383  shape_hessian_component(const unsigned int function_no,
2384  const unsigned int point_no,
2385  const unsigned int component) const;
2386 
2406  const Tensor<3, spacedim> &
2407  shape_3rd_derivative(const unsigned int function_no,
2408  const unsigned int point_no) const;
2409 
2427  shape_3rd_derivative_component(const unsigned int function_no,
2428  const unsigned int point_no,
2429  const unsigned int component) const;
2430 
2432 
2434 
2471  template <class InputVector>
2472  void
2473  get_function_values(
2474  const InputVector & fe_function,
2475  std::vector<typename InputVector::value_type> &values) const;
2476 
2490  template <class InputVector>
2491  void
2492  get_function_values(
2493  const InputVector & fe_function,
2494  std::vector<Vector<typename InputVector::value_type>> &values) const;
2495 
2552  template <class InputVector>
2553  void
2554  get_function_values(
2555  const InputVector & fe_function,
2557  std::vector<typename InputVector::value_type> & values) const;
2558 
2567  template <class InputVector>
2568  void
2569  get_function_values(
2570  const InputVector & fe_function,
2572  std::vector<Vector<typename InputVector::value_type>> &values) const;
2573 
2574 
2596  template <class InputVector>
2597  void
2598  get_function_values(
2599  const InputVector & fe_function,
2601  ArrayView<std::vector<typename InputVector::value_type>> values,
2602  const bool quadrature_points_fastest) const;
2603 
2605 
2607 
2644  template <class InputVector>
2645  void
2646  get_function_gradients(
2647  const InputVector &fe_function,
2649  &gradients) const;
2650 
2667  template <class InputVector>
2668  void
2669  get_function_gradients(
2670  const InputVector &fe_function,
2671  std::vector<
2673  &gradients) const;
2674 
2683  template <class InputVector>
2684  void
2685  get_function_gradients(
2686  const InputVector & fe_function,
2689  &gradients) const;
2690 
2699  template <class InputVector>
2700  void
2701  get_function_gradients(
2702  const InputVector & fe_function,
2704  ArrayView<
2706  gradients,
2707  const bool quadrature_points_fastest = false) const;
2708 
2710 
2714 
2752  template <class InputVector>
2753  void
2754  get_function_hessians(
2755  const InputVector &fe_function,
2757  &hessians) const;
2758 
2776  template <class InputVector>
2777  void
2778  get_function_hessians(
2779  const InputVector &fe_function,
2780  std::vector<
2782  & hessians,
2783  const bool quadrature_points_fastest = false) const;
2784 
2793  template <class InputVector>
2794  void
2795  get_function_hessians(
2796  const InputVector & fe_function,
2799  &hessians) const;
2800 
2809  template <class InputVector>
2810  void
2811  get_function_hessians(
2812  const InputVector & fe_function,
2814  ArrayView<
2816  hessians,
2817  const bool quadrature_points_fastest = false) const;
2818 
2859  template <class InputVector>
2860  void
2861  get_function_laplacians(
2862  const InputVector & fe_function,
2863  std::vector<typename InputVector::value_type> &laplacians) const;
2864 
2884  template <class InputVector>
2885  void
2886  get_function_laplacians(
2887  const InputVector & fe_function,
2888  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2889 
2898  template <class InputVector>
2899  void
2900  get_function_laplacians(
2901  const InputVector & fe_function,
2903  std::vector<typename InputVector::value_type> & laplacians) const;
2904 
2913  template <class InputVector>
2914  void
2915  get_function_laplacians(
2916  const InputVector & fe_function,
2918  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2919 
2928  template <class InputVector>
2929  void
2930  get_function_laplacians(
2931  const InputVector & fe_function,
2933  std::vector<std::vector<typename InputVector::value_type>> &laplacians,
2934  const bool quadrature_points_fastest = false) const;
2935 
2937 
2939 
2978  template <class InputVector>
2979  void
2980  get_function_third_derivatives(
2981  const InputVector &fe_function,
2983  &third_derivatives) const;
2984 
3003  template <class InputVector>
3004  void
3005  get_function_third_derivatives(
3006  const InputVector &fe_function,
3007  std::vector<
3009  & third_derivatives,
3010  const bool quadrature_points_fastest = false) const;
3011 
3020  template <class InputVector>
3021  void
3022  get_function_third_derivatives(
3023  const InputVector & fe_function,
3026  &third_derivatives) const;
3027 
3036  template <class InputVector>
3037  void
3038  get_function_third_derivatives(
3039  const InputVector & fe_function,
3041  ArrayView<
3043  third_derivatives,
3044  const bool quadrature_points_fastest = false) const;
3046 
3048 
3049 
3074  dof_indices() const;
3075 
3109  dof_indices_starting_at(const unsigned int start_dof_index) const;
3110 
3142  dof_indices_ending_at(const unsigned int end_dof_index) const;
3143 
3145 
3147 
3148 
3170  quadrature_point_indices() const;
3171 
3177  const Point<spacedim> &
3178  quadrature_point(const unsigned int q) const;
3179 
3185  const std::vector<Point<spacedim>> &
3186  get_quadrature_points() const;
3187 
3203  double
3204  JxW(const unsigned int quadrature_point) const;
3205 
3209  const std::vector<double> &
3210  get_JxW_values() const;
3211 
3219  jacobian(const unsigned int quadrature_point) const;
3220 
3227  const std::vector<DerivativeForm<1, dim, spacedim>> &
3228  get_jacobians() const;
3229 
3238  jacobian_grad(const unsigned int quadrature_point) const;
3239 
3246  const std::vector<DerivativeForm<2, dim, spacedim>> &
3247  get_jacobian_grads() const;
3248 
3257  const Tensor<3, spacedim> &
3258  jacobian_pushed_forward_grad(const unsigned int quadrature_point) const;
3259 
3266  const std::vector<Tensor<3, spacedim>> &
3267  get_jacobian_pushed_forward_grads() const;
3268 
3277  jacobian_2nd_derivative(const unsigned int quadrature_point) const;
3278 
3285  const std::vector<DerivativeForm<3, dim, spacedim>> &
3286  get_jacobian_2nd_derivatives() const;
3287 
3297  const Tensor<4, spacedim> &
3298  jacobian_pushed_forward_2nd_derivative(
3299  const unsigned int quadrature_point) const;
3300 
3307  const std::vector<Tensor<4, spacedim>> &
3308  get_jacobian_pushed_forward_2nd_derivatives() const;
3309 
3319  jacobian_3rd_derivative(const unsigned int quadrature_point) const;
3320 
3327  const std::vector<DerivativeForm<4, dim, spacedim>> &
3328  get_jacobian_3rd_derivatives() const;
3329 
3339  const Tensor<5, spacedim> &
3340  jacobian_pushed_forward_3rd_derivative(
3341  const unsigned int quadrature_point) const;
3342 
3349  const std::vector<Tensor<5, spacedim>> &
3350  get_jacobian_pushed_forward_3rd_derivatives() const;
3351 
3359  inverse_jacobian(const unsigned int quadrature_point) const;
3360 
3367  const std::vector<DerivativeForm<1, spacedim, dim>> &
3368  get_inverse_jacobians() const;
3369 
3389  const Tensor<1, spacedim> &
3390  normal_vector(const unsigned int i) const;
3391 
3399  const std::vector<Tensor<1, spacedim>> &
3400  get_normal_vectors() const;
3401 
3403 
3405 
3406 
3416  operator[](const FEValuesExtractors::Scalar &scalar) const;
3417 
3427  operator[](const FEValuesExtractors::Vector &vector) const;
3428 
3439  operator[](const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
3440 
3441 
3451  operator[](const FEValuesExtractors::Tensor<2> &tensor) const;
3452 
3454 
3456 
3457 
3461  const Mapping<dim, spacedim> &
3462  get_mapping() const;
3463 
3468  get_fe() const;
3469 
3473  UpdateFlags
3474  get_update_flags() const;
3475 
3480  get_cell() const;
3481 
3488  get_cell_similarity() const;
3489 
3494  std::size_t
3495  memory_consumption() const;
3497 
3498 
3507  std::string,
3508  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3509  << "object for which this kind of information has not been computed. What "
3510  << "information these objects compute is determined by the update_* flags you "
3511  << "pass to the constructor. Here, the operation you are attempting requires "
3512  << "the <" << arg1
3513  << "> flag to be set, but it was apparently not specified "
3514  << "upon construction.");
3515 
3523  ExcFEDontMatch,
3524  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3525  "to the DoFHandler that provided the cell iterator do not match.");
3531  DeclException1(ExcShapeFunctionNotPrimitive,
3532  int,
3533  << "The shape function with index " << arg1
3534  << " is not primitive, i.e. it is vector-valued and "
3535  << "has more than one non-zero vector component. This "
3536  << "function cannot be called for these shape functions. "
3537  << "Maybe you want to use the same function with the "
3538  << "_component suffix?");
3539 
3548  "The given FiniteElement is not a primitive element but the requested operation "
3549  "only works for those. See FiniteElement::is_primitive() for more information.");
3550 
3551 protected:
3582  class CellIteratorBase;
3583 
3588  template <typename CI>
3591 
3597  std::unique_ptr<const CellIteratorBase> present_cell;
3598 
3606  boost::signals2::connection tria_listener_refinement;
3607 
3615  boost::signals2::connection tria_listener_mesh_transform;
3616 
3622  void
3623  invalidate_present_cell();
3624 
3634  void
3635  maybe_invalidate_previous_present_cell(
3636  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3637 
3643 
3649  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3651 
3658 
3659 
3667 
3673  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3675 
3681  spacedim>
3683 
3684 
3689 
3698  UpdateFlags
3699  compute_update_flags(const UpdateFlags update_flags) const;
3700 
3707 
3713  void
3714  check_cell_similarity(
3715  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3716 
3717 private:
3722 
3723  // Make the view classes friends of this class, since they access internal
3724  // data.
3725  template <int, int>
3727  template <int, int>
3729  template <int, int, int>
3731  template <int, int, int>
3733 };
3734 
3735 
3736 
3746 template <int dim, int spacedim = dim>
3747 class FEValues : public FEValuesBase<dim, spacedim>
3748 {
3749 public:
3754  static const unsigned int integral_dimension = dim;
3755 
3760  FEValues(const Mapping<dim, spacedim> & mapping,
3761  const FiniteElement<dim, spacedim> &fe,
3762  const Quadrature<dim> & quadrature,
3763  const UpdateFlags update_flags);
3764 
3771  FEValues(const Mapping<dim, spacedim> & mapping,
3772  const FiniteElement<dim, spacedim> &fe,
3773  const hp::QCollection<dim> & quadrature,
3774  const UpdateFlags update_flags);
3775 
3782  const Quadrature<dim> & quadrature,
3783  const UpdateFlags update_flags);
3784 
3792  const hp::QCollection<dim> & quadrature,
3793  const UpdateFlags update_flags);
3794 
3801  template <bool level_dof_access>
3802  void
3803  reinit(
3805 
3819  void
3821 
3826  const Quadrature<dim> &
3827  get_quadrature() const;
3828 
3833  std::size_t
3834  memory_consumption() const;
3835 
3850  const FEValues<dim, spacedim> &
3851  get_present_fe_values() const;
3852 
3853 private:
3858 
3862  void
3863  initialize(const UpdateFlags update_flags);
3864 
3871  void
3872  do_reinit();
3873 };
3874 
3875 
3885 template <int dim, int spacedim = dim>
3886 class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
3887 {
3888 public:
3893  static const unsigned int integral_dimension = dim - 1;
3894 
3906  FEFaceValuesBase(const unsigned int dofs_per_cell,
3907  const UpdateFlags update_flags,
3908  const Mapping<dim, spacedim> & mapping,
3909  const FiniteElement<dim, spacedim> &fe,
3910  const Quadrature<dim - 1> & quadrature);
3911 
3918  FEFaceValuesBase(const unsigned int dofs_per_cell,
3919  const UpdateFlags update_flags,
3920  const Mapping<dim, spacedim> & mapping,
3921  const FiniteElement<dim, spacedim> &fe,
3922  const hp::QCollection<dim - 1> & quadrature);
3923 
3931  const Tensor<1, spacedim> &
3932  boundary_form(const unsigned int i) const;
3933 
3940  const std::vector<Tensor<1, spacedim>> &
3941  get_boundary_forms() const;
3942 
3947  unsigned int
3948  get_face_index() const;
3949 
3954  const Quadrature<dim - 1> &
3955  get_quadrature() const;
3956 
3961  std::size_t
3962  memory_consumption() const;
3963 
3964 protected:
3969  unsigned int present_face_no;
3970 
3975  unsigned int present_face_index;
3976 
3980  const hp::QCollection<dim - 1> quadrature;
3981 };
3982 
3983 
3984 
3998 template <int dim, int spacedim = dim>
3999 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
4000 {
4001 public:
4006  static const unsigned int dimension = dim;
4007 
4008  static const unsigned int space_dimension = spacedim;
4009 
4014  static const unsigned int integral_dimension = dim - 1;
4015 
4020  FEFaceValues(const Mapping<dim, spacedim> & mapping,
4021  const FiniteElement<dim, spacedim> &fe,
4022  const Quadrature<dim - 1> & quadrature,
4023  const UpdateFlags update_flags);
4024 
4031  FEFaceValues(const Mapping<dim, spacedim> & mapping,
4032  const FiniteElement<dim, spacedim> &fe,
4033  const hp::QCollection<dim - 1> & quadrature,
4034  const UpdateFlags update_flags);
4035 
4042  const Quadrature<dim - 1> & quadrature,
4043  const UpdateFlags update_flags);
4044 
4052  const hp::QCollection<dim - 1> & quadrature,
4053  const UpdateFlags update_flags);
4054 
4059  template <bool level_dof_access>
4060  void
4061  reinit(
4063  const unsigned int face_no);
4064 
4071  template <bool level_dof_access>
4072  void
4073  reinit(
4075  const typename Triangulation<dim, spacedim>::face_iterator & face);
4076 
4090  void
4092  const unsigned int face_no);
4093 
4094  /*
4095  * Reinitialize the gradients, Jacobi determinants, etc for the given face
4096  * on a given cell of type "iterator into a Triangulation object", and the
4097  * given finite element. Since iterators into a triangulation alone only
4098  * convey information about the geometry of a cell, but not about degrees of
4099  * freedom possibly associated with this cell, you will not be able to call
4100  * some functions of this class if they need information about degrees of
4101  * freedom. These functions are, above all, the
4102  * <tt>get_function_value/gradients/hessians/third_derivatives</tt>
4103  * functions. If you want to call these functions, you have to call the @p
4104  * reinit variants that take iterators into DoFHandler or other DoF handler
4105  * type objects.
4106  *
4107  * @note @p face must be one of @p cell's face iterators.
4108  */
4109  void
4111  const typename Triangulation<dim, spacedim>::face_iterator &face);
4112 
4128  get_present_fe_values() const;
4129 
4130 private:
4134  void
4135  initialize(const UpdateFlags update_flags);
4136 
4143  void
4144  do_reinit(const unsigned int face_no);
4145 };
4146 
4147 
4164 template <int dim, int spacedim = dim>
4165 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
4166 {
4167 public:
4171  static const unsigned int dimension = dim;
4172 
4176  static const unsigned int space_dimension = spacedim;
4177 
4182  static const unsigned int integral_dimension = dim - 1;
4183 
4188  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
4189  const FiniteElement<dim, spacedim> &fe,
4190  const Quadrature<dim - 1> & face_quadrature,
4191  const UpdateFlags update_flags);
4192 
4199  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
4200  const FiniteElement<dim, spacedim> &fe,
4201  const hp::QCollection<dim - 1> & face_quadrature,
4202  const UpdateFlags update_flags);
4203 
4210  const Quadrature<dim - 1> & face_quadrature,
4211  const UpdateFlags update_flags);
4212 
4220  const hp::QCollection<dim - 1> & face_quadrature,
4221  const UpdateFlags update_flags);
4222 
4229  template <bool level_dof_access>
4230  void
4231  reinit(
4233  const unsigned int face_no,
4234  const unsigned int subface_no);
4235 
4240  template <bool level_dof_access>
4241  void
4242  reinit(
4244  const typename Triangulation<dim, spacedim>::face_iterator & face,
4245  const typename Triangulation<dim, spacedim>::face_iterator &subface);
4246 
4260  void
4262  const unsigned int face_no,
4263  const unsigned int subface_no);
4264 
4284  void
4286  const typename Triangulation<dim, spacedim>::face_iterator &face,
4287  const typename Triangulation<dim, spacedim>::face_iterator &subface);
4288 
4304  get_present_fe_values() const;
4305 
4311  DeclException0(ExcReinitCalledWithBoundaryFace);
4312 
4318  DeclException0(ExcFaceHasNoSubfaces);
4319 
4320 private:
4324  void
4325  initialize(const UpdateFlags update_flags);
4326 
4333  void
4334  do_reinit(const unsigned int face_no, const unsigned int subface_no);
4335 };
4336 
4337 
4338 #ifndef DOXYGEN
4339 
4340 
4341 /*------------------------ Inline functions: namespace FEValuesViews --------*/
4342 
4343 namespace FEValuesViews
4344 {
4345  template <int dim, int spacedim>
4346  inline typename Scalar<dim, spacedim>::value_type
4347  Scalar<dim, spacedim>::value(const unsigned int shape_function,
4348  const unsigned int q_point) const
4349  {
4350  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4351  Assert(
4352  fe_values->update_flags & update_values,
4354  "update_values"))));
4355 
4356  // an adaptation of the FEValuesBase::shape_value_component function
4357  // except that here we know the component as fixed and we have
4358  // pre-computed and cached a bunch of information. See the comments there.
4359  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4360  return fe_values->finite_element_output.shape_values(
4361  shape_function_data[shape_function].row_index, q_point);
4362  else
4363  return 0;
4364  }
4365 
4366 
4367 
4368  template <int dim, int spacedim>
4369  inline typename Scalar<dim, spacedim>::gradient_type
4370  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
4371  const unsigned int q_point) const
4372  {
4373  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4374  Assert(fe_values->update_flags & update_gradients,
4376  "update_gradients")));
4377 
4378  // an adaptation of the FEValuesBase::shape_grad_component
4379  // function except that here we know the component as fixed and we have
4380  // pre-computed and cached a bunch of information. See the comments there.
4381  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4382  return fe_values->finite_element_output
4383  .shape_gradients[shape_function_data[shape_function].row_index]
4384  [q_point];
4385  else
4386  return gradient_type();
4387  }
4388 
4389 
4390 
4391  template <int dim, int spacedim>
4392  inline typename Scalar<dim, spacedim>::hessian_type
4393  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
4394  const unsigned int q_point) const
4395  {
4396  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4397  Assert(fe_values->update_flags & update_hessians,
4399  "update_hessians")));
4400 
4401  // an adaptation of the FEValuesBase::shape_hessian_component
4402  // function except that here we know the component as fixed and we have
4403  // pre-computed and cached a bunch of information. See the comments there.
4404  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4405  return fe_values->finite_element_output
4406  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
4407  else
4408  return hessian_type();
4409  }
4410 
4411 
4412 
4413  template <int dim, int spacedim>
4414  inline typename Scalar<dim, spacedim>::third_derivative_type
4415  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
4416  const unsigned int q_point) const
4417  {
4418  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4419  Assert(fe_values->update_flags & update_3rd_derivatives,
4421  "update_3rd_derivatives")));
4422 
4423  // an adaptation of the FEValuesBase::shape_3rdderivative_component
4424  // function except that here we know the component as fixed and we have
4425  // pre-computed and cached a bunch of information. See the comments there.
4426  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4427  return fe_values->finite_element_output
4428  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
4429  [q_point];
4430  else
4431  return third_derivative_type();
4432  }
4433 
4434 
4435 
4436  template <int dim, int spacedim>
4437  inline typename Vector<dim, spacedim>::value_type
4438  Vector<dim, spacedim>::value(const unsigned int shape_function,
4439  const unsigned int q_point) const
4440  {
4441  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4442  Assert(fe_values->update_flags & update_values,
4444  "update_values")));
4445 
4446  // same as for the scalar case except that we have one more index
4447  const int snc =
4448  shape_function_data[shape_function].single_nonzero_component;
4449  if (snc == -2)
4450  return value_type();
4451  else if (snc != -1)
4452  {
4453  value_type return_value;
4454  return_value[shape_function_data[shape_function]
4455  .single_nonzero_component_index] =
4456  fe_values->finite_element_output.shape_values(snc, q_point);
4457  return return_value;
4458  }
4459  else
4460  {
4461  value_type return_value;
4462  for (unsigned int d = 0; d < dim; ++d)
4463  if (shape_function_data[shape_function]
4464  .is_nonzero_shape_function_component[d])
4465  return_value[d] = fe_values->finite_element_output.shape_values(
4466  shape_function_data[shape_function].row_index[d], q_point);
4467 
4468  return return_value;
4469  }
4470  }
4471 
4472 
4473 
4474  template <int dim, int spacedim>
4475  inline typename Vector<dim, spacedim>::gradient_type
4476  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4477  const unsigned int q_point) const
4478  {
4479  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4480  Assert(fe_values->update_flags & update_gradients,
4482  "update_gradients")));
4483 
4484  // same as for the scalar case except that we have one more index
4485  const int snc =
4486  shape_function_data[shape_function].single_nonzero_component;
4487  if (snc == -2)
4488  return gradient_type();
4489  else if (snc != -1)
4490  {
4491  gradient_type return_value;
4492  return_value[shape_function_data[shape_function]
4493  .single_nonzero_component_index] =
4494  fe_values->finite_element_output.shape_gradients[snc][q_point];
4495  return return_value;
4496  }
4497  else
4498  {
4499  gradient_type return_value;
4500  for (unsigned int d = 0; d < dim; ++d)
4501  if (shape_function_data[shape_function]
4502  .is_nonzero_shape_function_component[d])
4503  return_value[d] =
4504  fe_values->finite_element_output.shape_gradients
4505  [shape_function_data[shape_function].row_index[d]][q_point];
4506 
4507  return return_value;
4508  }
4509  }
4510 
4511 
4512 
4513  template <int dim, int spacedim>
4514  inline typename Vector<dim, spacedim>::divergence_type
4515  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4516  const unsigned int q_point) const
4517  {
4518  // this function works like in the case above
4519  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4520  Assert(fe_values->update_flags & update_gradients,
4522  "update_gradients")));
4523 
4524  // same as for the scalar case except that we have one more index
4525  const int snc =
4526  shape_function_data[shape_function].single_nonzero_component;
4527  if (snc == -2)
4528  return divergence_type();
4529  else if (snc != -1)
4530  return fe_values->finite_element_output
4531  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4532  .single_nonzero_component_index];
4533  else
4534  {
4535  divergence_type return_value = 0;
4536  for (unsigned int d = 0; d < dim; ++d)
4537  if (shape_function_data[shape_function]
4538  .is_nonzero_shape_function_component[d])
4539  return_value +=
4540  fe_values->finite_element_output.shape_gradients
4541  [shape_function_data[shape_function].row_index[d]][q_point][d];
4542 
4543  return return_value;
4544  }
4545  }
4546 
4547 
4548 
4549  template <int dim, int spacedim>
4550  inline typename Vector<dim, spacedim>::curl_type
4551  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4552  const unsigned int q_point) const
4553  {
4554  // this function works like in the case above
4555 
4556  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4557  Assert(fe_values->update_flags & update_gradients,
4559  "update_gradients")));
4560  // same as for the scalar case except that we have one more index
4561  const int snc =
4562  shape_function_data[shape_function].single_nonzero_component;
4563 
4564  if (snc == -2)
4565  return curl_type();
4566 
4567  else
4568  switch (dim)
4569  {
4570  case 1:
4571  {
4572  Assert(false,
4573  ExcMessage(
4574  "Computing the curl in 1d is not a useful operation"));
4575  return curl_type();
4576  }
4577 
4578  case 2:
4579  {
4580  if (snc != -1)
4581  {
4582  curl_type return_value;
4583 
4584  // the single nonzero component can only be zero or one in 2d
4585  if (shape_function_data[shape_function]
4586  .single_nonzero_component_index == 0)
4587  return_value[0] =
4588  -1.0 * fe_values->finite_element_output
4589  .shape_gradients[snc][q_point][1];
4590  else
4591  return_value[0] = fe_values->finite_element_output
4592  .shape_gradients[snc][q_point][0];
4593 
4594  return return_value;
4595  }
4596 
4597  else
4598  {
4599  curl_type return_value;
4600 
4601  return_value[0] = 0.0;
4602 
4603  if (shape_function_data[shape_function]
4604  .is_nonzero_shape_function_component[0])
4605  return_value[0] -=
4606  fe_values->finite_element_output
4607  .shape_gradients[shape_function_data[shape_function]
4608  .row_index[0]][q_point][1];
4609 
4610  if (shape_function_data[shape_function]
4611  .is_nonzero_shape_function_component[1])
4612  return_value[0] +=
4613  fe_values->finite_element_output
4614  .shape_gradients[shape_function_data[shape_function]
4615  .row_index[1]][q_point][0];
4616 
4617  return return_value;
4618  }
4619  }
4620 
4621  case 3:
4622  {
4623  if (snc != -1)
4624  {
4625  curl_type return_value;
4626 
4627  switch (shape_function_data[shape_function]
4628  .single_nonzero_component_index)
4629  {
4630  case 0:
4631  {
4632  return_value[0] = 0;
4633  return_value[1] = fe_values->finite_element_output
4634  .shape_gradients[snc][q_point][2];
4635  return_value[2] =
4636  -1.0 * fe_values->finite_element_output
4637  .shape_gradients[snc][q_point][1];
4638  return return_value;
4639  }
4640 
4641  case 1:
4642  {
4643  return_value[0] =
4644  -1.0 * fe_values->finite_element_output
4645  .shape_gradients[snc][q_point][2];
4646  return_value[1] = 0;
4647  return_value[2] = fe_values->finite_element_output
4648  .shape_gradients[snc][q_point][0];
4649  return return_value;
4650  }
4651 
4652  default:
4653  {
4654  return_value[0] = fe_values->finite_element_output
4655  .shape_gradients[snc][q_point][1];
4656  return_value[1] =
4657  -1.0 * fe_values->finite_element_output
4658  .shape_gradients[snc][q_point][0];
4659  return_value[2] = 0;
4660  return return_value;
4661  }
4662  }
4663  }
4664 
4665  else
4666  {
4667  curl_type return_value;
4668 
4669  for (unsigned int i = 0; i < dim; ++i)
4670  return_value[i] = 0.0;
4671 
4672  if (shape_function_data[shape_function]
4673  .is_nonzero_shape_function_component[0])
4674  {
4675  return_value[1] +=
4676  fe_values->finite_element_output
4677  .shape_gradients[shape_function_data[shape_function]
4678  .row_index[0]][q_point][2];
4679  return_value[2] -=
4680  fe_values->finite_element_output
4681  .shape_gradients[shape_function_data[shape_function]
4682  .row_index[0]][q_point][1];
4683  }
4684 
4685  if (shape_function_data[shape_function]
4686  .is_nonzero_shape_function_component[1])
4687  {
4688  return_value[0] -=
4689  fe_values->finite_element_output
4690  .shape_gradients[shape_function_data[shape_function]
4691  .row_index[1]][q_point][2];
4692  return_value[2] +=
4693  fe_values->finite_element_output
4694  .shape_gradients[shape_function_data[shape_function]
4695  .row_index[1]][q_point][0];
4696  }
4697 
4698  if (shape_function_data[shape_function]
4699  .is_nonzero_shape_function_component[2])
4700  {
4701  return_value[0] +=
4702  fe_values->finite_element_output
4703  .shape_gradients[shape_function_data[shape_function]
4704  .row_index[2]][q_point][1];
4705  return_value[1] -=
4706  fe_values->finite_element_output
4707  .shape_gradients[shape_function_data[shape_function]
4708  .row_index[2]][q_point][0];
4709  }
4710 
4711  return return_value;
4712  }
4713  }
4714  }
4715  // should not end up here
4716  Assert(false, ExcInternalError());
4717  return curl_type();
4718  }
4719 
4720 
4721 
4722  template <int dim, int spacedim>
4723  inline typename Vector<dim, spacedim>::hessian_type
4724  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
4725  const unsigned int q_point) const
4726  {
4727  // this function works like in the case above
4728  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4729  Assert(fe_values->update_flags & update_hessians,
4731  "update_hessians")));
4732 
4733  // same as for the scalar case except that we have one more index
4734  const int snc =
4735  shape_function_data[shape_function].single_nonzero_component;
4736  if (snc == -2)
4737  return hessian_type();
4738  else if (snc != -1)
4739  {
4740  hessian_type return_value;
4741  return_value[shape_function_data[shape_function]
4742  .single_nonzero_component_index] =
4743  fe_values->finite_element_output.shape_hessians[snc][q_point];
4744  return return_value;
4745  }
4746  else
4747  {
4748  hessian_type return_value;
4749  for (unsigned int d = 0; d < dim; ++d)
4750  if (shape_function_data[shape_function]
4751  .is_nonzero_shape_function_component[d])
4752  return_value[d] =
4753  fe_values->finite_element_output.shape_hessians
4754  [shape_function_data[shape_function].row_index[d]][q_point];
4755 
4756  return return_value;
4757  }
4758  }
4759 
4760 
4761 
4762  template <int dim, int spacedim>
4763  inline typename Vector<dim, spacedim>::third_derivative_type
4764  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
4765  const unsigned int q_point) const
4766  {
4767  // this function works like in the case above
4768  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4769  Assert(fe_values->update_flags & update_3rd_derivatives,
4771  "update_3rd_derivatives")));
4772 
4773  // same as for the scalar case except that we have one more index
4774  const int snc =
4775  shape_function_data[shape_function].single_nonzero_component;
4776  if (snc == -2)
4777  return third_derivative_type();
4778  else if (snc != -1)
4779  {
4780  third_derivative_type return_value;
4781  return_value[shape_function_data[shape_function]
4782  .single_nonzero_component_index] =
4783  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
4784  return return_value;
4785  }
4786  else
4787  {
4788  third_derivative_type return_value;
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[d] =
4793  fe_values->finite_element_output.shape_3rd_derivatives
4794  [shape_function_data[shape_function].row_index[d]][q_point];
4795 
4796  return return_value;
4797  }
4798  }
4799 
4800 
4801 
4802  namespace internal
4803  {
4808  inline ::SymmetricTensor<2, 1>
4809  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
4810  {
4811  AssertIndexRange(n, 1);
4812  (void)n;
4813 
4814  return {{t[0]}};
4815  }
4816 
4817 
4818 
4819  inline ::SymmetricTensor<2, 2>
4820  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
4821  {
4822  switch (n)
4823  {
4824  case 0:
4825  {
4826  return {{t[0], 0, t[1] / 2}};
4827  }
4828  case 1:
4829  {
4830  return {{0, t[1], t[0] / 2}};
4831  }
4832  default:
4833  {
4834  AssertIndexRange(n, 2);
4835  return {};
4836  }
4837  }
4838  }
4839 
4840 
4841 
4842  inline ::SymmetricTensor<2, 3>
4843  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
4844  {
4845  switch (n)
4846  {
4847  case 0:
4848  {
4849  return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
4850  }
4851  case 1:
4852  {
4853  return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
4854  }
4855  case 2:
4856  {
4857  return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
4858  }
4859  default:
4860  {
4861  AssertIndexRange(n, 3);
4862  return {};
4863  }
4864  }
4865  }
4866  } // namespace internal
4867 
4868 
4869 
4870  template <int dim, int spacedim>
4871  inline typename Vector<dim, spacedim>::symmetric_gradient_type
4872  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
4873  const unsigned int q_point) const
4874  {
4875  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4876  Assert(fe_values->update_flags & update_gradients,
4878  "update_gradients")));
4879 
4880  // same as for the scalar case except that we have one more index
4881  const int snc =
4882  shape_function_data[shape_function].single_nonzero_component;
4883  if (snc == -2)
4884  return symmetric_gradient_type();
4885  else if (snc != -1)
4886  return internal::symmetrize_single_row(
4887  shape_function_data[shape_function].single_nonzero_component_index,
4888  fe_values->finite_element_output.shape_gradients[snc][q_point]);
4889  else
4890  {
4891  gradient_type return_value;
4892  for (unsigned int d = 0; d < dim; ++d)
4893  if (shape_function_data[shape_function]
4894  .is_nonzero_shape_function_component[d])
4895  return_value[d] =
4896  fe_values->finite_element_output.shape_gradients
4897  [shape_function_data[shape_function].row_index[d]][q_point];
4898 
4899  return symmetrize(return_value);
4900  }
4901  }
4902 
4903 
4904 
4905  template <int dim, int spacedim>
4907  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
4908  const unsigned int q_point) const
4909  {
4910  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4911  Assert(fe_values->update_flags & update_values,
4913  "update_values")));
4914 
4915  // similar to the vector case where we have more then one index and we need
4916  // to convert between unrolled and component indexing for tensors
4917  const int snc =
4918  shape_function_data[shape_function].single_nonzero_component;
4919 
4920  if (snc == -2)
4921  {
4922  // shape function is zero for the selected components
4923  return value_type();
4924  }
4925  else if (snc != -1)
4926  {
4927  value_type return_value;
4928  const unsigned int comp =
4929  shape_function_data[shape_function].single_nonzero_component_index;
4930  return_value[value_type::unrolled_to_component_indices(comp)] =
4931  fe_values->finite_element_output.shape_values(snc, q_point);
4932  return return_value;
4933  }
4934  else
4935  {
4936  value_type return_value;
4937  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
4938  if (shape_function_data[shape_function]
4939  .is_nonzero_shape_function_component[d])
4940  return_value[value_type::unrolled_to_component_indices(d)] =
4941  fe_values->finite_element_output.shape_values(
4942  shape_function_data[shape_function].row_index[d], q_point);
4943  return return_value;
4944  }
4945  }
4946 
4947 
4948 
4949  template <int dim, int spacedim>
4952  const unsigned int shape_function,
4953  const unsigned int q_point) const
4954  {
4955  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4956  Assert(fe_values->update_flags & update_gradients,
4958  "update_gradients")));
4959 
4960  const int snc =
4961  shape_function_data[shape_function].single_nonzero_component;
4962 
4963  if (snc == -2)
4964  {
4965  // shape function is zero for the selected components
4966  return divergence_type();
4967  }
4968  else if (snc != -1)
4969  {
4970  // we have a single non-zero component when the symmetric tensor is
4971  // represented in unrolled form. this implies we potentially have
4972  // two non-zero components when represented in component form! we
4973  // will only have one non-zero entry if the non-zero component lies on
4974  // the diagonal of the tensor.
4975  //
4976  // the divergence of a second-order tensor is a first order tensor.
4977  //
4978  // assume the second-order tensor is A with components A_{ij}. then
4979  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
4980  // entries in the tensorial representation. define the
4981  // divergence as:
4982  // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
4983  // (which is incidentally also
4984  // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
4985  // In both cases, a sum is implied.
4986  //
4987  // Now, we know the nonzero component in unrolled form: it is indicated
4988  // by 'snc'. we can figure out which tensor components belong to this:
4989  const unsigned int comp =
4990  shape_function_data[shape_function].single_nonzero_component_index;
4991  const unsigned int ii =
4992  value_type::unrolled_to_component_indices(comp)[0];
4993  const unsigned int jj =
4994  value_type::unrolled_to_component_indices(comp)[1];
4995 
4996  // given the form of the divergence above, if ii=jj there is only a
4997  // single nonzero component of the full tensor and the gradient
4998  // equals
4999  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
5000  // all other entries of 'b' are zero
5001  //
5002  // on the other hand, if ii!=jj, then there are two nonzero entries in
5003  // the full tensor and
5004  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
5005  // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
5006  // again, all other entries of 'b' are zero
5007  const ::Tensor<1, spacedim> &phi_grad =
5008  fe_values->finite_element_output.shape_gradients[snc][q_point];
5009 
5010  divergence_type return_value;
5011  return_value[ii] = phi_grad[jj];
5012 
5013  if (ii != jj)
5014  return_value[jj] = phi_grad[ii];
5015 
5016  return return_value;
5017  }
5018  else
5019  {
5020  Assert(false, ExcNotImplemented());
5021  divergence_type return_value;
5022  return return_value;
5023  }
5024  }
5025 
5026 
5027 
5028  template <int dim, int spacedim>
5029  inline typename Tensor<2, dim, spacedim>::value_type
5030  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
5031  const unsigned int q_point) const
5032  {
5033  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5034  Assert(fe_values->update_flags & update_values,
5036  "update_values")));
5037 
5038  // similar to the vector case where we have more then one index and we need
5039  // to convert between unrolled and component indexing for tensors
5040  const int snc =
5041  shape_function_data[shape_function].single_nonzero_component;
5042 
5043  if (snc == -2)
5044  {
5045  // shape function is zero for the selected components
5046  return value_type();
5047  }
5048  else if (snc != -1)
5049  {
5050  value_type return_value;
5051  const unsigned int comp =
5052  shape_function_data[shape_function].single_nonzero_component_index;
5053  const TableIndices<2> indices =
5055  return_value[indices] =
5056  fe_values->finite_element_output.shape_values(snc, q_point);
5057  return return_value;
5058  }
5059  else
5060  {
5061  value_type return_value;
5062  for (unsigned int d = 0; d < dim * dim; ++d)
5063  if (shape_function_data[shape_function]
5064  .is_nonzero_shape_function_component[d])
5065  {
5066  const TableIndices<2> indices =
5068  return_value[indices] =
5069  fe_values->finite_element_output.shape_values(
5070  shape_function_data[shape_function].row_index[d], q_point);
5071  }
5072  return return_value;
5073  }
5074  }
5075 
5076 
5077 
5078  template <int dim, int spacedim>
5080  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
5081  const unsigned int q_point) const
5082  {
5083  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5084  Assert(fe_values->update_flags & update_gradients,
5086  "update_gradients")));
5087 
5088  const int snc =
5089  shape_function_data[shape_function].single_nonzero_component;
5090 
5091  if (snc == -2)
5092  {
5093  // shape function is zero for the selected components
5094  return divergence_type();
5095  }
5096  else if (snc != -1)
5097  {
5098  // we have a single non-zero component when the tensor is
5099  // represented in unrolled form.
5100  //
5101  // the divergence of a second-order tensor is a first order tensor.
5102  //
5103  // assume the second-order tensor is A with components A_{ij},
5104  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
5105  //
5106  // Now, we know the nonzero component in unrolled form: it is indicated
5107  // by 'snc'. we can figure out which tensor components belong to this:
5108  const unsigned int comp =
5109  shape_function_data[shape_function].single_nonzero_component_index;
5110  const TableIndices<2> indices =
5112  const unsigned int ii = indices[0];
5113  const unsigned int jj = indices[1];
5114 
5115  const ::Tensor<1, spacedim> &phi_grad =
5116  fe_values->finite_element_output.shape_gradients[snc][q_point];
5117 
5118  divergence_type return_value;
5119  // note that we contract \nabla from the right
5120  return_value[ii] = phi_grad[jj];
5121 
5122  return return_value;
5123  }
5124  else
5125  {
5126  Assert(false, ExcNotImplemented());
5127  divergence_type return_value;
5128  return return_value;
5129  }
5130  }
5131 
5132 
5133 
5134  template <int dim, int spacedim>
5136  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
5137  const unsigned int q_point) const
5138  {
5139  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5140  Assert(fe_values->update_flags & update_gradients,
5142  "update_gradients")));
5143 
5144  const int snc =
5145  shape_function_data[shape_function].single_nonzero_component;
5146 
5147  if (snc == -2)
5148  {
5149  // shape function is zero for the selected components
5150  return gradient_type();
5151  }
5152  else if (snc != -1)
5153  {
5154  // we have a single non-zero component when the tensor is
5155  // represented in unrolled form.
5156  //
5157  // the gradient of a second-order tensor is a third order tensor.
5158  //
5159  // assume the second-order tensor is A with components A_{ij},
5160  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
5161  //
5162  // Now, we know the nonzero component in unrolled form: it is indicated
5163  // by 'snc'. we can figure out which tensor components belong to this:
5164  const unsigned int comp =
5165  shape_function_data[shape_function].single_nonzero_component_index;
5166  const TableIndices<2> indices =
5168  const unsigned int ii = indices[0];
5169  const unsigned int jj = indices[1];
5170 
5171  const ::Tensor<1, spacedim> &phi_grad =
5172  fe_values->finite_element_output.shape_gradients[snc][q_point];
5173 
5174  gradient_type return_value;
5175  return_value[ii][jj] = phi_grad;
5176 
5177  return return_value;
5178  }
5179  else
5180  {
5181  Assert(false, ExcNotImplemented());
5182  gradient_type return_value;
5183  return return_value;
5184  }
5185  }
5186 
5187 } // namespace FEValuesViews
5188 
5189 
5190 
5191 /*---------------------- Inline functions: FEValuesBase ---------------------*/
5192 
5193 
5194 
5195 template <int dim, int spacedim>
5197  operator[](const FEValuesExtractors::Scalar &scalar) const
5198 {
5199  AssertIndexRange(scalar.component, fe_values_views_cache.scalars.size());
5200 
5201  return fe_values_views_cache.scalars[scalar.component];
5202 }
5203 
5204 
5205 
5206 template <int dim, int spacedim>
5208  operator[](const FEValuesExtractors::Vector &vector) const
5209 {
5211  fe_values_views_cache.vectors.size());
5212 
5213  return fe_values_views_cache.vectors[vector.first_vector_component];
5214 }
5215 
5216 
5217 
5218 template <int dim, int spacedim>
5222 {
5223  Assert(
5224  tensor.first_tensor_component <
5225  fe_values_views_cache.symmetric_second_order_tensors.size(),
5227  0,
5228  fe_values_views_cache.symmetric_second_order_tensors.size()));
5229 
5230  return fe_values_views_cache
5231  .symmetric_second_order_tensors[tensor.first_tensor_component];
5232 }
5233 
5234 
5235 
5236 template <int dim, int spacedim>
5239  operator[](const FEValuesExtractors::Tensor<2> &tensor) const
5240 {
5242  fe_values_views_cache.second_order_tensors.size());
5243 
5244  return fe_values_views_cache
5245  .second_order_tensors[tensor.first_tensor_component];
5246 }
5247 
5248 
5249 
5250 template <int dim, int spacedim>
5251 inline const double &
5252 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
5253  const unsigned int j) const
5254 {
5255  AssertIndexRange(i, fe->n_dofs_per_cell());
5256  Assert(this->update_flags & update_values,
5257  ExcAccessToUninitializedField("update_values"));
5258  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5259  Assert(present_cell.get() != nullptr,
5260  ExcMessage("FEValues object is not reinit'ed to any cell"));
5261  // if the entire FE is primitive,
5262  // then we can take a short-cut:
5263  if (fe->is_primitive())
5264  return this->finite_element_output.shape_values(i, j);
5265  else
5266  {
5267  // otherwise, use the mapping
5268  // between shape function
5269  // numbers and rows. note that
5270  // by the assertions above, we
5271  // know that this particular
5272  // shape function is primitive,
5273  // so we can call
5274  // system_to_component_index
5275  const unsigned int row =
5276  this->finite_element_output
5277  .shape_function_to_row_table[i * fe->n_components() +
5278  fe->system_to_component_index(i).first];
5279  return this->finite_element_output.shape_values(row, j);
5280  }
5281 }
5282 
5283 
5284 
5285 template <int dim, int spacedim>
5286 inline double
5288  const unsigned int i,
5289  const unsigned int j,
5290  const unsigned int component) const
5291 {
5292  AssertIndexRange(i, fe->n_dofs_per_cell());
5293  Assert(this->update_flags & update_values,
5294  ExcAccessToUninitializedField("update_values"));
5295  AssertIndexRange(component, fe->n_components());
5296  Assert(present_cell.get() != nullptr,
5297  ExcMessage("FEValues object is not reinit'ed to any cell"));
5298 
5299  // check whether the shape function
5300  // is non-zero at all within
5301  // this component:
5302  if (fe->get_nonzero_components(i)[component] == false)
5303  return 0;
5304 
5305  // look up the right row in the
5306  // table and take the data from
5307  // there
5308  const unsigned int row =
5309  this->finite_element_output
5310  .shape_function_to_row_table[i * fe->n_components() + component];
5311  return this->finite_element_output.shape_values(row, j);
5312 }
5313 
5314 
5315 
5316 template <int dim, int spacedim>
5317 inline const Tensor<1, spacedim> &
5318 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
5319  const unsigned int j) const
5320 {
5321  AssertIndexRange(i, fe->n_dofs_per_cell());
5322  Assert(this->update_flags & update_gradients,
5323  ExcAccessToUninitializedField("update_gradients"));
5324  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5325  Assert(present_cell.get() != nullptr,
5326  ExcMessage("FEValues object is not reinit'ed to any cell"));
5327  // if the entire FE is primitive,
5328  // then we can take a short-cut:
5329  if (fe->is_primitive())
5330  return this->finite_element_output.shape_gradients[i][j];
5331  else
5332  {
5333  // otherwise, use the mapping
5334  // between shape function
5335  // numbers and rows. note that
5336  // by the assertions above, we
5337  // know that this particular
5338  // shape function is primitive,
5339  // so we can call
5340  // system_to_component_index
5341  const unsigned int row =
5342  this->finite_element_output
5343  .shape_function_to_row_table[i * fe->n_components() +
5344  fe->system_to_component_index(i).first];
5345  return this->finite_element_output.shape_gradients[row][j];
5346  }
5347 }
5348 
5349 
5350 
5351 template <int dim, int spacedim>
5352 inline Tensor<1, spacedim>
5354  const unsigned int i,
5355  const unsigned int j,
5356  const unsigned int component) const
5357 {
5358  AssertIndexRange(i, fe->n_dofs_per_cell());
5359  Assert(this->update_flags & update_gradients,
5360  ExcAccessToUninitializedField("update_gradients"));
5361  AssertIndexRange(component, fe->n_components());
5362  Assert(present_cell.get() != nullptr,
5363  ExcMessage("FEValues object is not reinit'ed to any cell"));
5364  // check whether the shape function
5365  // is non-zero at all within
5366  // this component:
5367  if (fe->get_nonzero_components(i)[component] == false)
5368  return Tensor<1, spacedim>();
5369 
5370  // look up the right row in the
5371  // table and take the data from
5372  // there
5373  const unsigned int row =
5374  this->finite_element_output
5375  .shape_function_to_row_table[i * fe->n_components() + component];
5376  return this->finite_element_output.shape_gradients[row][j];
5377 }
5378 
5379 
5380 
5381 template <int dim, int spacedim>
5382 inline const Tensor<2, spacedim> &
5383 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
5384  const unsigned int j) const
5385 {
5386  AssertIndexRange(i, fe->n_dofs_per_cell());
5387  Assert(this->update_flags & update_hessians,
5388  ExcAccessToUninitializedField("update_hessians"));
5389  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5390  Assert(present_cell.get() != nullptr,
5391  ExcMessage("FEValues object is not reinit'ed to any cell"));
5392  // if the entire FE is primitive,
5393  // then we can take a short-cut:
5394  if (fe->is_primitive())
5395  return this->finite_element_output.shape_hessians[i][j];
5396  else
5397  {
5398  // otherwise, use the mapping
5399  // between shape function
5400  // numbers and rows. note that
5401  // by the assertions above, we
5402  // know that this particular
5403  // shape function is primitive,
5404  // so we can call
5405  // system_to_component_index
5406  const unsigned int row =
5407  this->finite_element_output
5408  .shape_function_to_row_table[i * fe->n_components() +
5409  fe->system_to_component_index(i).first];
5410  return this->finite_element_output.shape_hessians[row][j];
5411  }
5412 }
5413 
5414 
5415 
5416 template <int dim, int spacedim>
5417 inline Tensor<2, spacedim>
5419  const unsigned int i,
5420  const unsigned int j,
5421  const unsigned int component) const
5422 {
5423  AssertIndexRange(i, fe->n_dofs_per_cell());
5424  Assert(this->update_flags & update_hessians,
5425  ExcAccessToUninitializedField("update_hessians"));
5426  AssertIndexRange(component, fe->n_components());
5427  Assert(present_cell.get() != nullptr,
5428  ExcMessage("FEValues object is not reinit'ed to any cell"));
5429  // check whether the shape function
5430  // is non-zero at all within
5431  // this component:
5432  if (fe->get_nonzero_components(i)[component] == false)
5433  return Tensor<2, spacedim>();
5434 
5435  // look up the right row in the
5436  // table and take the data from
5437  // there
5438  const unsigned int row =
5439  this->finite_element_output
5440  .shape_function_to_row_table[i * fe->n_components() + component];
5441  return this->finite_element_output.shape_hessians[row][j];
5442 }
5443 
5444 
5445 
5446 template <int dim, int spacedim>
5447 inline const Tensor<3, spacedim> &
5449  const unsigned int j) const
5450 {
5451  AssertIndexRange(i, fe->n_dofs_per_cell());
5452  Assert(this->update_flags & update_3rd_derivatives,
5453  ExcAccessToUninitializedField("update_3rd_derivatives"));
5454  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5455  Assert(present_cell.get() != nullptr,
5456  ExcMessage("FEValues object is not reinit'ed to any cell"));
5457  // if the entire FE is primitive,
5458  // then we can take a short-cut:
5459  if (fe->is_primitive())
5460  return this->finite_element_output.shape_3rd_derivatives[i][j];
5461  else
5462  {
5463  // otherwise, use the mapping
5464  // between shape function
5465  // numbers and rows. note that
5466  // by the assertions above, we
5467  // know that this particular
5468  // shape function is primitive,
5469  // so we can call
5470  // system_to_component_index
5471  const unsigned int row =
5472  this->finite_element_output
5473  .shape_function_to_row_table[i * fe->n_components() +
5474  fe->system_to_component_index(i).first];
5475  return this->finite_element_output.shape_3rd_derivatives[row][j];
5476  }
5477 }
5478 
5479 
5480 
5481 template <int dim, int spacedim>
5482 inline Tensor<3, spacedim>
5484  const unsigned int i,
5485  const unsigned int j,
5486  const unsigned int component) const
5487 {
5488  AssertIndexRange(i, fe->n_dofs_per_cell());
5489  Assert(this->update_flags & update_3rd_derivatives,
5490  ExcAccessToUninitializedField("update_3rd_derivatives"));
5491  AssertIndexRange(component, fe->n_components());
5492  Assert(present_cell.get() != nullptr,
5493  ExcMessage("FEValues object is not reinit'ed to any cell"));
5494  // check whether the shape function
5495  // is non-zero at all within
5496  // this component:
5497  if (fe->get_nonzero_components(i)[component] == false)
5498  return Tensor<3, spacedim>();
5499 
5500  // look up the right row in the
5501  // table and take the data from
5502  // there
5503  const unsigned int row =
5504  this->finite_element_output
5505  .shape_function_to_row_table[i * fe->n_components() + component];
5506  return this->finite_element_output.shape_3rd_derivatives[row][j];
5507 }
5508 
5509 
5510 
5511 template <int dim, int spacedim>
5512 inline const FiniteElement<dim, spacedim> &
5514 {
5515  return *fe;
5516 }
5517 
5518 
5519 
5520 template <int dim, int spacedim>
5521 inline const Mapping<dim, spacedim> &
5523 {
5524  return *mapping;
5525 }
5526 
5527 
5528 
5529 template <int dim, int spacedim>
5530 inline UpdateFlags
5532 {
5533  return this->update_flags;
5534 }
5535 
5536 
5537 
5538 template <int dim, int spacedim>
5539 inline const std::vector<Point<spacedim>> &
5541 {
5542  Assert(this->update_flags & update_quadrature_points,
5543  ExcAccessToUninitializedField("update_quadrature_points"));
5544  Assert(present_cell.get() != nullptr,
5545  ExcMessage("FEValues object is not reinit'ed to any cell"));
5546  return this->mapping_output.quadrature_points;
5547 }
5548 
5549 
5550 
5551 template <int dim, int spacedim>
5552 inline const std::vector<double> &
5554 {
5555  Assert(this->update_flags & update_JxW_values,
5556  ExcAccessToUninitializedField("update_JxW_values"));
5557  Assert(present_cell.get() != nullptr,
5558  ExcMessage("FEValues object is not reinit'ed to any cell"));
5559  return this->mapping_output.JxW_values;
5560 }
5561 
5562 
5563 
5564 template <int dim, int spacedim>
5565 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5567 {
5568  Assert(this->update_flags & update_jacobians,
5569  ExcAccessToUninitializedField("update_jacobians"));
5570  Assert(present_cell.get() != nullptr,
5571  ExcMessage("FEValues object is not reinit'ed to any cell"));
5572  return this->mapping_output.jacobians;
5573 }
5574 
5575 
5576 
5577 template <int dim, int spacedim>
5578 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5580 {
5581  Assert(this->update_flags & update_jacobian_grads,
5582  ExcAccessToUninitializedField("update_jacobians_grads"));
5583  Assert(present_cell.get() != nullptr,
5584  ExcMessage("FEValues object is not reinit'ed to any cell"));
5585  return this->mapping_output.jacobian_grads;
5586 }
5587 
5588 
5589 
5590 template <int dim, int spacedim>
5591 inline const Tensor<3, spacedim> &
5593  const unsigned int i) const
5594 {
5595  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5596  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5597  Assert(present_cell.get() != nullptr,
5598  ExcMessage("FEValues object is not reinit'ed to any cell"));
5599  return this->mapping_output.jacobian_pushed_forward_grads[i];
5600 }
5601 
5602 
5603 
5604 template <int dim, int spacedim>
5605 inline const std::vector<Tensor<3, spacedim>> &
5607 {
5608  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5609  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5610  Assert(present_cell.get() != nullptr,
5611  ExcMessage("FEValues object is not reinit'ed to any cell"));
5612  return this->mapping_output.jacobian_pushed_forward_grads;
5613 }
5614 
5615 
5616 
5617 template <int dim, int spacedim>
5618 inline const DerivativeForm<3, dim, spacedim> &
5619 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5620 {
5621  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5622  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5623  Assert(present_cell.get() != nullptr,
5624  ExcMessage("FEValues object is not reinit'ed to any cell"));
5625  return this->mapping_output.jacobian_2nd_derivatives[i];
5626 }
5627 
5628 
5629 
5630 template <int dim, int spacedim>
5631 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5633 {
5634  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5635  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5636  Assert(present_cell.get() != nullptr,
5637  ExcMessage("FEValues object is not reinit'ed to any cell"));
5638  return this->mapping_output.jacobian_2nd_derivatives;
5639 }
5640 
5641 
5642 
5643 template <int dim, int spacedim>
5644 inline const Tensor<4, spacedim> &
5646  const unsigned int i) const
5647 {
5650  "update_jacobian_pushed_forward_2nd_derivatives"));
5651  Assert(present_cell.get() != nullptr,
5652  ExcMessage("FEValues object is not reinit'ed to any cell"));
5653  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5654 }
5655 
5656 
5657 
5658 template <int dim, int spacedim>
5659 inline const std::vector<Tensor<4, spacedim>> &
5661 {
5662  Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
5664  "update_jacobian_pushed_forward_2nd_derivatives"));
5665  Assert(present_cell.get() != nullptr,
5666  ExcMessage("FEValues object is not reinit'ed to any cell"));
5667  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5668 }
5669 
5670 
5671 
5672 template <int dim, int spacedim>
5673 inline const DerivativeForm<4, dim, spacedim> &
5674 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5675 {
5676  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5677  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5678  Assert(present_cell.get() != nullptr,
5679  ExcMessage("FEValues object is not reinit'ed to any cell"));
5680  return this->mapping_output.jacobian_3rd_derivatives[i];
5681 }
5682 
5683 
5684 
5685 template <int dim, int spacedim>
5686 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5688 {
5689  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5690  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5691  Assert(present_cell.get() != nullptr,
5692  ExcMessage("FEValues object is not reinit'ed to any cell"));
5693  return this->mapping_output.jacobian_3rd_derivatives;
5694 }
5695 
5696 
5697 
5698 template <int dim, int spacedim>
5699 inline const Tensor<5, spacedim> &
5701  const unsigned int i) const
5702 {
5705  "update_jacobian_pushed_forward_3rd_derivatives"));
5706  Assert(present_cell.get() != nullptr,
5707  ExcMessage("FEValues object is not reinit'ed to any cell"));
5708  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
5709 }
5710 
5711 
5712 
5713 template <int dim, int spacedim>
5714 inline const std::vector<Tensor<5, spacedim>> &
5716 {
5717  Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
5719  "update_jacobian_pushed_forward_3rd_derivatives"));
5720  Assert(present_cell.get() != nullptr,
5721  ExcMessage("FEValues object is not reinit'ed to any cell"));
5722  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
5723 }
5724 
5725 
5726 
5727 template <int dim, int spacedim>
5728 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
5730 {
5731  Assert(this->update_flags & update_inverse_jacobians,
5732  ExcAccessToUninitializedField("update_inverse_jacobians"));
5733  Assert(present_cell.get() != nullptr,
5734  ExcMessage("FEValues object is not reinit'ed to any cell"));
5735  return this->mapping_output.inverse_jacobians;
5736 }
5737 
5738 
5739 
5740 template <int dim, int spacedim>
5743 {
5744  return {0U, dofs_per_cell};
5745 }
5746 
5747 
5748 
5749 template <int dim, int spacedim>
5752  const unsigned int start_dof_index) const
5753 {
5754  Assert(start_dof_index <= dofs_per_cell,
5755  ExcIndexRange(start_dof_index, 0, dofs_per_cell + 1));
5756  return {start_dof_index, dofs_per_cell};
5757 }
5758 
5759 
5760 
5761 template <int dim, int spacedim>
5764  const unsigned int end_dof_index) const
5765 {
5766  Assert(end_dof_index < dofs_per_cell,
5767  ExcIndexRange(end_dof_index, 0, dofs_per_cell));
5768  return {0U, end_dof_index + 1};
5769 }
5770 
5771 
5772 
5773 template <int dim, int spacedim>
5776 {
5777  return {0U, n_quadrature_points};
5778 }
5779 
5780 
5781 
5782 template <int dim, int spacedim>
5783 inline const Point<spacedim> &
5784 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
5785 {
5786  Assert(this->update_flags & update_quadrature_points,
5787  ExcAccessToUninitializedField("update_quadrature_points"));
5788  AssertIndexRange(i, this->mapping_output.quadrature_points.size());
5789  Assert(present_cell.get() != nullptr,
5790  ExcMessage("FEValues object is not reinit'ed to any cell"));
5791 
5792  return this->mapping_output.quadrature_points[i];
5793 }
5794 
5795 
5796 
5797 template <int dim, int spacedim>
5798 inline double
5799 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
5800 {
5801  Assert(this->update_flags & update_JxW_values,
5802  ExcAccessToUninitializedField("update_JxW_values"));
5803  AssertIndexRange(i, this->mapping_output.JxW_values.size());
5804  Assert(present_cell.get() != nullptr,
5805  ExcMessage("FEValues object is not reinit'ed to any cell"));
5806 
5807  return this->mapping_output.JxW_values[i];
5808 }
5809 
5810 
5811 
5812 template <int dim, int spacedim>
5813 inline const DerivativeForm<1, dim, spacedim> &
5814 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
5815 {
5816  Assert(this->update_flags & update_jacobians,
5817  ExcAccessToUninitializedField("update_jacobians"));
5818  AssertIndexRange(i, this->mapping_output.jacobians.size());
5819  Assert(present_cell.get() != nullptr,
5820  ExcMessage("FEValues object is not reinit'ed to any cell"));
5821 
5822  return this->mapping_output.jacobians[i];
5823 }
5824 
5825 
5826 
5827 template <int dim, int spacedim>
5828 inline const DerivativeForm<2, dim, spacedim> &
5829 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
5830 {
5831  Assert(this->update_flags & update_jacobian_grads,
5832  ExcAccessToUninitializedField("update_jacobians_grads"));
5833  AssertIndexRange(i, this->mapping_output.jacobian_grads.size());
5834  Assert(present_cell.get() != nullptr,
5835  ExcMessage("FEValues object is not reinit'ed to any cell"));
5836 
5837  return this->mapping_output.jacobian_grads[i];
5838 }
5839 
5840 
5841 
5842 template <int dim, int spacedim>
5843 inline const DerivativeForm<1, spacedim, dim> &
5844 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
5845 {
5846  Assert(this->update_flags & update_inverse_jacobians,
5847  ExcAccessToUninitializedField("update_inverse_jacobians"));
5848  AssertIndexRange(i, this->mapping_output.inverse_jacobians.size());
5849  Assert(present_cell.get() != nullptr,
5850  ExcMessage("FEValues object is not reinit'ed to any cell"));
5851 
5852  return this->mapping_output.inverse_jacobians[i];
5853 }
5854 
5855 
5856 
5857 template <int dim, int spacedim>
5858 inline const Tensor<1, spacedim> &
5859 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
5860 {
5861  Assert(this->update_flags & update_normal_vectors,
5863  "update_normal_vectors")));
5864  AssertIndexRange(i, this->mapping_output.normal_vectors.size());
5865  Assert(present_cell.get() != nullptr,
5866  ExcMessage("FEValues object is not reinit'ed to any cell"));
5867 
5868  return this->mapping_output.normal_vectors[i];
5869 }
5870 
5871 
5872 
5873 /*--------------------- Inline functions: FEValues --------------------------*/
5874 
5875 
5876 template <int dim, int spacedim>
5877 inline const Quadrature<dim> &
5879 {
5880  return quadrature;
5881 }
5882 
5883 
5884 
5885 template <int dim, int spacedim>
5886 inline const FEValues<dim, spacedim> &
5888 {
5889  return *this;
5890 }
5891 
5892 
5893 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
5894 
5895 
5896 template <int dim, int spacedim>
5897 inline unsigned int
5899 {
5900  return present_face_index;
5901 }
5902 
5903 
5904 /*----------------------- Inline functions: FE*FaceValues -------------------*/
5905 
5906 template <int dim, int spacedim>
5907 inline const Quadrature<dim - 1> &
5909 {
5910  return quadrature[quadrature.size() == 1 ? 0 : present_face_no];
5911 }
5912 
5913 
5914 
5915 template <int dim, int spacedim>
5916 inline const FEFaceValues<dim, spacedim> &
5918 {
5919  return *this;
5920 }
5921 
5922 
5923 
5924 template <int dim, int spacedim>
5925 inline const FESubfaceValues<dim, spacedim> &
5927 {
5928  return *this;
5929 }
5930 
5931 
5932 
5933 template <int dim, int spacedim>
5934 inline const Tensor<1, spacedim> &
5935 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
5936 {
5937  AssertIndexRange(i, this->mapping_output.boundary_forms.size());
5938  Assert(this->update_flags & update_boundary_forms,
5940  "update_boundary_forms")));
5941 
5942  return this->mapping_output.boundary_forms[i];
5943 }
5944 
5945 #endif // DOXYGEN
5946 
5948 
5949 #endif
Transformed quadrature weights.
Shape function values.
typename ProductType< Number, typename Vector< dim, spacedim >::curl_type >::type curl_type
Definition: fe_values.h:721
typename FEValuesViews::View< dim, spacedim, Extractor >::template OutputType< NumberType > OutputType
Definition: scratch_data.h:47
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:3674
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:3706
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:1666
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1350
unsigned int present_face_no
Definition: fe_values.h:3969
unsigned int present_face_index
Definition: fe_values.h:3975
std::vector<::FEValuesViews::Vector< dim, spacedim > > vectors
Definition: fe_values.h:2040
typename ::internal::FEValuesViews::ViewType< dim, spacedim, Extractor >::type View
Definition: fe_values.h:2062
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1591
constexpr SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
static ::ExceptionBase & ExcAccessToUninitializedField()
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:561
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:2204
typename ::internal::CurlType< spacedim >::type curl_type
Definition: fe_values.h:652
const unsigned int component
Definition: fe_values.h:567
Volume element.
typename ::FEValuesViews::SymmetricTensor< rank, dim, spacedim > type
Definition: fe_values.h:2022
const ReferenceCell::internal::Info::Base & get_cell(const ReferenceCell::Type &type)
#define AssertIndexRange(index, range)
Definition: exceptions.h:1691
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:214
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:3590
UpdateFlags get_update_flags() const
static const char U
Transformed quadrature points.
std::vector< double > get_quadrature_points(const unsigned int n)
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:190
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:3642
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int quadrature_point) const
const DerivativeForm< 3, dim, spacedim > & jacobian_2nd_derivative(const unsigned int quadrature_point) const
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition: fe_values.h:3721
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()
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:2042
typename ProductType< Number, typename Vector< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:689
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:198
const hp::QCollection< dim - 1 > quadrature
Definition: fe_values.h:3980
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1650
static ::ExceptionBase & ExcMessage(std::string arg1)
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:737
static constexpr TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:515
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_ending_at(const unsigned int end_dof_index) const
Tensor< 3, spacedim > shape_3rd_derivative_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
Third derivatives of shape functions.
const FEValues< dim, spacedim > & get_present_fe_values() const
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1954
unsigned int get_face_index() const
#define Assert(cond, exc)
Definition: exceptions.h:1466
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:1285
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
const Quadrature< dim > quadrature
Definition: fe_values.h:3857
const unsigned int first_vector_component
Definition: fe_values.h:1280
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:380
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition: fe_values.h:3650
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int quadrature_point) const
typename Tensor< rank_ - 1, dim, Number >::tensor_type value_type
Definition: tensor.h:483
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:1358
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:2039
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:2186
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
boost::signals2::connection tria_listener_mesh_transform
Definition: fe_values.h:3615
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:681
typename ::FEValuesViews::Tensor< rank, dim, spacedim > type
Definition: fe_values.h:2008
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
Definition: tensor.h:449
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:697
typename ProductType< Number, typename Vector< dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:705
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:379
typename ProductType< Number, typename Vector< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:729
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:1658
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:1943
Definition: fe.h:38
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1580
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:1274
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:3606
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3682
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:713
::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:3657
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:2197
std::vector<::FEValuesViews::Tensor< 2, dim, spacedim > > second_order_tensors
Definition: fe_values.h:2044
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:206
typename ProductType< Number, typename Scalar< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:222
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:572
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
UpdateFlags update_flags
Definition: fe_values.h:3688
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:3666
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