Reference documentation for deal.II version Git 2e98021fd4 2020-10-20 17:25:20 -0400
\(\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 
41 #include <algorithm>
42 #include <memory>
43 #include <type_traits>
44 
45 
46 // dummy include in order to have the
47 // definition of PetscScalar available
48 // without including other PETSc stuff
49 #ifdef DEAL_II_WITH_PETSC
50 # include <petsc.h>
51 #endif
52 
54 
55 // Forward declaration
56 #ifndef DOXYGEN
57 template <int dim, int spacedim = dim>
58 class FEValuesBase;
59 #endif
60 
61 namespace internal
62 {
67  template <int dim, class NumberType = double>
68  struct CurlType;
69 
76  template <class NumberType>
77  struct CurlType<1, NumberType>
78  {
80  };
81 
88  template <class NumberType>
89  struct CurlType<2, NumberType>
90  {
92  };
93 
100  template <class NumberType>
101  struct CurlType<3, NumberType>
102  {
104  };
105 } // namespace internal
106 
107 
108 
130 namespace FEValuesViews
131 {
143  template <int dim, int spacedim = dim>
144  class Scalar
145  {
146  public:
153 
160 
167 
174 
179  template <typename Number>
180  struct OutputType
181  {
186  using value_type =
187  typename ProductType<Number,
189 
194  using gradient_type = typename ProductType<
195  Number,
197 
202  using laplacian_type =
203  typename ProductType<Number,
205 
210  using hessian_type = typename ProductType<
211  Number,
213 
218  using third_derivative_type = typename ProductType<
219  Number,
221  };
222 
228  {
238 
247  unsigned int row_index;
248  };
249 
253  Scalar();
254 
260  Scalar(const FEValuesBase<dim, spacedim> &fe_values_base,
261  const unsigned int component);
262 
267  Scalar &
268  operator=(const Scalar<dim, spacedim> &) = delete;
269 
283  value_type
284  value(const unsigned int shape_function, const unsigned int q_point) const;
285 
297  gradient(const unsigned int shape_function,
298  const unsigned int q_point) const;
299 
311  hessian(const unsigned int shape_function,
312  const unsigned int q_point) const;
313 
325  third_derivative(const unsigned int shape_function,
326  const unsigned int q_point) const;
327 
345  template <class InputVector>
346  void
347  get_function_values(
348  const InputVector &fe_function,
349  std::vector<typename ProductType<value_type,
350  typename InputVector::value_type>::type>
351  &values) const;
352 
375  template <class InputVector>
376  void
377  get_function_values_from_local_dof_values(
378  const InputVector &dof_values,
379  std::vector<
381  &values) const;
382 
400  template <class InputVector>
401  void
402  get_function_gradients(
403  const InputVector &fe_function,
404  std::vector<typename ProductType<gradient_type,
405  typename InputVector::value_type>::type>
406  &gradients) const;
407 
411  template <class InputVector>
412  void
413  get_function_gradients_from_local_dof_values(
414  const InputVector &dof_values,
415  std::vector<
417  &gradients) const;
418 
436  template <class InputVector>
437  void
438  get_function_hessians(
439  const InputVector &fe_function,
440  std::vector<typename ProductType<hessian_type,
441  typename InputVector::value_type>::type>
442  &hessians) const;
443 
447  template <class InputVector>
448  void
449  get_function_hessians_from_local_dof_values(
450  const InputVector &dof_values,
451  std::vector<
453  &hessians) const;
454 
455 
474  template <class InputVector>
475  void
476  get_function_laplacians(
477  const InputVector &fe_function,
478  std::vector<typename ProductType<value_type,
479  typename InputVector::value_type>::type>
480  &laplacians) const;
481 
485  template <class InputVector>
486  void
487  get_function_laplacians_from_local_dof_values(
488  const InputVector &dof_values,
489  std::vector<
491  &laplacians) const;
492 
493 
512  template <class InputVector>
513  void
514  get_function_third_derivatives(
515  const InputVector &fe_function,
516  std::vector<typename ProductType<third_derivative_type,
517  typename InputVector::value_type>::type>
518  &third_derivatives) const;
519 
523  template <class InputVector>
524  void
525  get_function_third_derivatives_from_local_dof_values(
526  const InputVector & dof_values,
527  std::vector<typename OutputType<typename InputVector::value_type>::
528  third_derivative_type> &third_derivatives) const;
529 
530 
531  private:
536 
541  const unsigned int component;
542 
546  std::vector<ShapeFunctionData> shape_function_data;
547  };
548 
549 
550 
580  template <int dim, int spacedim = dim>
581  class Vector
582  {
583  public:
590 
600 
612 
619 
626  using curl_type = typename ::internal::CurlType<spacedim>::type;
627 
634 
641 
646  template <typename Number>
647  struct OutputType
648  {
653  using value_type =
654  typename ProductType<Number,
656 
661  using gradient_type = typename ProductType<
662  Number,
664 
669  using symmetric_gradient_type = typename ProductType<
670  Number,
672 
677  using divergence_type = typename ProductType<
678  Number,
680 
685  using laplacian_type =
686  typename ProductType<Number,
688 
693  using curl_type =
694  typename ProductType<Number,
696 
701  using hessian_type = typename ProductType<
702  Number,
704 
709  using third_derivative_type = typename ProductType<
710  Number,
712  };
713 
719  {
728  bool is_nonzero_shape_function_component[spacedim];
729 
739  unsigned int row_index[spacedim];
740 
751  };
752 
756  Vector();
757 
766  Vector(const FEValuesBase<dim, spacedim> &fe_values_base,
767  const unsigned int first_vector_component);
768 
773  Vector &
774  operator=(const Vector<dim, spacedim> &) = delete;
775 
792  value_type
793  value(const unsigned int shape_function, const unsigned int q_point) const;
794 
809  gradient(const unsigned int shape_function,
810  const unsigned int q_point) const;
811 
828  symmetric_gradient(const unsigned int shape_function,
829  const unsigned int q_point) const;
830 
842  divergence(const unsigned int shape_function,
843  const unsigned int q_point) const;
844 
865  curl_type
866  curl(const unsigned int shape_function, const unsigned int q_point) const;
867 
879  hessian(const unsigned int shape_function,
880  const unsigned int q_point) const;
881 
893  third_derivative(const unsigned int shape_function,
894  const unsigned int q_point) const;
895 
913  template <class InputVector>
914  void
915  get_function_values(
916  const InputVector &fe_function,
917  std::vector<typename ProductType<value_type,
918  typename InputVector::value_type>::type>
919  &values) const;
920 
943  template <class InputVector>
944  void
945  get_function_values_from_local_dof_values(
946  const InputVector &dof_values,
947  std::vector<
949  &values) const;
950 
968  template <class InputVector>
969  void
970  get_function_gradients(
971  const InputVector &fe_function,
972  std::vector<typename ProductType<gradient_type,
973  typename InputVector::value_type>::type>
974  &gradients) const;
975 
979  template <class InputVector>
980  void
981  get_function_gradients_from_local_dof_values(
982  const InputVector &dof_values,
983  std::vector<
985  &gradients) const;
986 
1010  template <class InputVector>
1011  void
1012  get_function_symmetric_gradients(
1013  const InputVector &fe_function,
1014  std::vector<typename ProductType<symmetric_gradient_type,
1015  typename InputVector::value_type>::type>
1016  &symmetric_gradients) const;
1017 
1021  template <class InputVector>
1022  void
1023  get_function_symmetric_gradients_from_local_dof_values(
1024  const InputVector & dof_values,
1025  std::vector<typename OutputType<typename InputVector::value_type>::
1026  symmetric_gradient_type> &symmetric_gradients) const;
1027 
1046  template <class InputVector>
1047  void
1048  get_function_divergences(
1049  const InputVector &fe_function,
1050  std::vector<typename ProductType<divergence_type,
1051  typename InputVector::value_type>::type>
1052  &divergences) const;
1053 
1057  template <class InputVector>
1058  void
1059  get_function_divergences_from_local_dof_values(
1060  const InputVector &dof_values,
1061  std::vector<
1063  &divergences) const;
1064 
1083  template <class InputVector>
1084  void
1085  get_function_curls(
1086  const InputVector &fe_function,
1087  std::vector<
1089  &curls) const;
1090 
1094  template <class InputVector>
1095  void
1096  get_function_curls_from_local_dof_values(
1097  const InputVector &dof_values,
1098  std::vector<
1100  &curls) const;
1101 
1119  template <class InputVector>
1120  void
1121  get_function_hessians(
1122  const InputVector &fe_function,
1123  std::vector<typename ProductType<hessian_type,
1124  typename InputVector::value_type>::type>
1125  &hessians) const;
1126 
1130  template <class InputVector>
1131  void
1132  get_function_hessians_from_local_dof_values(
1133  const InputVector &dof_values,
1134  std::vector<
1136  &hessians) const;
1137 
1156  template <class InputVector>
1157  void
1158  get_function_laplacians(
1159  const InputVector &fe_function,
1160  std::vector<typename ProductType<value_type,
1161  typename InputVector::value_type>::type>
1162  &laplacians) const;
1163 
1167  template <class InputVector>
1168  void
1169  get_function_laplacians_from_local_dof_values(
1170  const InputVector &dof_values,
1171  std::vector<
1173  &laplacians) const;
1174 
1193  template <class InputVector>
1194  void
1195  get_function_third_derivatives(
1196  const InputVector &fe_function,
1197  std::vector<typename ProductType<third_derivative_type,
1198  typename InputVector::value_type>::type>
1199  &third_derivatives) const;
1200 
1204  template <class InputVector>
1205  void
1206  get_function_third_derivatives_from_local_dof_values(
1207  const InputVector & dof_values,
1208  std::vector<typename OutputType<typename InputVector::value_type>::
1209  third_derivative_type> &third_derivatives) const;
1210 
1211  private:
1216 
1221  const unsigned int first_vector_component;
1222 
1226  std::vector<ShapeFunctionData> shape_function_data;
1227  };
1228 
1229 
1230  template <int rank, int dim, int spacedim = dim>
1232 
1255  template <int dim, int spacedim>
1256  class SymmetricTensor<2, dim, spacedim>
1257  {
1258  public:
1266 
1277 
1282  template <typename Number>
1283  struct OutputType
1284  {
1289  using value_type = typename ProductType<
1290  Number,
1292 
1297  using divergence_type = typename ProductType<
1298  Number,
1300  };
1301 
1306  struct ShapeFunctionData
1307  {
1316  bool is_nonzero_shape_function_component
1317  [value_type::n_independent_components];
1318 
1328  unsigned int row_index[value_type::n_independent_components];
1329 
1339 
1344  };
1345 
1349  SymmetricTensor();
1350 
1360  SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1361  const unsigned int first_tensor_component);
1362 
1367  SymmetricTensor &
1368  operator=(const SymmetricTensor<2, dim, spacedim> &) = delete;
1369 
1387  value_type
1388  value(const unsigned int shape_function, const unsigned int q_point) const;
1389 
1404  divergence(const unsigned int shape_function,
1405  const unsigned int q_point) const;
1406 
1424  template <class InputVector>
1425  void
1426  get_function_values(
1427  const InputVector &fe_function,
1428  std::vector<typename ProductType<value_type,
1429  typename InputVector::value_type>::type>
1430  &values) const;
1431 
1454  template <class InputVector>
1455  void
1456  get_function_values_from_local_dof_values(
1457  const InputVector &dof_values,
1458  std::vector<
1460  &values) const;
1461 
1483  template <class InputVector>
1484  void
1485  get_function_divergences(
1486  const InputVector &fe_function,
1487  std::vector<typename ProductType<divergence_type,
1488  typename InputVector::value_type>::type>
1489  &divergences) const;
1490 
1494  template <class InputVector>
1495  void
1496  get_function_divergences_from_local_dof_values(
1497  const InputVector &dof_values,
1498  std::vector<
1500  &divergences) const;
1501 
1502  private:
1507 
1512  const unsigned int first_tensor_component;
1513 
1517  std::vector<ShapeFunctionData> shape_function_data;
1518  };
1519 
1520 
1521  template <int rank, int dim, int spacedim = dim>
1522  class Tensor;
1523 
1542  template <int dim, int spacedim>
1543  class Tensor<2, dim, spacedim>
1544  {
1545  public:
1551 
1556 
1562 
1567  template <typename Number>
1568  struct OutputType
1569  {
1574  using value_type = typename ProductType<
1575  Number,
1577 
1582  using divergence_type = typename ProductType<
1583  Number,
1585 
1590  using gradient_type = typename ProductType<
1591  Number,
1593  };
1594 
1599  struct ShapeFunctionData
1600  {
1609  bool is_nonzero_shape_function_component
1610  [value_type::n_independent_components];
1611 
1621  unsigned int row_index[value_type::n_independent_components];
1622 
1632 
1637  };
1638 
1642  Tensor();
1643 
1653  Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1654  const unsigned int first_tensor_component);
1655 
1656 
1661  Tensor &
1662  operator=(const Tensor<2, dim, spacedim> &) = delete;
1663 
1680  value_type
1681  value(const unsigned int shape_function, const unsigned int q_point) const;
1682 
1697  divergence(const unsigned int shape_function,
1698  const unsigned int q_point) const;
1699 
1714  gradient(const unsigned int shape_function,
1715  const unsigned int q_point) const;
1716 
1734  template <class InputVector>
1735  void
1736  get_function_values(
1737  const InputVector &fe_function,
1738  std::vector<typename ProductType<value_type,
1739  typename InputVector::value_type>::type>
1740  &values) const;
1741 
1764  template <class InputVector>
1765  void
1766  get_function_values_from_local_dof_values(
1767  const InputVector &dof_values,
1768  std::vector<
1770  &values) const;
1771 
1793  template <class InputVector>
1794  void
1795  get_function_divergences(
1796  const InputVector &fe_function,
1797  std::vector<typename ProductType<divergence_type,
1798  typename InputVector::value_type>::type>
1799  &divergences) const;
1800 
1804  template <class InputVector>
1805  void
1806  get_function_divergences_from_local_dof_values(
1807  const InputVector &dof_values,
1808  std::vector<
1810  &divergences) const;
1811 
1828  template <class InputVector>
1829  void
1830  get_function_gradients(
1831  const InputVector &fe_function,
1832  std::vector<typename ProductType<gradient_type,
1833  typename InputVector::value_type>::type>
1834  &gradients) const;
1835 
1839  template <class InputVector>
1840  void
1841  get_function_gradients_from_local_dof_values(
1842  const InputVector &dof_values,
1843  std::vector<
1845  &gradients) const;
1846 
1847  private:
1852 
1857  const unsigned int first_tensor_component;
1858 
1862  std::vector<ShapeFunctionData> shape_function_data;
1863  };
1864 
1865 } // namespace FEValuesViews
1866 
1867 
1868 namespace internal
1869 {
1870  namespace FEValuesViews
1871  {
1876  template <int dim, int spacedim, typename Extractor>
1877  struct ViewType
1878  {};
1879 
1887  template <int dim, int spacedim>
1888  struct ViewType<dim, spacedim, FEValuesExtractors::Scalar>
1889  {
1890  using type = typename ::FEValuesViews::Scalar<dim, spacedim>;
1891  };
1892 
1900  template <int dim, int spacedim>
1901  struct ViewType<dim, spacedim, FEValuesExtractors::Vector>
1902  {
1903  using type = typename ::FEValuesViews::Vector<dim, spacedim>;
1904  };
1905 
1913  template <int dim, int spacedim, int rank>
1914  struct ViewType<dim, spacedim, FEValuesExtractors::Tensor<rank>>
1915  {
1916  using type = typename ::FEValuesViews::Tensor<rank, dim, spacedim>;
1917  };
1918 
1926  template <int dim, int spacedim, int rank>
1927  struct ViewType<dim, spacedim, FEValuesExtractors::SymmetricTensor<rank>>
1928  {
1929  using type =
1930  typename ::FEValuesViews::SymmetricTensor<rank, dim, spacedim>;
1931  };
1932 
1940  template <int dim, int spacedim>
1941  struct Cache
1942  {
1947  std::vector<::FEValuesViews::Scalar<dim, spacedim>> scalars;
1948  std::vector<::FEValuesViews::Vector<dim, spacedim>> vectors;
1949  std::vector<::FEValuesViews::SymmetricTensor<2, dim, spacedim>>
1951  std::vector<::FEValuesViews::Tensor<2, dim, spacedim>>
1953 
1957  Cache(const FEValuesBase<dim, spacedim> &fe_values);
1958  };
1959  } // namespace FEValuesViews
1960 } // namespace internal
1961 
1962 namespace FEValuesViews
1963 {
1968  template <int dim, int spacedim, typename Extractor>
1969  using View = typename ::internal::FEValuesViews::
1970  ViewType<dim, spacedim, Extractor>::type;
1971 } // namespace FEValuesViews
1972 
1973 
2073 template <int dim, int spacedim>
2074 class FEValuesBase : public Subscriptor
2075 {
2076 public:
2080  static const unsigned int dimension = dim;
2081 
2085  static const unsigned int space_dimension = spacedim;
2086 
2090  const unsigned int n_quadrature_points;
2091 
2097  const unsigned int dofs_per_cell;
2098 
2099 
2107  FEValuesBase(const unsigned int n_q_points,
2108  const unsigned int dofs_per_cell,
2109  const UpdateFlags update_flags,
2110  const Mapping<dim, spacedim> & mapping,
2111  const FiniteElement<dim, spacedim> &fe);
2112 
2117  FEValuesBase &
2118  operator=(const FEValuesBase &) = delete;
2119 
2124  FEValuesBase(const FEValuesBase &) = delete;
2125 
2129  virtual ~FEValuesBase() override;
2130 
2131 
2135 
2136 
2157  const double &
2158  shape_value(const unsigned int function_no,
2159  const unsigned int point_no) const;
2160 
2181  double
2182  shape_value_component(const unsigned int function_no,
2183  const unsigned int point_no,
2184  const unsigned int component) const;
2185 
2211  const Tensor<1, spacedim> &
2212  shape_grad(const unsigned int function_no,
2213  const unsigned int quadrature_point) const;
2214 
2232  shape_grad_component(const unsigned int function_no,
2233  const unsigned int point_no,
2234  const unsigned int component) const;
2235 
2255  const Tensor<2, spacedim> &
2256  shape_hessian(const unsigned int function_no,
2257  const unsigned int point_no) const;
2258 
2276  shape_hessian_component(const unsigned int function_no,
2277  const unsigned int point_no,
2278  const unsigned int component) const;
2279 
2299  const Tensor<3, spacedim> &
2300  shape_3rd_derivative(const unsigned int function_no,
2301  const unsigned int point_no) const;
2302 
2320  shape_3rd_derivative_component(const unsigned int function_no,
2321  const unsigned int point_no,
2322  const unsigned int component) const;
2323 
2325 
2327 
2364  template <class InputVector>
2365  void
2366  get_function_values(
2367  const InputVector & fe_function,
2368  std::vector<typename InputVector::value_type> &values) const;
2369 
2383  template <class InputVector>
2384  void
2385  get_function_values(
2386  const InputVector & fe_function,
2387  std::vector<Vector<typename InputVector::value_type>> &values) const;
2388 
2407  template <class InputVector>
2408  void
2409  get_function_values(
2410  const InputVector & fe_function,
2412  std::vector<typename InputVector::value_type> & values) const;
2413 
2435  template <class InputVector>
2436  void
2437  get_function_values(
2438  const InputVector & fe_function,
2440  std::vector<Vector<typename InputVector::value_type>> &values) const;
2441 
2442 
2473  template <class InputVector>
2474  void
2475  get_function_values(
2476  const InputVector & fe_function,
2478  ArrayView<std::vector<typename InputVector::value_type>> values,
2479  const bool quadrature_points_fastest) const;
2480 
2482 
2484 
2521  template <class InputVector>
2522  void
2523  get_function_gradients(
2524  const InputVector &fe_function,
2526  &gradients) const;
2527 
2544  template <class InputVector>
2545  void
2546  get_function_gradients(
2547  const InputVector &fe_function,
2548  std::vector<
2550  &gradients) const;
2551 
2558  template <class InputVector>
2559  void
2560  get_function_gradients(
2561  const InputVector & fe_function,
2564  &gradients) const;
2565 
2572  template <class InputVector>
2573  void
2574  get_function_gradients(
2575  const InputVector & fe_function,
2577  ArrayView<
2579  gradients,
2580  bool quadrature_points_fastest = false) const;
2581 
2583 
2587 
2625  template <class InputVector>
2626  void
2627  get_function_hessians(
2628  const InputVector &fe_function,
2630  &hessians) const;
2631 
2649  template <class InputVector>
2650  void
2651  get_function_hessians(
2652  const InputVector &fe_function,
2653  std::vector<
2655  & hessians,
2656  bool quadrature_points_fastest = false) const;
2657 
2662  template <class InputVector>
2663  void
2664  get_function_hessians(
2665  const InputVector & fe_function,
2668  &hessians) const;
2669 
2676  template <class InputVector>
2677  void
2678  get_function_hessians(
2679  const InputVector & fe_function,
2681  ArrayView<
2683  hessians,
2684  bool quadrature_points_fastest = false) const;
2685 
2726  template <class InputVector>
2727  void
2728  get_function_laplacians(
2729  const InputVector & fe_function,
2730  std::vector<typename InputVector::value_type> &laplacians) const;
2731 
2751  template <class InputVector>
2752  void
2753  get_function_laplacians(
2754  const InputVector & fe_function,
2755  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2756 
2763  template <class InputVector>
2764  void
2765  get_function_laplacians(
2766  const InputVector & fe_function,
2768  std::vector<typename InputVector::value_type> & laplacians) const;
2769 
2776  template <class InputVector>
2777  void
2778  get_function_laplacians(
2779  const InputVector & fe_function,
2781  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2782 
2789  template <class InputVector>
2790  void
2791  get_function_laplacians(
2792  const InputVector & fe_function,
2794  std::vector<std::vector<typename InputVector::value_type>> &laplacians,
2795  bool quadrature_points_fastest = false) const;
2796 
2798 
2800 
2839  template <class InputVector>
2840  void
2841  get_function_third_derivatives(
2842  const InputVector &fe_function,
2844  &third_derivatives) const;
2845 
2864  template <class InputVector>
2865  void
2866  get_function_third_derivatives(
2867  const InputVector &fe_function,
2868  std::vector<
2870  & third_derivatives,
2871  bool quadrature_points_fastest = false) const;
2872 
2877  template <class InputVector>
2878  void
2879  get_function_third_derivatives(
2880  const InputVector & fe_function,
2883  &third_derivatives) const;
2884 
2891  template <class InputVector>
2892  void
2893  get_function_third_derivatives(
2894  const InputVector & fe_function,
2896  ArrayView<
2898  third_derivatives,
2899  bool quadrature_points_fastest = false) const;
2901 
2903 
2904 
2929  dof_indices() const;
2930 
2964  dof_indices_starting_at(const unsigned int start_dof_index) const;
2965 
2997  dof_indices_ending_at(const unsigned int end_dof_index) const;
2998 
3000 
3002 
3003 
3025  quadrature_point_indices() const;
3026 
3032  const Point<spacedim> &
3033  quadrature_point(const unsigned int q) const;
3034 
3040  const std::vector<Point<spacedim>> &
3041  get_quadrature_points() const;
3042 
3058  double
3059  JxW(const unsigned int quadrature_point) const;
3060 
3064  const std::vector<double> &
3065  get_JxW_values() const;
3066 
3074  jacobian(const unsigned int quadrature_point) const;
3075 
3082  const std::vector<DerivativeForm<1, dim, spacedim>> &
3083  get_jacobians() const;
3084 
3093  jacobian_grad(const unsigned int quadrature_point) const;
3094 
3101  const std::vector<DerivativeForm<2, dim, spacedim>> &
3102  get_jacobian_grads() const;
3103 
3112  const Tensor<3, spacedim> &
3113  jacobian_pushed_forward_grad(const unsigned int quadrature_point) const;
3114 
3121  const std::vector<Tensor<3, spacedim>> &
3122  get_jacobian_pushed_forward_grads() const;
3123 
3132  jacobian_2nd_derivative(const unsigned int quadrature_point) const;
3133 
3140  const std::vector<DerivativeForm<3, dim, spacedim>> &
3141  get_jacobian_2nd_derivatives() const;
3142 
3152  const Tensor<4, spacedim> &
3153  jacobian_pushed_forward_2nd_derivative(
3154  const unsigned int quadrature_point) const;
3155 
3162  const std::vector<Tensor<4, spacedim>> &
3163  get_jacobian_pushed_forward_2nd_derivatives() const;
3164 
3174  jacobian_3rd_derivative(const unsigned int quadrature_point) const;
3175 
3182  const std::vector<DerivativeForm<4, dim, spacedim>> &
3183  get_jacobian_3rd_derivatives() const;
3184 
3194  const Tensor<5, spacedim> &
3195  jacobian_pushed_forward_3rd_derivative(
3196  const unsigned int quadrature_point) const;
3197 
3204  const std::vector<Tensor<5, spacedim>> &
3205  get_jacobian_pushed_forward_3rd_derivatives() const;
3206 
3214  inverse_jacobian(const unsigned int quadrature_point) const;
3215 
3222  const std::vector<DerivativeForm<1, spacedim, dim>> &
3223  get_inverse_jacobians() const;
3224 
3244  const Tensor<1, spacedim> &
3245  normal_vector(const unsigned int i) const;
3246 
3254  const std::vector<Tensor<1, spacedim>> &
3255  get_normal_vectors() const;
3256 
3258 
3260 
3261 
3271  operator[](const FEValuesExtractors::Scalar &scalar) const;
3272 
3282  operator[](const FEValuesExtractors::Vector &vector) const;
3283 
3294  operator[](const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
3295 
3296 
3306  operator[](const FEValuesExtractors::Tensor<2> &tensor) const;
3307 
3309 
3311 
3312 
3316  const Mapping<dim, spacedim> &
3317  get_mapping() const;
3318 
3323  get_fe() const;
3324 
3328  UpdateFlags
3329  get_update_flags() const;
3330 
3335  get_cell() const;
3336 
3343  get_cell_similarity() const;
3344 
3349  std::size_t
3350  memory_consumption() const;
3352 
3353 
3362  std::string,
3363  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3364  << "object for which this kind of information has not been computed. What "
3365  << "information these objects compute is determined by the update_* flags you "
3366  << "pass to the constructor. Here, the operation you are attempting requires "
3367  << "the <" << arg1
3368  << "> flag to be set, but it was apparently not specified "
3369  << "upon construction.");
3370 
3378  ExcFEDontMatch,
3379  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3380  "to the DoFHandler that provided the cell iterator do not match.");
3386  DeclException1(ExcShapeFunctionNotPrimitive,
3387  int,
3388  << "The shape function with index " << arg1
3389  << " is not primitive, i.e. it is vector-valued and "
3390  << "has more than one non-zero vector component. This "
3391  << "function cannot be called for these shape functions. "
3392  << "Maybe you want to use the same function with the "
3393  << "_component suffix?");
3394 
3403  "The given FiniteElement is not a primitive element but the requested operation "
3404  "only works for those. See FiniteElement::is_primitive() for more information.");
3405 
3406 protected:
3437  class CellIteratorBase;
3438 
3443  template <typename CI>
3446 
3452  std::unique_ptr<const CellIteratorBase> present_cell;
3453 
3461  boost::signals2::connection tria_listener_refinement;
3462 
3470  boost::signals2::connection tria_listener_mesh_transform;
3471 
3477  void
3478  invalidate_present_cell();
3479 
3489  void
3490  maybe_invalidate_previous_present_cell(
3491  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3492 
3498 
3504  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3506 
3513 
3514 
3522 
3528  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3530 
3536  spacedim>
3538 
3539 
3544 
3553  UpdateFlags
3554  compute_update_flags(const UpdateFlags update_flags) const;
3555 
3562 
3568  void
3569  check_cell_similarity(
3570  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3571 
3572 private:
3577 
3578  // Make the view classes friends of this class, since they access internal
3579  // data.
3580  template <int, int>
3582  template <int, int>
3584  template <int, int, int>
3586  template <int, int, int>
3588 };
3589 
3590 
3591 
3601 template <int dim, int spacedim = dim>
3602 class FEValues : public FEValuesBase<dim, spacedim>
3603 {
3604 public:
3609  static const unsigned int integral_dimension = dim;
3610 
3615  FEValues(const Mapping<dim, spacedim> & mapping,
3616  const FiniteElement<dim, spacedim> &fe,
3617  const Quadrature<dim> & quadrature,
3618  const UpdateFlags update_flags);
3619 
3626  const Quadrature<dim> & quadrature,
3627  const UpdateFlags update_flags);
3628 
3635  template <bool level_dof_access>
3636  void
3637  reinit(
3639 
3653  void
3655 
3660  const Quadrature<dim> &
3661  get_quadrature() const;
3662 
3667  std::size_t
3668  memory_consumption() const;
3669 
3684  const FEValues<dim, spacedim> &
3685  get_present_fe_values() const;
3686 
3687 private:
3692 
3696  void
3697  initialize(const UpdateFlags update_flags);
3698 
3705  void
3706  do_reinit();
3707 };
3708 
3709 
3719 template <int dim, int spacedim = dim>
3720 class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
3721 {
3722 public:
3727  static const unsigned int integral_dimension = dim - 1;
3728 
3740  FEFaceValuesBase(const unsigned int n_q_points,
3741  const unsigned int dofs_per_cell,
3742  const UpdateFlags update_flags,
3743  const Mapping<dim, spacedim> & mapping,
3744  const FiniteElement<dim, spacedim> &fe,
3745  const Quadrature<dim - 1> & quadrature);
3746 
3754  const Tensor<1, spacedim> &
3755  boundary_form(const unsigned int i) const;
3756 
3763  const std::vector<Tensor<1, spacedim>> &
3764  get_boundary_forms() const;
3765 
3770  unsigned int
3771  get_face_index() const;
3772 
3777  const Quadrature<dim - 1> &
3778  get_quadrature() const;
3779 
3784  std::size_t
3785  memory_consumption() const;
3786 
3787 protected:
3792  unsigned int present_face_index;
3793 
3797  const Quadrature<dim - 1> quadrature;
3798 };
3799 
3800 
3801 
3815 template <int dim, int spacedim = dim>
3816 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
3817 {
3818 public:
3823  static const unsigned int dimension = dim;
3824 
3825  static const unsigned int space_dimension = spacedim;
3826 
3831  static const unsigned int integral_dimension = dim - 1;
3832 
3837  FEFaceValues(const Mapping<dim, spacedim> & mapping,
3838  const FiniteElement<dim, spacedim> &fe,
3839  const Quadrature<dim - 1> & quadrature,
3840  const UpdateFlags update_flags);
3841 
3848  const Quadrature<dim - 1> & quadrature,
3849  const UpdateFlags update_flags);
3850 
3855  template <bool level_dof_access>
3856  void
3857  reinit(
3859  const unsigned int face_no);
3860 
3867  template <bool level_dof_access>
3868  void
3869  reinit(
3871  const typename Triangulation<dim, spacedim>::face_iterator & face);
3872 
3886  void
3888  const unsigned int face_no);
3889 
3890  /*
3891  * Reinitialize the gradients, Jacobi determinants, etc for the given face
3892  * on a given cell of type "iterator into a Triangulation object", and the
3893  * given finite element. Since iterators into a triangulation alone only
3894  * convey information about the geometry of a cell, but not about degrees of
3895  * freedom possibly associated with this cell, you will not be able to call
3896  * some functions of this class if they need information about degrees of
3897  * freedom. These functions are, above all, the
3898  * <tt>get_function_value/gradients/hessians/third_derivatives</tt>
3899  * functions. If you want to call these functions, you have to call the @p
3900  * reinit variants that take iterators into DoFHandler or other DoF handler
3901  * type objects.
3902  *
3903  * @note @p face must be one of @p cell's face iterators.
3904  */
3905  void
3907  const typename Triangulation<dim, spacedim>::face_iterator &face);
3908 
3924  get_present_fe_values() const;
3925 
3926 private:
3930  void
3931  initialize(const UpdateFlags update_flags);
3932 
3939  void
3940  do_reinit(const unsigned int face_no);
3941 };
3942 
3943 
3960 template <int dim, int spacedim = dim>
3961 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
3962 {
3963 public:
3967  static const unsigned int dimension = dim;
3968 
3972  static const unsigned int space_dimension = spacedim;
3973 
3978  static const unsigned int integral_dimension = dim - 1;
3979 
3984  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
3985  const FiniteElement<dim, spacedim> &fe,
3986  const Quadrature<dim - 1> & face_quadrature,
3987  const UpdateFlags update_flags);
3988 
3995  const Quadrature<dim - 1> & face_quadrature,
3996  const UpdateFlags update_flags);
3997 
4004  template <bool level_dof_access>
4005  void
4006  reinit(
4008  const unsigned int face_no,
4009  const unsigned int subface_no);
4010 
4015  template <bool level_dof_access>
4016  void
4017  reinit(
4019  const typename Triangulation<dim, spacedim>::face_iterator & face,
4020  const typename Triangulation<dim, spacedim>::face_iterator &subface);
4021 
4035  void
4037  const unsigned int face_no,
4038  const unsigned int subface_no);
4039 
4059  void
4061  const typename Triangulation<dim, spacedim>::face_iterator &face,
4062  const typename Triangulation<dim, spacedim>::face_iterator &subface);
4063 
4079  get_present_fe_values() const;
4080 
4086  DeclException0(ExcReinitCalledWithBoundaryFace);
4087 
4093  DeclException0(ExcFaceHasNoSubfaces);
4094 
4095 private:
4099  void
4100  initialize(const UpdateFlags update_flags);
4101 
4108  void
4109  do_reinit(const unsigned int face_no, const unsigned int subface_no);
4110 };
4111 
4112 
4113 #ifndef DOXYGEN
4114 
4115 
4116 /*------------------------ Inline functions: namespace FEValuesViews --------*/
4117 
4118 namespace FEValuesViews
4119 {
4120  template <int dim, int spacedim>
4121  inline typename Scalar<dim, spacedim>::value_type
4122  Scalar<dim, spacedim>::value(const unsigned int shape_function,
4123  const unsigned int q_point) const
4124  {
4125  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4126  Assert(
4127  fe_values->update_flags & update_values,
4129  "update_values"))));
4130 
4131  // an adaptation of the FEValuesBase::shape_value_component function
4132  // except that here we know the component as fixed and we have
4133  // pre-computed and cached a bunch of information. See the comments there.
4134  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4135  return fe_values->finite_element_output.shape_values(
4136  shape_function_data[shape_function].row_index, q_point);
4137  else
4138  return 0;
4139  }
4140 
4141 
4142 
4143  template <int dim, int spacedim>
4144  inline typename Scalar<dim, spacedim>::gradient_type
4145  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
4146  const unsigned int q_point) const
4147  {
4148  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4149  Assert(fe_values->update_flags & update_gradients,
4151  "update_gradients")));
4152 
4153  // an adaptation of the FEValuesBase::shape_grad_component
4154  // function except that here we know the component as fixed and we have
4155  // pre-computed and cached a bunch of information. See the comments there.
4156  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4157  return fe_values->finite_element_output
4158  .shape_gradients[shape_function_data[shape_function].row_index]
4159  [q_point];
4160  else
4161  return gradient_type();
4162  }
4163 
4164 
4165 
4166  template <int dim, int spacedim>
4167  inline typename Scalar<dim, spacedim>::hessian_type
4168  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
4169  const unsigned int q_point) const
4170  {
4171  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4172  Assert(fe_values->update_flags & update_hessians,
4174  "update_hessians")));
4175 
4176  // an adaptation of the FEValuesBase::shape_hessian_component
4177  // function except that here we know the component as fixed and we have
4178  // pre-computed and cached a bunch of information. See the comments there.
4179  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4180  return fe_values->finite_element_output
4181  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
4182  else
4183  return hessian_type();
4184  }
4185 
4186 
4187 
4188  template <int dim, int spacedim>
4189  inline typename Scalar<dim, spacedim>::third_derivative_type
4190  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
4191  const unsigned int q_point) const
4192  {
4193  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4194  Assert(fe_values->update_flags & update_3rd_derivatives,
4196  "update_3rd_derivatives")));
4197 
4198  // an adaptation of the FEValuesBase::shape_3rdderivative_component
4199  // function except that here we know the component as fixed and we have
4200  // pre-computed and cached a bunch of information. See the comments there.
4201  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4202  return fe_values->finite_element_output
4203  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
4204  [q_point];
4205  else
4206  return third_derivative_type();
4207  }
4208 
4209 
4210 
4211  template <int dim, int spacedim>
4212  inline typename Vector<dim, spacedim>::value_type
4213  Vector<dim, spacedim>::value(const unsigned int shape_function,
4214  const unsigned int q_point) const
4215  {
4216  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4217  Assert(fe_values->update_flags & update_values,
4219  "update_values")));
4220 
4221  // same as for the scalar case except that we have one more index
4222  const int snc =
4223  shape_function_data[shape_function].single_nonzero_component;
4224  if (snc == -2)
4225  return value_type();
4226  else if (snc != -1)
4227  {
4228  value_type return_value;
4229  return_value[shape_function_data[shape_function]
4230  .single_nonzero_component_index] =
4231  fe_values->finite_element_output.shape_values(snc, q_point);
4232  return return_value;
4233  }
4234  else
4235  {
4236  value_type return_value;
4237  for (unsigned int d = 0; d < dim; ++d)
4238  if (shape_function_data[shape_function]
4239  .is_nonzero_shape_function_component[d])
4240  return_value[d] = fe_values->finite_element_output.shape_values(
4241  shape_function_data[shape_function].row_index[d], q_point);
4242 
4243  return return_value;
4244  }
4245  }
4246 
4247 
4248 
4249  template <int dim, int spacedim>
4250  inline typename Vector<dim, spacedim>::gradient_type
4251  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4252  const unsigned int q_point) const
4253  {
4254  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4255  Assert(fe_values->update_flags & update_gradients,
4257  "update_gradients")));
4258 
4259  // same as for the scalar case except that we have one more index
4260  const int snc =
4261  shape_function_data[shape_function].single_nonzero_component;
4262  if (snc == -2)
4263  return gradient_type();
4264  else if (snc != -1)
4265  {
4266  gradient_type return_value;
4267  return_value[shape_function_data[shape_function]
4268  .single_nonzero_component_index] =
4269  fe_values->finite_element_output.shape_gradients[snc][q_point];
4270  return return_value;
4271  }
4272  else
4273  {
4274  gradient_type return_value;
4275  for (unsigned int d = 0; d < dim; ++d)
4276  if (shape_function_data[shape_function]
4277  .is_nonzero_shape_function_component[d])
4278  return_value[d] =
4279  fe_values->finite_element_output.shape_gradients
4280  [shape_function_data[shape_function].row_index[d]][q_point];
4281 
4282  return return_value;
4283  }
4284  }
4285 
4286 
4287 
4288  template <int dim, int spacedim>
4289  inline typename Vector<dim, spacedim>::divergence_type
4290  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4291  const unsigned int q_point) const
4292  {
4293  // this function works like in the case above
4294  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4295  Assert(fe_values->update_flags & update_gradients,
4297  "update_gradients")));
4298 
4299  // same as for the scalar case except that we have one more index
4300  const int snc =
4301  shape_function_data[shape_function].single_nonzero_component;
4302  if (snc == -2)
4303  return divergence_type();
4304  else if (snc != -1)
4305  return fe_values->finite_element_output
4306  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4307  .single_nonzero_component_index];
4308  else
4309  {
4310  divergence_type return_value = 0;
4311  for (unsigned int d = 0; d < dim; ++d)
4312  if (shape_function_data[shape_function]
4313  .is_nonzero_shape_function_component[d])
4314  return_value +=
4315  fe_values->finite_element_output.shape_gradients
4316  [shape_function_data[shape_function].row_index[d]][q_point][d];
4317 
4318  return return_value;
4319  }
4320  }
4321 
4322 
4323 
4324  template <int dim, int spacedim>
4325  inline typename Vector<dim, spacedim>::curl_type
4326  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4327  const unsigned int q_point) const
4328  {
4329  // this function works like in the case above
4330 
4331  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4332  Assert(fe_values->update_flags & update_gradients,
4334  "update_gradients")));
4335  // same as for the scalar case except that we have one more index
4336  const int snc =
4337  shape_function_data[shape_function].single_nonzero_component;
4338 
4339  if (snc == -2)
4340  return curl_type();
4341 
4342  else
4343  switch (dim)
4344  {
4345  case 1:
4346  {
4347  Assert(false,
4348  ExcMessage(
4349  "Computing the curl in 1d is not a useful operation"));
4350  return curl_type();
4351  }
4352 
4353  case 2:
4354  {
4355  if (snc != -1)
4356  {
4357  curl_type return_value;
4358 
4359  // the single nonzero component can only be zero or one in 2d
4360  if (shape_function_data[shape_function]
4361  .single_nonzero_component_index == 0)
4362  return_value[0] =
4363  -1.0 * fe_values->finite_element_output
4364  .shape_gradients[snc][q_point][1];
4365  else
4366  return_value[0] = fe_values->finite_element_output
4367  .shape_gradients[snc][q_point][0];
4368 
4369  return return_value;
4370  }
4371 
4372  else
4373  {
4374  curl_type return_value;
4375 
4376  return_value[0] = 0.0;
4377 
4378  if (shape_function_data[shape_function]
4379  .is_nonzero_shape_function_component[0])
4380  return_value[0] -=
4381  fe_values->finite_element_output
4382  .shape_gradients[shape_function_data[shape_function]
4383  .row_index[0]][q_point][1];
4384 
4385  if (shape_function_data[shape_function]
4386  .is_nonzero_shape_function_component[1])
4387  return_value[0] +=
4388  fe_values->finite_element_output
4389  .shape_gradients[shape_function_data[shape_function]
4390  .row_index[1]][q_point][0];
4391 
4392  return return_value;
4393  }
4394  }
4395 
4396  case 3:
4397  {
4398  if (snc != -1)
4399  {
4400  curl_type return_value;
4401 
4402  switch (shape_function_data[shape_function]
4403  .single_nonzero_component_index)
4404  {
4405  case 0:
4406  {
4407  return_value[0] = 0;
4408  return_value[1] = fe_values->finite_element_output
4409  .shape_gradients[snc][q_point][2];
4410  return_value[2] =
4411  -1.0 * fe_values->finite_element_output
4412  .shape_gradients[snc][q_point][1];
4413  return return_value;
4414  }
4415 
4416  case 1:
4417  {
4418  return_value[0] =
4419  -1.0 * fe_values->finite_element_output
4420  .shape_gradients[snc][q_point][2];
4421  return_value[1] = 0;
4422  return_value[2] = fe_values->finite_element_output
4423  .shape_gradients[snc][q_point][0];
4424  return return_value;
4425  }
4426 
4427  default:
4428  {
4429  return_value[0] = fe_values->finite_element_output
4430  .shape_gradients[snc][q_point][1];
4431  return_value[1] =
4432  -1.0 * fe_values->finite_element_output
4433  .shape_gradients[snc][q_point][0];
4434  return_value[2] = 0;
4435  return return_value;
4436  }
4437  }
4438  }
4439 
4440  else
4441  {
4442  curl_type return_value;
4443 
4444  for (unsigned int i = 0; i < dim; ++i)
4445  return_value[i] = 0.0;
4446 
4447  if (shape_function_data[shape_function]
4448  .is_nonzero_shape_function_component[0])
4449  {
4450  return_value[1] +=
4451  fe_values->finite_element_output
4452  .shape_gradients[shape_function_data[shape_function]
4453  .row_index[0]][q_point][2];
4454  return_value[2] -=
4455  fe_values->finite_element_output
4456  .shape_gradients[shape_function_data[shape_function]
4457  .row_index[0]][q_point][1];
4458  }
4459 
4460  if (shape_function_data[shape_function]
4461  .is_nonzero_shape_function_component[1])
4462  {
4463  return_value[0] -=
4464  fe_values->finite_element_output
4465  .shape_gradients[shape_function_data[shape_function]
4466  .row_index[1]][q_point][2];
4467  return_value[2] +=
4468  fe_values->finite_element_output
4469  .shape_gradients[shape_function_data[shape_function]
4470  .row_index[1]][q_point][0];
4471  }
4472 
4473  if (shape_function_data[shape_function]
4474  .is_nonzero_shape_function_component[2])
4475  {
4476  return_value[0] +=
4477  fe_values->finite_element_output
4478  .shape_gradients[shape_function_data[shape_function]
4479  .row_index[2]][q_point][1];
4480  return_value[1] -=
4481  fe_values->finite_element_output
4482  .shape_gradients[shape_function_data[shape_function]
4483  .row_index[2]][q_point][0];
4484  }
4485 
4486  return return_value;
4487  }
4488  }
4489  }
4490  // should not end up here
4491  Assert(false, ExcInternalError());
4492  return curl_type();
4493  }
4494 
4495 
4496 
4497  template <int dim, int spacedim>
4498  inline typename Vector<dim, spacedim>::hessian_type
4499  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
4500  const unsigned int q_point) const
4501  {
4502  // this function works like in the case above
4503  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4504  Assert(fe_values->update_flags & update_hessians,
4506  "update_hessians")));
4507 
4508  // same as for the scalar case except that we have one more index
4509  const int snc =
4510  shape_function_data[shape_function].single_nonzero_component;
4511  if (snc == -2)
4512  return hessian_type();
4513  else if (snc != -1)
4514  {
4515  hessian_type return_value;
4516  return_value[shape_function_data[shape_function]
4517  .single_nonzero_component_index] =
4518  fe_values->finite_element_output.shape_hessians[snc][q_point];
4519  return return_value;
4520  }
4521  else
4522  {
4523  hessian_type return_value;
4524  for (unsigned int d = 0; d < dim; ++d)
4525  if (shape_function_data[shape_function]
4526  .is_nonzero_shape_function_component[d])
4527  return_value[d] =
4528  fe_values->finite_element_output.shape_hessians
4529  [shape_function_data[shape_function].row_index[d]][q_point];
4530 
4531  return return_value;
4532  }
4533  }
4534 
4535 
4536 
4537  template <int dim, int spacedim>
4538  inline typename Vector<dim, spacedim>::third_derivative_type
4539  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
4540  const unsigned int q_point) const
4541  {
4542  // this function works like in the case above
4543  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4544  Assert(fe_values->update_flags & update_3rd_derivatives,
4546  "update_3rd_derivatives")));
4547 
4548  // same as for the scalar case except that we have one more index
4549  const int snc =
4550  shape_function_data[shape_function].single_nonzero_component;
4551  if (snc == -2)
4552  return third_derivative_type();
4553  else if (snc != -1)
4554  {
4555  third_derivative_type return_value;
4556  return_value[shape_function_data[shape_function]
4557  .single_nonzero_component_index] =
4558  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
4559  return return_value;
4560  }
4561  else
4562  {
4563  third_derivative_type return_value;
4564  for (unsigned int d = 0; d < dim; ++d)
4565  if (shape_function_data[shape_function]
4566  .is_nonzero_shape_function_component[d])
4567  return_value[d] =
4568  fe_values->finite_element_output.shape_3rd_derivatives
4569  [shape_function_data[shape_function].row_index[d]][q_point];
4570 
4571  return return_value;
4572  }
4573  }
4574 
4575 
4576 
4577  namespace internal
4578  {
4583  inline ::SymmetricTensor<2, 1>
4584  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
4585  {
4586  AssertIndexRange(n, 1);
4587  (void)n;
4588 
4589  return {{t[0]}};
4590  }
4591 
4592 
4593 
4594  inline ::SymmetricTensor<2, 2>
4595  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
4596  {
4597  switch (n)
4598  {
4599  case 0:
4600  {
4601  return {{t[0], 0, t[1] / 2}};
4602  }
4603  case 1:
4604  {
4605  return {{0, t[1], t[0] / 2}};
4606  }
4607  default:
4608  {
4609  AssertIndexRange(n, 2);
4610  return {};
4611  }
4612  }
4613  }
4614 
4615 
4616 
4617  inline ::SymmetricTensor<2, 3>
4618  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
4619  {
4620  switch (n)
4621  {
4622  case 0:
4623  {
4624  return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
4625  }
4626  case 1:
4627  {
4628  return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
4629  }
4630  case 2:
4631  {
4632  return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
4633  }
4634  default:
4635  {
4636  AssertIndexRange(n, 3);
4637  return {};
4638  }
4639  }
4640  }
4641  } // namespace internal
4642 
4643 
4644 
4645  template <int dim, int spacedim>
4646  inline typename Vector<dim, spacedim>::symmetric_gradient_type
4647  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
4648  const unsigned int q_point) const
4649  {
4650  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4651  Assert(fe_values->update_flags & update_gradients,
4653  "update_gradients")));
4654 
4655  // same as for the scalar case except that we have one more index
4656  const int snc =
4657  shape_function_data[shape_function].single_nonzero_component;
4658  if (snc == -2)
4659  return symmetric_gradient_type();
4660  else if (snc != -1)
4661  return internal::symmetrize_single_row(
4662  shape_function_data[shape_function].single_nonzero_component_index,
4663  fe_values->finite_element_output.shape_gradients[snc][q_point]);
4664  else
4665  {
4666  gradient_type return_value;
4667  for (unsigned int d = 0; d < dim; ++d)
4668  if (shape_function_data[shape_function]
4669  .is_nonzero_shape_function_component[d])
4670  return_value[d] =
4671  fe_values->finite_element_output.shape_gradients
4672  [shape_function_data[shape_function].row_index[d]][q_point];
4673 
4674  return symmetrize(return_value);
4675  }
4676  }
4677 
4678 
4679 
4680  template <int dim, int spacedim>
4682  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
4683  const unsigned int q_point) const
4684  {
4685  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4686  Assert(fe_values->update_flags & update_values,
4688  "update_values")));
4689 
4690  // similar to the vector case where we have more then one index and we need
4691  // to convert between unrolled and component indexing for tensors
4692  const int snc =
4693  shape_function_data[shape_function].single_nonzero_component;
4694 
4695  if (snc == -2)
4696  {
4697  // shape function is zero for the selected components
4698  return value_type();
4699  }
4700  else if (snc != -1)
4701  {
4702  value_type return_value;
4703  const unsigned int comp =
4704  shape_function_data[shape_function].single_nonzero_component_index;
4705  return_value[value_type::unrolled_to_component_indices(comp)] =
4706  fe_values->finite_element_output.shape_values(snc, q_point);
4707  return return_value;
4708  }
4709  else
4710  {
4711  value_type return_value;
4712  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
4713  if (shape_function_data[shape_function]
4714  .is_nonzero_shape_function_component[d])
4715  return_value[value_type::unrolled_to_component_indices(d)] =
4716  fe_values->finite_element_output.shape_values(
4717  shape_function_data[shape_function].row_index[d], q_point);
4718  return return_value;
4719  }
4720  }
4721 
4722 
4723 
4724  template <int dim, int spacedim>
4727  const unsigned int shape_function,
4728  const unsigned int q_point) const
4729  {
4730  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4731  Assert(fe_values->update_flags & update_gradients,
4733  "update_gradients")));
4734 
4735  const int snc =
4736  shape_function_data[shape_function].single_nonzero_component;
4737 
4738  if (snc == -2)
4739  {
4740  // shape function is zero for the selected components
4741  return divergence_type();
4742  }
4743  else if (snc != -1)
4744  {
4745  // we have a single non-zero component when the symmetric tensor is
4746  // represented in unrolled form. this implies we potentially have
4747  // two non-zero components when represented in component form! we
4748  // will only have one non-zero entry if the non-zero component lies on
4749  // the diagonal of the tensor.
4750  //
4751  // the divergence of a second-order tensor is a first order tensor.
4752  //
4753  // assume the second-order tensor is A with components A_{ij}. then
4754  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
4755  // entries in the tensorial representation. define the
4756  // divergence as:
4757  // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
4758  // (which is incidentally also
4759  // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
4760  // In both cases, a sum is implied.
4761  //
4762  // Now, we know the nonzero component in unrolled form: it is indicated
4763  // by 'snc'. we can figure out which tensor components belong to this:
4764  const unsigned int comp =
4765  shape_function_data[shape_function].single_nonzero_component_index;
4766  const unsigned int ii =
4767  value_type::unrolled_to_component_indices(comp)[0];
4768  const unsigned int jj =
4769  value_type::unrolled_to_component_indices(comp)[1];
4770 
4771  // given the form of the divergence above, if ii=jj there is only a
4772  // single nonzero component of the full tensor and the gradient
4773  // equals
4774  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
4775  // all other entries of 'b' are zero
4776  //
4777  // on the other hand, if ii!=jj, then there are two nonzero entries in
4778  // the full tensor and
4779  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
4780  // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
4781  // again, all other entries of 'b' are zero
4782  const ::Tensor<1, spacedim> &phi_grad =
4783  fe_values->finite_element_output.shape_gradients[snc][q_point];
4784 
4785  divergence_type return_value;
4786  return_value[ii] = phi_grad[jj];
4787 
4788  if (ii != jj)
4789  return_value[jj] = phi_grad[ii];
4790 
4791  return return_value;
4792  }
4793  else
4794  {
4795  Assert(false, ExcNotImplemented());
4796  divergence_type return_value;
4797  return return_value;
4798  }
4799  }
4800 
4801 
4802 
4803  template <int dim, int spacedim>
4804  inline typename Tensor<2, dim, spacedim>::value_type
4805  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
4806  const unsigned int q_point) const
4807  {
4808  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4809  Assert(fe_values->update_flags & update_values,
4811  "update_values")));
4812 
4813  // similar to the vector case where we have more then one index and we need
4814  // to convert between unrolled and component indexing for tensors
4815  const int snc =
4816  shape_function_data[shape_function].single_nonzero_component;
4817 
4818  if (snc == -2)
4819  {
4820  // shape function is zero for the selected components
4821  return value_type();
4822  }
4823  else if (snc != -1)
4824  {
4825  value_type return_value;
4826  const unsigned int comp =
4827  shape_function_data[shape_function].single_nonzero_component_index;
4828  const TableIndices<2> indices =
4830  return_value[indices] =
4831  fe_values->finite_element_output.shape_values(snc, q_point);
4832  return return_value;
4833  }
4834  else
4835  {
4836  value_type return_value;
4837  for (unsigned int d = 0; d < dim * dim; ++d)
4838  if (shape_function_data[shape_function]
4839  .is_nonzero_shape_function_component[d])
4840  {
4841  const TableIndices<2> indices =
4843  return_value[indices] =
4844  fe_values->finite_element_output.shape_values(
4845  shape_function_data[shape_function].row_index[d], q_point);
4846  }
4847  return return_value;
4848  }
4849  }
4850 
4851 
4852 
4853  template <int dim, int spacedim>
4855  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
4856  const unsigned int q_point) const
4857  {
4858  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4859  Assert(fe_values->update_flags & update_gradients,
4861  "update_gradients")));
4862 
4863  const int snc =
4864  shape_function_data[shape_function].single_nonzero_component;
4865 
4866  if (snc == -2)
4867  {
4868  // shape function is zero for the selected components
4869  return divergence_type();
4870  }
4871  else if (snc != -1)
4872  {
4873  // we have a single non-zero component when the tensor is
4874  // represented in unrolled form.
4875  //
4876  // the divergence of a second-order tensor is a first order tensor.
4877  //
4878  // assume the second-order tensor is A with components A_{ij},
4879  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
4880  //
4881  // Now, we know the nonzero component in unrolled form: it is indicated
4882  // by 'snc'. we can figure out which tensor components belong to this:
4883  const unsigned int comp =
4884  shape_function_data[shape_function].single_nonzero_component_index;
4885  const TableIndices<2> indices =
4887  const unsigned int ii = indices[0];
4888  const unsigned int jj = indices[1];
4889 
4890  const ::Tensor<1, spacedim> &phi_grad =
4891  fe_values->finite_element_output.shape_gradients[snc][q_point];
4892 
4893  divergence_type return_value;
4894  // note that we contract \nabla from the right
4895  return_value[ii] = phi_grad[jj];
4896 
4897  return return_value;
4898  }
4899  else
4900  {
4901  Assert(false, ExcNotImplemented());
4902  divergence_type return_value;
4903  return return_value;
4904  }
4905  }
4906 
4907 
4908 
4909  template <int dim, int spacedim>
4911  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
4912  const unsigned int q_point) const
4913  {
4914  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4915  Assert(fe_values->update_flags & update_gradients,
4917  "update_gradients")));
4918 
4919  const int snc =
4920  shape_function_data[shape_function].single_nonzero_component;
4921 
4922  if (snc == -2)
4923  {
4924  // shape function is zero for the selected components
4925  return gradient_type();
4926  }
4927  else if (snc != -1)
4928  {
4929  // we have a single non-zero component when the tensor is
4930  // represented in unrolled form.
4931  //
4932  // the gradient of a second-order tensor is a third order tensor.
4933  //
4934  // assume the second-order tensor is A with components A_{ij},
4935  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
4936  //
4937  // Now, we know the nonzero component in unrolled form: it is indicated
4938  // by 'snc'. we can figure out which tensor components belong to this:
4939  const unsigned int comp =
4940  shape_function_data[shape_function].single_nonzero_component_index;
4941  const TableIndices<2> indices =
4943  const unsigned int ii = indices[0];
4944  const unsigned int jj = indices[1];
4945 
4946  const ::Tensor<1, spacedim> &phi_grad =
4947  fe_values->finite_element_output.shape_gradients[snc][q_point];
4948 
4949  gradient_type return_value;
4950  return_value[ii][jj] = phi_grad;
4951 
4952  return return_value;
4953  }
4954  else
4955  {
4956  Assert(false, ExcNotImplemented());
4957  gradient_type return_value;
4958  return return_value;
4959  }
4960  }
4961 
4962 } // namespace FEValuesViews
4963 
4964 
4965 
4966 /*---------------------- Inline functions: FEValuesBase ---------------------*/
4967 
4968 
4969 
4970 template <int dim, int spacedim>
4972  operator[](const FEValuesExtractors::Scalar &scalar) const
4973 {
4974  AssertIndexRange(scalar.component, fe_values_views_cache.scalars.size());
4975 
4976  return fe_values_views_cache.scalars[scalar.component];
4977 }
4978 
4979 
4980 
4981 template <int dim, int spacedim>
4983  operator[](const FEValuesExtractors::Vector &vector) const
4984 {
4986  fe_values_views_cache.vectors.size());
4987 
4988  return fe_values_views_cache.vectors[vector.first_vector_component];
4989 }
4990 
4991 
4992 
4993 template <int dim, int spacedim>
4997 {
4998  Assert(
4999  tensor.first_tensor_component <
5000  fe_values_views_cache.symmetric_second_order_tensors.size(),
5002  0,
5003  fe_values_views_cache.symmetric_second_order_tensors.size()));
5004 
5005  return fe_values_views_cache
5006  .symmetric_second_order_tensors[tensor.first_tensor_component];
5007 }
5008 
5009 
5010 
5011 template <int dim, int spacedim>
5014  operator[](const FEValuesExtractors::Tensor<2> &tensor) const
5015 {
5017  fe_values_views_cache.second_order_tensors.size());
5018 
5019  return fe_values_views_cache
5020  .second_order_tensors[tensor.first_tensor_component];
5021 }
5022 
5023 
5024 
5025 template <int dim, int spacedim>
5026 inline const double &
5027 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
5028  const unsigned int j) const
5029 {
5030  AssertIndexRange(i, fe->n_dofs_per_cell());
5031  Assert(this->update_flags & update_values,
5032  ExcAccessToUninitializedField("update_values"));
5033  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5034  Assert(present_cell.get() != nullptr,
5035  ExcMessage("FEValues object is not reinit'ed to any cell"));
5036  // if the entire FE is primitive,
5037  // then we can take a short-cut:
5038  if (fe->is_primitive())
5039  return this->finite_element_output.shape_values(i, j);
5040  else
5041  {
5042  // otherwise, use the mapping
5043  // between shape function
5044  // numbers and rows. note that
5045  // by the assertions above, we
5046  // know that this particular
5047  // shape function is primitive,
5048  // so we can call
5049  // system_to_component_index
5050  const unsigned int row =
5051  this->finite_element_output
5052  .shape_function_to_row_table[i * fe->n_components() +
5053  fe->system_to_component_index(i).first];
5054  return this->finite_element_output.shape_values(row, j);
5055  }
5056 }
5057 
5058 
5059 
5060 template <int dim, int spacedim>
5061 inline double
5063  const unsigned int i,
5064  const unsigned int j,
5065  const unsigned int component) const
5066 {
5067  AssertIndexRange(i, fe->n_dofs_per_cell());
5068  Assert(this->update_flags & update_values,
5069  ExcAccessToUninitializedField("update_values"));
5070  AssertIndexRange(component, fe->n_components());
5071  Assert(present_cell.get() != nullptr,
5072  ExcMessage("FEValues object is not reinit'ed to any cell"));
5073 
5074  // check whether the shape function
5075  // is non-zero at all within
5076  // this component:
5077  if (fe->get_nonzero_components(i)[component] == false)
5078  return 0;
5079 
5080  // look up the right row in the
5081  // table and take the data from
5082  // there
5083  const unsigned int row =
5084  this->finite_element_output
5085  .shape_function_to_row_table[i * fe->n_components() + component];
5086  return this->finite_element_output.shape_values(row, j);
5087 }
5088 
5089 
5090 
5091 template <int dim, int spacedim>
5092 inline const Tensor<1, spacedim> &
5093 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
5094  const unsigned int j) const
5095 {
5096  AssertIndexRange(i, fe->n_dofs_per_cell());
5097  Assert(this->update_flags & update_gradients,
5098  ExcAccessToUninitializedField("update_gradients"));
5099  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5100  Assert(present_cell.get() != nullptr,
5101  ExcMessage("FEValues object is not reinit'ed to any cell"));
5102  // if the entire FE is primitive,
5103  // then we can take a short-cut:
5104  if (fe->is_primitive())
5105  return this->finite_element_output.shape_gradients[i][j];
5106  else
5107  {
5108  // otherwise, use the mapping
5109  // between shape function
5110  // numbers and rows. note that
5111  // by the assertions above, we
5112  // know that this particular
5113  // shape function is primitive,
5114  // so we can call
5115  // system_to_component_index
5116  const unsigned int row =
5117  this->finite_element_output
5118  .shape_function_to_row_table[i * fe->n_components() +
5119  fe->system_to_component_index(i).first];
5120  return this->finite_element_output.shape_gradients[row][j];
5121  }
5122 }
5123 
5124 
5125 
5126 template <int dim, int spacedim>
5127 inline Tensor<1, spacedim>
5129  const unsigned int i,
5130  const unsigned int j,
5131  const unsigned int component) const
5132 {
5133  AssertIndexRange(i, fe->n_dofs_per_cell());
5134  Assert(this->update_flags & update_gradients,
5135  ExcAccessToUninitializedField("update_gradients"));
5136  AssertIndexRange(component, fe->n_components());
5137  Assert(present_cell.get() != nullptr,
5138  ExcMessage("FEValues object is not reinit'ed to any cell"));
5139  // check whether the shape function
5140  // is non-zero at all within
5141  // this component:
5142  if (fe->get_nonzero_components(i)[component] == false)
5143  return Tensor<1, spacedim>();
5144 
5145  // look up the right row in the
5146  // table and take the data from
5147  // there
5148  const unsigned int row =
5149  this->finite_element_output
5150  .shape_function_to_row_table[i * fe->n_components() + component];
5151  return this->finite_element_output.shape_gradients[row][j];
5152 }
5153 
5154 
5155 
5156 template <int dim, int spacedim>
5157 inline const Tensor<2, spacedim> &
5158 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
5159  const unsigned int j) const
5160 {
5161  AssertIndexRange(i, fe->n_dofs_per_cell());
5162  Assert(this->update_flags & update_hessians,
5163  ExcAccessToUninitializedField("update_hessians"));
5164  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5165  Assert(present_cell.get() != nullptr,
5166  ExcMessage("FEValues object is not reinit'ed to any cell"));
5167  // if the entire FE is primitive,
5168  // then we can take a short-cut:
5169  if (fe->is_primitive())
5170  return this->finite_element_output.shape_hessians[i][j];
5171  else
5172  {
5173  // otherwise, use the mapping
5174  // between shape function
5175  // numbers and rows. note that
5176  // by the assertions above, we
5177  // know that this particular
5178  // shape function is primitive,
5179  // so we can call
5180  // system_to_component_index
5181  const unsigned int row =
5182  this->finite_element_output
5183  .shape_function_to_row_table[i * fe->n_components() +
5184  fe->system_to_component_index(i).first];
5185  return this->finite_element_output.shape_hessians[row][j];
5186  }
5187 }
5188 
5189 
5190 
5191 template <int dim, int spacedim>
5192 inline Tensor<2, spacedim>
5194  const unsigned int i,
5195  const unsigned int j,
5196  const unsigned int component) const
5197 {
5198  AssertIndexRange(i, fe->n_dofs_per_cell());
5199  Assert(this->update_flags & update_hessians,
5200  ExcAccessToUninitializedField("update_hessians"));
5201  AssertIndexRange(component, fe->n_components());
5202  Assert(present_cell.get() != nullptr,
5203  ExcMessage("FEValues object is not reinit'ed to any cell"));
5204  // check whether the shape function
5205  // is non-zero at all within
5206  // this component:
5207  if (fe->get_nonzero_components(i)[component] == false)
5208  return Tensor<2, spacedim>();
5209 
5210  // look up the right row in the
5211  // table and take the data from
5212  // there
5213  const unsigned int row =
5214  this->finite_element_output
5215  .shape_function_to_row_table[i * fe->n_components() + component];
5216  return this->finite_element_output.shape_hessians[row][j];
5217 }
5218 
5219 
5220 
5221 template <int dim, int spacedim>
5222 inline const Tensor<3, spacedim> &
5224  const unsigned int j) const
5225 {
5226  AssertIndexRange(i, fe->n_dofs_per_cell());
5227  Assert(this->update_flags & update_3rd_derivatives,
5228  ExcAccessToUninitializedField("update_3rd_derivatives"));
5229  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5230  Assert(present_cell.get() != nullptr,
5231  ExcMessage("FEValues object is not reinit'ed to any cell"));
5232  // if the entire FE is primitive,
5233  // then we can take a short-cut:
5234  if (fe->is_primitive())
5235  return this->finite_element_output.shape_3rd_derivatives[i][j];
5236  else
5237  {
5238  // otherwise, use the mapping
5239  // between shape function
5240  // numbers and rows. note that
5241  // by the assertions above, we
5242  // know that this particular
5243  // shape function is primitive,
5244  // so we can call
5245  // system_to_component_index
5246  const unsigned int row =
5247  this->finite_element_output
5248  .shape_function_to_row_table[i * fe->n_components() +
5249  fe->system_to_component_index(i).first];
5250  return this->finite_element_output.shape_3rd_derivatives[row][j];
5251  }
5252 }
5253 
5254 
5255 
5256 template <int dim, int spacedim>
5257 inline Tensor<3, spacedim>
5259  const unsigned int i,
5260  const unsigned int j,
5261  const unsigned int component) const
5262 {
5263  AssertIndexRange(i, fe->n_dofs_per_cell());
5264  Assert(this->update_flags & update_3rd_derivatives,
5265  ExcAccessToUninitializedField("update_3rd_derivatives"));
5266  AssertIndexRange(component, fe->n_components());
5267  Assert(present_cell.get() != nullptr,
5268  ExcMessage("FEValues object is not reinit'ed to any cell"));
5269  // check whether the shape function
5270  // is non-zero at all within
5271  // this component:
5272  if (fe->get_nonzero_components(i)[component] == false)
5273  return Tensor<3, spacedim>();
5274 
5275  // look up the right row in the
5276  // table and take the data from
5277  // there
5278  const unsigned int row =
5279  this->finite_element_output
5280  .shape_function_to_row_table[i * fe->n_components() + component];
5281  return this->finite_element_output.shape_3rd_derivatives[row][j];
5282 }
5283 
5284 
5285 
5286 template <int dim, int spacedim>
5287 inline const FiniteElement<dim, spacedim> &
5289 {
5290  return *fe;
5291 }
5292 
5293 
5294 
5295 template <int dim, int spacedim>
5296 inline const Mapping<dim, spacedim> &
5298 {
5299  return *mapping;
5300 }
5301 
5302 
5303 
5304 template <int dim, int spacedim>
5305 inline UpdateFlags
5307 {
5308  return this->update_flags;
5309 }
5310 
5311 
5312 
5313 template <int dim, int spacedim>
5314 inline const std::vector<Point<spacedim>> &
5316 {
5317  Assert(this->update_flags & update_quadrature_points,
5318  ExcAccessToUninitializedField("update_quadrature_points"));
5319  Assert(present_cell.get() != nullptr,
5320  ExcMessage("FEValues object is not reinit'ed to any cell"));
5321  return this->mapping_output.quadrature_points;
5322 }
5323 
5324 
5325 
5326 template <int dim, int spacedim>
5327 inline const std::vector<double> &
5329 {
5330  Assert(this->update_flags & update_JxW_values,
5331  ExcAccessToUninitializedField("update_JxW_values"));
5332  Assert(present_cell.get() != nullptr,
5333  ExcMessage("FEValues object is not reinit'ed to any cell"));
5334  return this->mapping_output.JxW_values;
5335 }
5336 
5337 
5338 
5339 template <int dim, int spacedim>
5340 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5342 {
5343  Assert(this->update_flags & update_jacobians,
5344  ExcAccessToUninitializedField("update_jacobians"));
5345  Assert(present_cell.get() != nullptr,
5346  ExcMessage("FEValues object is not reinit'ed to any cell"));
5347  return this->mapping_output.jacobians;
5348 }
5349 
5350 
5351 
5352 template <int dim, int spacedim>
5353 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5355 {
5356  Assert(this->update_flags & update_jacobian_grads,
5357  ExcAccessToUninitializedField("update_jacobians_grads"));
5358  Assert(present_cell.get() != nullptr,
5359  ExcMessage("FEValues object is not reinit'ed to any cell"));
5360  return this->mapping_output.jacobian_grads;
5361 }
5362 
5363 
5364 
5365 template <int dim, int spacedim>
5366 inline const Tensor<3, spacedim> &
5368  const unsigned int i) const
5369 {
5370  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5371  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5372  Assert(present_cell.get() != nullptr,
5373  ExcMessage("FEValues object is not reinit'ed to any cell"));
5374  return this->mapping_output.jacobian_pushed_forward_grads[i];
5375 }
5376 
5377 
5378 
5379 template <int dim, int spacedim>
5380 inline const std::vector<Tensor<3, spacedim>> &
5382 {
5383  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5384  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5385  Assert(present_cell.get() != nullptr,
5386  ExcMessage("FEValues object is not reinit'ed to any cell"));
5387  return this->mapping_output.jacobian_pushed_forward_grads;
5388 }
5389 
5390 
5391 
5392 template <int dim, int spacedim>
5393 inline const DerivativeForm<3, dim, spacedim> &
5394 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5395 {
5396  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5397  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5398  Assert(present_cell.get() != nullptr,
5399  ExcMessage("FEValues object is not reinit'ed to any cell"));
5400  return this->mapping_output.jacobian_2nd_derivatives[i];
5401 }
5402 
5403 
5404 
5405 template <int dim, int spacedim>
5406 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5408 {
5409  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5410  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5411  Assert(present_cell.get() != nullptr,
5412  ExcMessage("FEValues object is not reinit'ed to any cell"));
5413  return this->mapping_output.jacobian_2nd_derivatives;
5414 }
5415 
5416 
5417 
5418 template <int dim, int spacedim>
5419 inline const Tensor<4, spacedim> &
5421  const unsigned int i) const
5422 {
5425  "update_jacobian_pushed_forward_2nd_derivatives"));
5426  Assert(present_cell.get() != nullptr,
5427  ExcMessage("FEValues object is not reinit'ed to any cell"));
5428  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5429 }
5430 
5431 
5432 
5433 template <int dim, int spacedim>
5434 inline const std::vector<Tensor<4, spacedim>> &
5436 {
5437  Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
5439  "update_jacobian_pushed_forward_2nd_derivatives"));
5440  Assert(present_cell.get() != nullptr,
5441  ExcMessage("FEValues object is not reinit'ed to any cell"));
5442  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5443 }
5444 
5445 
5446 
5447 template <int dim, int spacedim>
5448 inline const DerivativeForm<4, dim, spacedim> &
5449 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5450 {
5451  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5452  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5453  Assert(present_cell.get() != nullptr,
5454  ExcMessage("FEValues object is not reinit'ed to any cell"));
5455  return this->mapping_output.jacobian_3rd_derivatives[i];
5456 }
5457 
5458 
5459 
5460 template <int dim, int spacedim>
5461 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5463 {
5464  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5465  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5466  Assert(present_cell.get() != nullptr,
5467  ExcMessage("FEValues object is not reinit'ed to any cell"));
5468  return this->mapping_output.jacobian_3rd_derivatives;
5469 }
5470 
5471 
5472 
5473 template <int dim, int spacedim>
5474 inline const Tensor<5, spacedim> &
5476  const unsigned int i) const
5477 {
5480  "update_jacobian_pushed_forward_3rd_derivatives"));
5481  Assert(present_cell.get() != nullptr,
5482  ExcMessage("FEValues object is not reinit'ed to any cell"));
5483  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
5484 }
5485 
5486 
5487 
5488 template <int dim, int spacedim>
5489 inline const std::vector<Tensor<5, spacedim>> &
5491 {
5492  Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
5494  "update_jacobian_pushed_forward_3rd_derivatives"));
5495  Assert(present_cell.get() != nullptr,
5496  ExcMessage("FEValues object is not reinit'ed to any cell"));
5497  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
5498 }
5499 
5500 
5501 
5502 template <int dim, int spacedim>
5503 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
5505 {
5506  Assert(this->update_flags & update_inverse_jacobians,
5507  ExcAccessToUninitializedField("update_inverse_jacobians"));
5508  Assert(present_cell.get() != nullptr,
5509  ExcMessage("FEValues object is not reinit'ed to any cell"));
5510  return this->mapping_output.inverse_jacobians;
5511 }
5512 
5513 
5514 
5515 template <int dim, int spacedim>
5518 {
5519  return {0U, dofs_per_cell};
5520 }
5521 
5522 
5523 
5524 template <int dim, int spacedim>
5527  const unsigned int start_dof_index) const
5528 {
5529  Assert(start_dof_index <= dofs_per_cell,
5530  ExcIndexRange(start_dof_index, 0, dofs_per_cell + 1));
5531  return {start_dof_index, dofs_per_cell};
5532 }
5533 
5534 
5535 
5536 template <int dim, int spacedim>
5539  const unsigned int end_dof_index) const
5540 {
5541  Assert(end_dof_index < dofs_per_cell,
5542  ExcIndexRange(end_dof_index, 0, dofs_per_cell));
5543  return {0U, end_dof_index + 1};
5544 }
5545 
5546 
5547 
5548 template <int dim, int spacedim>
5551 {
5552  return {0U, n_quadrature_points};
5553 }
5554 
5555 
5556 
5557 template <int dim, int spacedim>
5558 inline const Point<spacedim> &
5559 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
5560 {
5561  Assert(this->update_flags & update_quadrature_points,
5562  ExcAccessToUninitializedField("update_quadrature_points"));
5563  AssertIndexRange(i, this->mapping_output.quadrature_points.size());
5564  Assert(present_cell.get() != nullptr,
5565  ExcMessage("FEValues object is not reinit'ed to any cell"));
5566 
5567  return this->mapping_output.quadrature_points[i];
5568 }
5569 
5570 
5571 
5572 template <int dim, int spacedim>
5573 inline double
5574 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
5575 {
5576  Assert(this->update_flags & update_JxW_values,
5577  ExcAccessToUninitializedField("update_JxW_values"));
5578  AssertIndexRange(i, this->mapping_output.JxW_values.size());
5579  Assert(present_cell.get() != nullptr,
5580  ExcMessage("FEValues object is not reinit'ed to any cell"));
5581 
5582  return this->mapping_output.JxW_values[i];
5583 }
5584 
5585 
5586 
5587 template <int dim, int spacedim>
5588 inline const DerivativeForm<1, dim, spacedim> &
5589 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
5590 {
5591  Assert(this->update_flags & update_jacobians,
5592  ExcAccessToUninitializedField("update_jacobians"));
5593  AssertIndexRange(i, this->mapping_output.jacobians.size());
5594  Assert(present_cell.get() != nullptr,
5595  ExcMessage("FEValues object is not reinit'ed to any cell"));
5596 
5597  return this->mapping_output.jacobians[i];
5598 }
5599 
5600 
5601 
5602 template <int dim, int spacedim>
5603 inline const DerivativeForm<2, dim, spacedim> &
5604 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
5605 {
5606  Assert(this->update_flags & update_jacobian_grads,
5607  ExcAccessToUninitializedField("update_jacobians_grads"));
5608  AssertIndexRange(i, this->mapping_output.jacobian_grads.size());
5609  Assert(present_cell.get() != nullptr,
5610  ExcMessage("FEValues object is not reinit'ed to any cell"));
5611 
5612  return this->mapping_output.jacobian_grads[i];
5613 }
5614 
5615 
5616 
5617 template <int dim, int spacedim>
5618 inline const DerivativeForm<1, spacedim, dim> &
5619 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
5620 {
5621  Assert(this->update_flags & update_inverse_jacobians,
5622  ExcAccessToUninitializedField("update_inverse_jacobians"));
5623  AssertIndexRange(i, this->mapping_output.inverse_jacobians.size());
5624  Assert(present_cell.get() != nullptr,
5625  ExcMessage("FEValues object is not reinit'ed to any cell"));
5626 
5627  return this->mapping_output.inverse_jacobians[i];
5628 }
5629 
5630 
5631 
5632 template <int dim, int spacedim>
5633 inline const Tensor<1, spacedim> &
5634 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
5635 {
5636  Assert(this->update_flags & update_normal_vectors,
5638  "update_normal_vectors")));
5639  AssertIndexRange(i, this->mapping_output.normal_vectors.size());
5640  Assert(present_cell.get() != nullptr,
5641  ExcMessage("FEValues object is not reinit'ed to any cell"));
5642 
5643  return this->mapping_output.normal_vectors[i];
5644 }
5645 
5646 
5647 
5648 /*--------------------- Inline functions: FEValues --------------------------*/
5649 
5650 
5651 template <int dim, int spacedim>
5652 inline const Quadrature<dim> &
5654 {
5655  return quadrature;
5656 }
5657 
5658 
5659 
5660 template <int dim, int spacedim>
5661 inline const FEValues<dim, spacedim> &
5663 {
5664  return *this;
5665 }
5666 
5667 
5668 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
5669 
5670 
5671 template <int dim, int spacedim>
5672 inline unsigned int
5674 {
5675  return present_face_index;
5676 }
5677 
5678 
5679 /*----------------------- Inline functions: FE*FaceValues -------------------*/
5680 
5681 template <int dim, int spacedim>
5682 inline const Quadrature<dim - 1> &
5684 {
5685  return quadrature;
5686 }
5687 
5688 
5689 
5690 template <int dim, int spacedim>
5691 inline const FEFaceValues<dim, spacedim> &
5693 {
5694  return *this;
5695 }
5696 
5697 
5698 
5699 template <int dim, int spacedim>
5700 inline const FESubfaceValues<dim, spacedim> &
5702 {
5703  return *this;
5704 }
5705 
5706 
5707 
5708 template <int dim, int spacedim>
5709 inline const Tensor<1, spacedim> &
5710 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
5711 {
5712  AssertIndexRange(i, this->mapping_output.boundary_forms.size());
5713  Assert(this->update_flags & update_boundary_forms,
5715  "update_boundary_forms")));
5716 
5717  return this->mapping_output.boundary_forms[i];
5718 }
5719 
5720 #endif // DOXYGEN
5721 
5723 
5724 #endif
Transformed quadrature weights.
Shape function values.
typename ProductType< Number, typename Vector< dim, spacedim >::curl_type >::type curl_type
Definition: fe_values.h:695
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:3529
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:3561
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:1592
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1291
unsigned int present_face_index
Definition: fe_values.h:3792
std::vector<::FEValuesViews::Vector< dim, spacedim > > vectors
Definition: fe_values.h:1948
typename ::internal::FEValuesViews::ViewType< dim, spacedim, Extractor >::type View
Definition: fe_values.h:1970
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1517
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:535
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:2097
typename ::internal::CurlType< spacedim >::type curl_type
Definition: fe_values.h:626
const unsigned int component
Definition: fe_values.h:541
const Quadrature< dim - 1 > quadrature
Definition: fe_values.h:3797
Volume element.
typename ::FEValuesViews::SymmetricTensor< rank, dim, spacedim > type
Definition: fe_values.h:1930
const ReferenceCell::internal::Info::Base & get_cell(const ReferenceCell::Type &type)
#define AssertIndexRange(index, range)
Definition: exceptions.h:1648
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:212
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:3445
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:188
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:3497
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:3576
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:1950
typename ProductType< Number, typename Vector< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:663
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:196
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1576
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:711
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:1862
unsigned int get_face_index() const
#define Assert(cond, exc)
Definition: exceptions.h:1423
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:301
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1226
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
const Quadrature< dim > quadrature
Definition: fe_values.h:3691
const unsigned int first_vector_component
Definition: fe_values.h:1221
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:369
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition: fe_values.h:3505
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:482
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:1299
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:1947
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:2090
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:3470
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:655
typename ::FEValuesViews::Tensor< rank, dim, spacedim > type
Definition: fe_values.h:1916
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
Definition: tensor.h:448
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:671
typename ProductType< Number, typename Vector< dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:679
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:368
typename ProductType< Number, typename Vector< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:703
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:1584
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:1851
Definition: fe.h:38
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1506
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:1215
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:3461
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3537
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:687
::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:3512
const std::vector< Tensor< 5, spacedim > > & get_jacobian_pushed_forward_3rd_derivatives() const
const FESubfaceValues< dim, spacedim > & get_present_fe_values() const
std::vector<::FEValuesViews::Tensor< 2, dim, spacedim > > second_order_tensors
Definition: fe_values.h:1952
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:204
typename ProductType< Number, typename Scalar< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:220
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:546
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
UpdateFlags update_flags
Definition: fe_values.h:3543
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:3521
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