Reference documentation for deal.II version Git 0b5e0ec63c 2020-01-23 00:52:06 -0500
\(\newcommand{\vcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\vcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
fe_values.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2019 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 
22 #include <deal.II/base/derivative_form.h>
23 #include <deal.II/base/exceptions.h>
24 #include <deal.II/base/point.h>
25 #include <deal.II/base/quadrature.h>
26 #include <deal.II/base/subscriptor.h>
28 
29 #include <deal.II/dofs/dof_accessor.h>
30 #include <deal.II/dofs/dof_handler.h>
31 
32 #include <deal.II/fe/fe.h>
33 #include <deal.II/fe/fe_update_flags.h>
34 #include <deal.II/fe/fe_values_extractors.h>
35 #include <deal.II/fe/mapping.h>
36 
37 #include <deal.II/grid/tria.h>
38 #include <deal.II/grid/tria_iterator.h>
39 
40 #include <deal.II/hp/dof_handler.h>
41 
42 #include <algorithm>
43 #include <memory>
44 #include <type_traits>
45 
46 
47 // dummy include in order to have the
48 // definition of PetscScalar available
49 // without including other PETSc stuff
50 #ifdef DEAL_II_WITH_PETSC
51 # include <petsc.h>
52 #endif
53 
54 DEAL_II_NAMESPACE_OPEN
55 
56 // Forward declaration
57 #ifndef DOXYGEN
58 template <int dim, int spacedim = dim>
59 class FEValuesBase;
60 #endif
61 
62 namespace internal
63 {
68  template <int dim, class NumberType = double>
69  struct CurlType;
70 
77  template <class NumberType>
78  struct CurlType<1, NumberType>
79  {
81  };
82 
89  template <class NumberType>
90  struct CurlType<2, NumberType>
91  {
93  };
94 
101  template <class NumberType>
102  struct CurlType<3, NumberType>
103  {
105  };
106 } // namespace internal
107 
108 
109 
131 namespace FEValuesViews
132 {
144  template <int dim, int spacedim = dim>
145  class Scalar
146  {
147  public:
153  using value_type = double;
154 
161 
168 
175 
180  template <typename Number>
181  struct OutputType
182  {
187  using value_type =
188  typename ProductType<Number,
190 
195  using gradient_type = typename ProductType<
196  Number,
198 
203  using laplacian_type =
204  typename ProductType<Number,
206 
211  using hessian_type = typename ProductType<
212  Number,
214 
219  using third_derivative_type = typename ProductType<
220  Number,
222  };
223 
229  {
239 
248  unsigned int row_index;
249  };
250 
254  Scalar();
255 
261  Scalar(const FEValuesBase<dim, spacedim> &fe_values_base,
262  const unsigned int component);
263 
268  Scalar &
269  operator=(const Scalar<dim, spacedim> &) = delete;
270 
284  value_type
285  value(const unsigned int shape_function, const unsigned int q_point) const;
286 
298  gradient(const unsigned int shape_function,
299  const unsigned int q_point) const;
300 
312  hessian(const unsigned int shape_function,
313  const unsigned int q_point) const;
314 
326  third_derivative(const unsigned int shape_function,
327  const unsigned int q_point) const;
328 
346  template <class InputVector>
347  void
348  get_function_values(
349  const InputVector &fe_function,
350  std::vector<typename ProductType<value_type,
351  typename InputVector::value_type>::type>
352  &values) const;
353 
376  template <class InputVector>
377  void
378  get_function_values_from_local_dof_values(
379  const InputVector &dof_values,
380  std::vector<
382  &values) const;
383 
401  template <class InputVector>
402  void
403  get_function_gradients(
404  const InputVector &fe_function,
405  std::vector<typename ProductType<gradient_type,
406  typename InputVector::value_type>::type>
407  &gradients) const;
408 
412  template <class InputVector>
413  void
414  get_function_gradients_from_local_dof_values(
415  const InputVector &dof_values,
416  std::vector<
418  &gradients) const;
419 
437  template <class InputVector>
438  void
439  get_function_hessians(
440  const InputVector &fe_function,
441  std::vector<typename ProductType<hessian_type,
442  typename InputVector::value_type>::type>
443  &hessians) const;
444 
448  template <class InputVector>
449  void
450  get_function_hessians_from_local_dof_values(
451  const InputVector &dof_values,
452  std::vector<
454  &hessians) const;
455 
456 
475  template <class InputVector>
476  void
477  get_function_laplacians(
478  const InputVector &fe_function,
479  std::vector<typename ProductType<value_type,
480  typename InputVector::value_type>::type>
481  &laplacians) const;
482 
486  template <class InputVector>
487  void
488  get_function_laplacians_from_local_dof_values(
489  const InputVector &dof_values,
490  std::vector<
492  &laplacians) const;
493 
494 
513  template <class InputVector>
514  void
515  get_function_third_derivatives(
516  const InputVector &fe_function,
517  std::vector<typename ProductType<third_derivative_type,
518  typename InputVector::value_type>::type>
519  &third_derivatives) const;
520 
524  template <class InputVector>
525  void
526  get_function_third_derivatives_from_local_dof_values(
527  const InputVector & dof_values,
528  std::vector<typename OutputType<typename InputVector::value_type>::
529  third_derivative_type> &third_derivatives) const;
530 
531 
532  private:
537 
542  const unsigned int component;
543 
547  std::vector<ShapeFunctionData> shape_function_data;
548  };
549 
550 
551 
581  template <int dim, int spacedim = dim>
582  class Vector
583  {
584  public:
591 
601 
613 
619  using divergence_type = double;
620 
627  using curl_type = typename ::internal::CurlType<spacedim>::type;
628 
635 
642 
647  template <typename Number>
648  struct OutputType
649  {
654  using value_type =
655  typename ProductType<Number,
657 
662  using gradient_type = typename ProductType<
663  Number,
665 
670  using symmetric_gradient_type = typename ProductType<
671  Number,
673 
678  using divergence_type = typename ProductType<
679  Number,
681 
686  using laplacian_type =
687  typename ProductType<Number,
689 
694  using curl_type =
695  typename ProductType<Number,
697 
702  using hessian_type = typename ProductType<
703  Number,
705 
710  using third_derivative_type = typename ProductType<
711  Number,
713  };
714 
720  {
729  bool is_nonzero_shape_function_component[spacedim];
730 
740  unsigned int row_index[spacedim];
741 
751  unsigned int single_nonzero_component_index;
752  };
753 
757  Vector();
758 
767  Vector(const FEValuesBase<dim, spacedim> &fe_values_base,
768  const unsigned int first_vector_component);
769 
774  Vector &
775  operator=(const Vector<dim, spacedim> &) = delete;
776 
793  value_type
794  value(const unsigned int shape_function, const unsigned int q_point) const;
795 
810  gradient(const unsigned int shape_function,
811  const unsigned int q_point) const;
812 
829  symmetric_gradient(const unsigned int shape_function,
830  const unsigned int q_point) const;
831 
843  divergence(const unsigned int shape_function,
844  const unsigned int q_point) const;
845 
866  curl_type
867  curl(const unsigned int shape_function, const unsigned int q_point) const;
868 
880  hessian(const unsigned int shape_function,
881  const unsigned int q_point) const;
882 
894  third_derivative(const unsigned int shape_function,
895  const unsigned int q_point) const;
896 
914  template <class InputVector>
915  void
916  get_function_values(
917  const InputVector &fe_function,
918  std::vector<typename ProductType<value_type,
919  typename InputVector::value_type>::type>
920  &values) const;
921 
944  template <class InputVector>
945  void
946  get_function_values_from_local_dof_values(
947  const InputVector &dof_values,
948  std::vector<
950  &values) const;
951 
969  template <class InputVector>
970  void
971  get_function_gradients(
972  const InputVector &fe_function,
973  std::vector<typename ProductType<gradient_type,
974  typename InputVector::value_type>::type>
975  &gradients) const;
976 
980  template <class InputVector>
981  void
982  get_function_gradients_from_local_dof_values(
983  const InputVector &dof_values,
984  std::vector<
986  &gradients) const;
987 
1011  template <class InputVector>
1012  void
1013  get_function_symmetric_gradients(
1014  const InputVector &fe_function,
1015  std::vector<typename ProductType<symmetric_gradient_type,
1016  typename InputVector::value_type>::type>
1017  &symmetric_gradients) const;
1018 
1022  template <class InputVector>
1023  void
1024  get_function_symmetric_gradients_from_local_dof_values(
1025  const InputVector & dof_values,
1026  std::vector<typename OutputType<typename InputVector::value_type>::
1027  symmetric_gradient_type> &symmetric_gradients) const;
1028 
1047  template <class InputVector>
1048  void
1049  get_function_divergences(
1050  const InputVector &fe_function,
1051  std::vector<typename ProductType<divergence_type,
1052  typename InputVector::value_type>::type>
1053  &divergences) const;
1054 
1058  template <class InputVector>
1059  void
1060  get_function_divergences_from_local_dof_values(
1061  const InputVector &dof_values,
1062  std::vector<
1064  &divergences) const;
1065 
1084  template <class InputVector>
1085  void
1086  get_function_curls(
1087  const InputVector &fe_function,
1088  std::vector<
1089  typename ProductType<curl_type, typename InputVector::value_type>::type>
1090  &curls) const;
1091 
1095  template <class InputVector>
1096  void
1097  get_function_curls_from_local_dof_values(
1098  const InputVector &dof_values,
1099  std::vector<
1101  &curls) const;
1102 
1120  template <class InputVector>
1121  void
1122  get_function_hessians(
1123  const InputVector &fe_function,
1124  std::vector<typename ProductType<hessian_type,
1125  typename InputVector::value_type>::type>
1126  &hessians) const;
1127 
1131  template <class InputVector>
1132  void
1133  get_function_hessians_from_local_dof_values(
1134  const InputVector &dof_values,
1135  std::vector<
1137  &hessians) const;
1138 
1157  template <class InputVector>
1158  void
1159  get_function_laplacians(
1160  const InputVector &fe_function,
1161  std::vector<typename ProductType<value_type,
1162  typename InputVector::value_type>::type>
1163  &laplacians) const;
1164 
1168  template <class InputVector>
1169  void
1170  get_function_laplacians_from_local_dof_values(
1171  const InputVector &dof_values,
1172  std::vector<
1174  &laplacians) const;
1175 
1194  template <class InputVector>
1195  void
1196  get_function_third_derivatives(
1197  const InputVector &fe_function,
1198  std::vector<typename ProductType<third_derivative_type,
1199  typename InputVector::value_type>::type>
1200  &third_derivatives) const;
1201 
1205  template <class InputVector>
1206  void
1207  get_function_third_derivatives_from_local_dof_values(
1208  const InputVector & dof_values,
1209  std::vector<typename OutputType<typename InputVector::value_type>::
1210  third_derivative_type> &third_derivatives) const;
1211 
1212  private:
1217 
1222  const unsigned int first_vector_component;
1223 
1227  std::vector<ShapeFunctionData> shape_function_data;
1228  };
1229 
1230 
1231  template <int rank, int dim, int spacedim = dim>
1232  class SymmetricTensor;
1233 
1258  template <int dim, int spacedim>
1259  class SymmetricTensor<2, dim, spacedim>
1260  {
1261  public:
1269 
1280 
1285  template <typename Number>
1286  struct OutputType
1287  {
1292  using value_type = typename ProductType<
1293  Number,
1294  typename SymmetricTensor<2, dim, spacedim>::value_type>::type;
1295 
1300  using divergence_type = typename ProductType<
1301  Number,
1302  typename SymmetricTensor<2, dim, spacedim>::divergence_type>::type;
1303  };
1304 
1309  struct ShapeFunctionData
1310  {
1319  bool is_nonzero_shape_function_component
1320  [value_type::n_independent_components];
1321 
1331  unsigned int row_index[value_type::n_independent_components];
1332 
1342 
1347  };
1348 
1352  SymmetricTensor();
1353 
1363  SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1364  const unsigned int first_tensor_component);
1365 
1370  SymmetricTensor &
1371  operator=(const SymmetricTensor<2, dim, spacedim> &) = delete;
1372 
1390  value_type
1391  value(const unsigned int shape_function, const unsigned int q_point) const;
1392 
1407  divergence(const unsigned int shape_function,
1408  const unsigned int q_point) const;
1409 
1427  template <class InputVector>
1428  void
1429  get_function_values(
1430  const InputVector &fe_function,
1431  std::vector<typename ProductType<value_type,
1432  typename InputVector::value_type>::type>
1433  &values) const;
1434 
1457  template <class InputVector>
1458  void
1459  get_function_values_from_local_dof_values(
1460  const InputVector &dof_values,
1461  std::vector<
1463  &values) const;
1464 
1486  template <class InputVector>
1487  void
1488  get_function_divergences(
1489  const InputVector &fe_function,
1490  std::vector<typename ProductType<divergence_type,
1491  typename InputVector::value_type>::type>
1492  &divergences) const;
1493 
1497  template <class InputVector>
1498  void
1499  get_function_divergences_from_local_dof_values(
1500  const InputVector &dof_values,
1501  std::vector<
1503  &divergences) const;
1504 
1505  private:
1510 
1515  const unsigned int first_tensor_component;
1516 
1520  std::vector<ShapeFunctionData> shape_function_data;
1521  };
1522 
1523 
1524  template <int rank, int dim, int spacedim = dim>
1525  class Tensor;
1526 
1547  template <int dim, int spacedim>
1548  class Tensor<2, dim, spacedim>
1549  {
1550  public:
1556 
1561 
1567 
1572  template <typename Number>
1573  struct OutputType
1574  {
1579  using value_type = typename ProductType<
1580  Number,
1582 
1587  using divergence_type = typename ProductType<
1588  Number,
1589  typename Tensor<2, dim, spacedim>::divergence_type>::type;
1590 
1595  using gradient_type = typename ProductType<
1596  Number,
1597  typename Tensor<2, dim, spacedim>::gradient_type>::type;
1598  };
1599 
1604  struct ShapeFunctionData
1605  {
1614  bool is_nonzero_shape_function_component
1615  [value_type::n_independent_components];
1616 
1626  unsigned int row_index[value_type::n_independent_components];
1627 
1637 
1642  };
1643 
1647  Tensor();
1648 
1658  Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1659  const unsigned int first_tensor_component);
1660 
1661 
1666  Tensor &
1667  operator=(const Tensor<2, dim, spacedim> &) = delete;
1668 
1685  value_type
1686  value(const unsigned int shape_function, const unsigned int q_point) const;
1687 
1702  divergence(const unsigned int shape_function,
1703  const unsigned int q_point) const;
1704 
1719  gradient(const unsigned int shape_function,
1720  const unsigned int q_point) const;
1721 
1739  template <class InputVector>
1740  void
1741  get_function_values(
1742  const InputVector &fe_function,
1743  std::vector<typename ProductType<value_type,
1744  typename InputVector::value_type>::type>
1745  &values) const;
1746 
1769  template <class InputVector>
1770  void
1771  get_function_values_from_local_dof_values(
1772  const InputVector &dof_values,
1773  std::vector<
1775  &values) const;
1776 
1798  template <class InputVector>
1799  void
1800  get_function_divergences(
1801  const InputVector &fe_function,
1802  std::vector<typename ProductType<divergence_type,
1803  typename InputVector::value_type>::type>
1804  &divergences) const;
1805 
1809  template <class InputVector>
1810  void
1811  get_function_divergences_from_local_dof_values(
1812  const InputVector &dof_values,
1813  std::vector<
1815  &divergences) const;
1816 
1833  template <class InputVector>
1834  void
1835  get_function_gradients(
1836  const InputVector &fe_function,
1837  std::vector<typename ProductType<gradient_type,
1838  typename InputVector::value_type>::type>
1839  &gradients) const;
1840 
1844  template <class InputVector>
1845  void
1846  get_function_gradients_from_local_dof_values(
1847  const InputVector &dof_values,
1848  std::vector<
1850  &gradients) const;
1851 
1852  private:
1857 
1862  const unsigned int first_tensor_component;
1863 
1867  std::vector<ShapeFunctionData> shape_function_data;
1868  };
1869 
1870 } // namespace FEValuesViews
1871 
1872 
1873 namespace internal
1874 {
1875  namespace FEValuesViews
1876  {
1883  template <int dim, int spacedim, typename Extractor>
1884  struct ViewType
1885  {};
1886 
1894  template <int dim, int spacedim>
1895  struct ViewType<dim, spacedim, FEValuesExtractors::Scalar>
1896  {
1897  using type = typename ::FEValuesViews::Scalar<dim, spacedim>;
1898  };
1899 
1907  template <int dim, int spacedim>
1908  struct ViewType<dim, spacedim, FEValuesExtractors::Vector>
1909  {
1910  using type = typename ::FEValuesViews::Vector<dim, spacedim>;
1911  };
1912 
1920  template <int dim, int spacedim, int rank>
1921  struct ViewType<dim, spacedim, FEValuesExtractors::Tensor<rank>>
1922  {
1923  using type = typename ::FEValuesViews::Tensor<rank, dim, spacedim>;
1924  };
1925 
1933  template <int dim, int spacedim, int rank>
1934  struct ViewType<dim, spacedim, FEValuesExtractors::SymmetricTensor<rank>>
1935  {
1936  using type =
1937  typename ::FEValuesViews::SymmetricTensor<rank, dim, spacedim>;
1938  };
1939 
1947  template <int dim, int spacedim>
1948  struct Cache
1949  {
1954  std::vector<::FEValuesViews::Scalar<dim, spacedim>> scalars;
1955  std::vector<::FEValuesViews::Vector<dim, spacedim>> vectors;
1956  std::vector<::FEValuesViews::SymmetricTensor<2, dim, spacedim>>
1957  symmetric_second_order_tensors;
1958  std::vector<::FEValuesViews::Tensor<2, dim, spacedim>>
1959  second_order_tensors;
1960 
1964  Cache(const FEValuesBase<dim, spacedim> &fe_values);
1965  };
1966  } // namespace FEValuesViews
1967 } // namespace internal
1968 
1969 namespace FEValuesViews
1970 {
1977  template <int dim, int spacedim, typename Extractor>
1978  using View = typename ::internal::FEValuesViews::
1979  ViewType<dim, spacedim, Extractor>::type;
1980 } // namespace FEValuesViews
1981 
1982 
2083 template <int dim, int spacedim>
2084 class FEValuesBase : public Subscriptor
2085 {
2086 public:
2090  static const unsigned int dimension = dim;
2091 
2095  static const unsigned int space_dimension = spacedim;
2096 
2100  const unsigned int n_quadrature_points;
2101 
2107  const unsigned int dofs_per_cell;
2108 
2109 
2117  FEValuesBase(const unsigned int n_q_points,
2118  const unsigned int dofs_per_cell,
2119  const UpdateFlags update_flags,
2120  const Mapping<dim, spacedim> & mapping,
2121  const FiniteElement<dim, spacedim> &fe);
2122 
2127  FEValuesBase &
2128  operator=(const FEValuesBase &) = delete;
2129 
2134  FEValuesBase(const FEValuesBase &) = delete;
2135 
2139  virtual ~FEValuesBase() override;
2140 
2141 
2145 
2146 
2167  const double &
2168  shape_value(const unsigned int function_no,
2169  const unsigned int point_no) const;
2170 
2191  double
2192  shape_value_component(const unsigned int function_no,
2193  const unsigned int point_no,
2194  const unsigned int component) const;
2195 
2221  const Tensor<1, spacedim> &
2222  shape_grad(const unsigned int function_no,
2223  const unsigned int quadrature_point) const;
2224 
2242  shape_grad_component(const unsigned int function_no,
2243  const unsigned int point_no,
2244  const unsigned int component) const;
2245 
2265  const Tensor<2, spacedim> &
2266  shape_hessian(const unsigned int function_no,
2267  const unsigned int point_no) const;
2268 
2286  shape_hessian_component(const unsigned int function_no,
2287  const unsigned int point_no,
2288  const unsigned int component) const;
2289 
2309  const Tensor<3, spacedim> &
2310  shape_3rd_derivative(const unsigned int function_no,
2311  const unsigned int point_no) const;
2312 
2330  shape_3rd_derivative_component(const unsigned int function_no,
2331  const unsigned int point_no,
2332  const unsigned int component) const;
2333 
2335 
2337 
2374  template <class InputVector>
2375  void
2376  get_function_values(
2377  const InputVector & fe_function,
2378  std::vector<typename InputVector::value_type> &values) const;
2379 
2393  template <class InputVector>
2394  void
2395  get_function_values(
2396  const InputVector & fe_function,
2397  std::vector<Vector<typename InputVector::value_type>> &values) const;
2398 
2417  template <class InputVector>
2418  void
2419  get_function_values(
2420  const InputVector & fe_function,
2422  std::vector<typename InputVector::value_type> & values) const;
2423 
2445  template <class InputVector>
2446  void
2447  get_function_values(
2448  const InputVector & fe_function,
2450  std::vector<Vector<typename InputVector::value_type>> &values) const;
2451 
2452 
2483  template <class InputVector>
2484  void
2485  get_function_values(
2486  const InputVector & fe_function,
2488  ArrayView<std::vector<typename InputVector::value_type>> values,
2489  const bool quadrature_points_fastest) const;
2490 
2492 
2494 
2531  template <class InputVector>
2532  void
2533  get_function_gradients(
2534  const InputVector &fe_function,
2536  &gradients) const;
2537 
2554  template <class InputVector>
2555  void
2556  get_function_gradients(
2557  const InputVector &fe_function,
2558  std::vector<
2560  &gradients) const;
2561 
2568  template <class InputVector>
2569  void
2570  get_function_gradients(
2571  const InputVector & fe_function,
2574  &gradients) const;
2575 
2582  template <class InputVector>
2583  void
2584  get_function_gradients(
2585  const InputVector & fe_function,
2587  ArrayView<
2589  gradients,
2590  bool quadrature_points_fastest = false) const;
2591 
2593 
2597 
2635  template <class InputVector>
2636  void
2637  get_function_hessians(
2638  const InputVector &fe_function,
2640  &hessians) const;
2641 
2659  template <class InputVector>
2660  void
2661  get_function_hessians(
2662  const InputVector &fe_function,
2663  std::vector<
2665  & hessians,
2666  bool quadrature_points_fastest = false) const;
2667 
2672  template <class InputVector>
2673  void
2674  get_function_hessians(
2675  const InputVector & fe_function,
2678  &hessians) const;
2679 
2686  template <class InputVector>
2687  void
2688  get_function_hessians(
2689  const InputVector & fe_function,
2691  ArrayView<
2693  hessians,
2694  bool quadrature_points_fastest = false) const;
2695 
2736  template <class InputVector>
2737  void
2738  get_function_laplacians(
2739  const InputVector & fe_function,
2740  std::vector<typename InputVector::value_type> &laplacians) const;
2741 
2761  template <class InputVector>
2762  void
2763  get_function_laplacians(
2764  const InputVector & fe_function,
2765  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2766 
2773  template <class InputVector>
2774  void
2775  get_function_laplacians(
2776  const InputVector & fe_function,
2778  std::vector<typename InputVector::value_type> & laplacians) const;
2779 
2786  template <class InputVector>
2787  void
2788  get_function_laplacians(
2789  const InputVector & fe_function,
2791  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2792 
2799  template <class InputVector>
2800  void
2801  get_function_laplacians(
2802  const InputVector & fe_function,
2804  std::vector<std::vector<typename InputVector::value_type>> &laplacians,
2805  bool quadrature_points_fastest = false) const;
2806 
2808 
2810 
2849  template <class InputVector>
2850  void
2851  get_function_third_derivatives(
2852  const InputVector &fe_function,
2854  &third_derivatives) const;
2855 
2874  template <class InputVector>
2875  void
2876  get_function_third_derivatives(
2877  const InputVector &fe_function,
2878  std::vector<
2880  & third_derivatives,
2881  bool quadrature_points_fastest = false) const;
2882 
2887  template <class InputVector>
2888  void
2889  get_function_third_derivatives(
2890  const InputVector & fe_function,
2893  &third_derivatives) const;
2894 
2901  template <class InputVector>
2902  void
2903  get_function_third_derivatives(
2904  const InputVector & fe_function,
2906  ArrayView<
2908  third_derivatives,
2909  bool quadrature_points_fastest = false) const;
2911 
2913 
2914 
2920  const Point<spacedim> &
2921  quadrature_point(const unsigned int q) const;
2922 
2928  const std::vector<Point<spacedim>> &
2929  get_quadrature_points() const;
2930 
2946  double
2947  JxW(const unsigned int quadrature_point) const;
2948 
2952  const std::vector<double> &
2953  get_JxW_values() const;
2954 
2962  jacobian(const unsigned int quadrature_point) const;
2963 
2970  const std::vector<DerivativeForm<1, dim, spacedim>> &
2971  get_jacobians() const;
2972 
2981  jacobian_grad(const unsigned int quadrature_point) const;
2982 
2989  const std::vector<DerivativeForm<2, dim, spacedim>> &
2990  get_jacobian_grads() const;
2991 
3000  const Tensor<3, spacedim> &
3001  jacobian_pushed_forward_grad(const unsigned int quadrature_point) const;
3002 
3009  const std::vector<Tensor<3, spacedim>> &
3010  get_jacobian_pushed_forward_grads() const;
3011 
3020  jacobian_2nd_derivative(const unsigned int quadrature_point) const;
3021 
3028  const std::vector<DerivativeForm<3, dim, spacedim>> &
3029  get_jacobian_2nd_derivatives() const;
3030 
3040  const Tensor<4, spacedim> &
3041  jacobian_pushed_forward_2nd_derivative(
3042  const unsigned int quadrature_point) const;
3043 
3050  const std::vector<Tensor<4, spacedim>> &
3051  get_jacobian_pushed_forward_2nd_derivatives() const;
3052 
3062  jacobian_3rd_derivative(const unsigned int quadrature_point) const;
3063 
3070  const std::vector<DerivativeForm<4, dim, spacedim>> &
3071  get_jacobian_3rd_derivatives() const;
3072 
3082  const Tensor<5, spacedim> &
3083  jacobian_pushed_forward_3rd_derivative(
3084  const unsigned int quadrature_point) const;
3085 
3092  const std::vector<Tensor<5, spacedim>> &
3093  get_jacobian_pushed_forward_3rd_derivatives() const;
3094 
3102  inverse_jacobian(const unsigned int quadrature_point) const;
3103 
3110  const std::vector<DerivativeForm<1, spacedim, dim>> &
3111  get_inverse_jacobians() const;
3112 
3132  const Tensor<1, spacedim> &
3133  normal_vector(const unsigned int i) const;
3134 
3145  DEAL_II_DEPRECATED
3146  const std::vector<Tensor<1, spacedim>> &
3147  get_all_normal_vectors() const;
3148 
3156  const std::vector<Tensor<1, spacedim>> &
3157  get_normal_vectors() const;
3158 
3160 
3162 
3163 
3173  operator[](const FEValuesExtractors::Scalar &scalar) const;
3174 
3184  operator[](const FEValuesExtractors::Vector &vector) const;
3185 
3196  operator[](const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
3197 
3198 
3208  operator[](const FEValuesExtractors::Tensor<2> &tensor) const;
3209 
3211 
3213 
3214 
3218  const Mapping<dim, spacedim> &
3219  get_mapping() const;
3220 
3225  get_fe() const;
3226 
3230  UpdateFlags
3231  get_update_flags() const;
3232 
3237  get_cell() const;
3238 
3245  get_cell_similarity() const;
3246 
3251  std::size_t
3252  memory_consumption() const;
3254 
3255 
3264  std::string,
3265  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3266  << "object for which this kind of information has not been computed. What "
3267  << "information these objects compute is determined by the update_* flags you "
3268  << "pass to the constructor. Here, the operation you are attempting requires "
3269  << "the <" << arg1
3270  << "> flag to be set, but it was apparently not specified "
3271  << "upon construction.");
3272 
3280  ExcFEDontMatch,
3281  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3282  "to the DoFHandler that provided the cell iterator do not match.");
3288  DeclException1(ExcShapeFunctionNotPrimitive,
3289  int,
3290  << "The shape function with index " << arg1
3291  << " is not primitive, i.e. it is vector-valued and "
3292  << "has more than one non-zero vector component. This "
3293  << "function cannot be called for these shape functions. "
3294  << "Maybe you want to use the same function with the "
3295  << "_component suffix?");
3296 
3305  "The given FiniteElement is not a primitive element but the requested operation "
3306  "only works for those. See FiniteElement::is_primitive() for more information.");
3307 
3308 protected:
3337  class CellIteratorBase;
3338 
3343  template <typename CI>
3345  class TriaCellIterator;
3346 
3352  std::unique_ptr<const CellIteratorBase> present_cell;
3353 
3361  boost::signals2::connection tria_listener_refinement;
3362 
3370  boost::signals2::connection tria_listener_mesh_transform;
3371 
3377  void
3378  invalidate_present_cell();
3379 
3389  void
3390  maybe_invalidate_previous_present_cell(
3391  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3392 
3398 
3404  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3406 
3413 
3414 
3422 
3428  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3430 
3436  spacedim>
3438 
3439 
3444 
3453  UpdateFlags
3454  compute_update_flags(const UpdateFlags update_flags) const;
3455 
3462 
3468  void
3469  check_cell_similarity(
3470  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3471 
3472 private:
3477 
3478  // Make the view classes friends of this class, since they access internal
3479  // data.
3480  template <int, int>
3481  friend class FEValuesViews::Scalar;
3482  template <int, int>
3483  friend class FEValuesViews::Vector;
3484  template <int, int, int>
3485  friend class FEValuesViews::SymmetricTensor;
3486  template <int, int, int>
3487  friend class FEValuesViews::Tensor;
3488 };
3489 
3490 
3491 
3502 template <int dim, int spacedim = dim>
3503 class FEValues : public FEValuesBase<dim, spacedim>
3504 {
3505 public:
3510  static const unsigned int integral_dimension = dim;
3511 
3516  FEValues(const Mapping<dim, spacedim> & mapping,
3517  const FiniteElement<dim, spacedim> &fe,
3518  const Quadrature<dim> & quadrature,
3519  const UpdateFlags update_flags);
3520 
3527  const Quadrature<dim> & quadrature,
3528  const UpdateFlags update_flags);
3529 
3536  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3537  void
3538  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3539  level_dof_access>> &cell);
3540 
3554  void
3555  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3556 
3561  const Quadrature<dim> &
3562  get_quadrature() const;
3563 
3568  std::size_t
3569  memory_consumption() const;
3570 
3585  const FEValues<dim, spacedim> &
3586  get_present_fe_values() const;
3587 
3588 private:
3593 
3597  void
3598  initialize(const UpdateFlags update_flags);
3599 
3606  void
3607  do_reinit();
3608 };
3609 
3610 
3621 template <int dim, int spacedim = dim>
3622 class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
3623 {
3624 public:
3629  static const unsigned int integral_dimension = dim - 1;
3630 
3642  FEFaceValuesBase(const unsigned int n_q_points,
3643  const unsigned int dofs_per_cell,
3644  const UpdateFlags update_flags,
3645  const Mapping<dim, spacedim> & mapping,
3646  const FiniteElement<dim, spacedim> &fe,
3647  const Quadrature<dim - 1> & quadrature);
3648 
3656  const Tensor<1, spacedim> &
3657  boundary_form(const unsigned int i) const;
3658 
3665  const std::vector<Tensor<1, spacedim>> &
3666  get_boundary_forms() const;
3667 
3672  unsigned int
3673  get_face_index() const;
3674 
3679  const Quadrature<dim - 1> &
3680  get_quadrature() const;
3681 
3686  std::size_t
3687  memory_consumption() const;
3688 
3689 protected:
3694  unsigned int present_face_index;
3695 
3699  const Quadrature<dim - 1> quadrature;
3700 };
3701 
3702 
3703 
3718 template <int dim, int spacedim = dim>
3719 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
3720 {
3721 public:
3726  static const unsigned int dimension = dim;
3727 
3728  static const unsigned int space_dimension = spacedim;
3729 
3734  static const unsigned int integral_dimension = dim - 1;
3735 
3740  FEFaceValues(const Mapping<dim, spacedim> & mapping,
3741  const FiniteElement<dim, spacedim> &fe,
3742  const Quadrature<dim - 1> & quadrature,
3743  const UpdateFlags update_flags);
3744 
3751  const Quadrature<dim - 1> & quadrature,
3752  const UpdateFlags update_flags);
3753 
3758  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3759  void
3760  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3761  level_dof_access>> &cell,
3762  const unsigned int face_no);
3763 
3770  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3771  void
3772  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3773  level_dof_access>> & cell,
3774  const typename Triangulation<dim, spacedim>::face_iterator &face);
3775 
3789  void
3790  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3791  const unsigned int face_no);
3792 
3793  /*
3794  * Reinitialize the gradients, Jacobi determinants, etc for the given face
3795  * on a given cell of type "iterator into a Triangulation object", and the
3796  * given finite element. Since iterators into a triangulation alone only
3797  * convey information about the geometry of a cell, but not about degrees of
3798  * freedom possibly associated with this cell, you will not be able to call
3799  * some functions of this class if they need information about degrees of
3800  * freedom. These functions are, above all, the
3801  * <tt>get_function_value/gradients/hessians/third_derivatives</tt>
3802  * functions. If you want to call these functions, you have to call the @p
3803  * reinit variants that take iterators into DoFHandler or other DoF handler
3804  * type objects.
3805  *
3806  * @note @p face must be one of @p cell's face iterators.
3807  */
3808  void
3809  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3810  const typename Triangulation<dim, spacedim>::face_iterator &face);
3811 
3827  get_present_fe_values() const;
3828 
3829 private:
3833  void
3834  initialize(const UpdateFlags update_flags);
3835 
3842  void
3843  do_reinit(const unsigned int face_no);
3844 };
3845 
3846 
3864 template <int dim, int spacedim = dim>
3865 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
3866 {
3867 public:
3871  static const unsigned int dimension = dim;
3872 
3876  static const unsigned int space_dimension = spacedim;
3877 
3882  static const unsigned int integral_dimension = dim - 1;
3883 
3888  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
3889  const FiniteElement<dim, spacedim> &fe,
3890  const Quadrature<dim - 1> & face_quadrature,
3891  const UpdateFlags update_flags);
3892 
3899  const Quadrature<dim - 1> & face_quadrature,
3900  const UpdateFlags update_flags);
3901 
3908  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3909  void
3910  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3911  level_dof_access>> &cell,
3912  const unsigned int face_no,
3913  const unsigned int subface_no);
3914 
3919  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3920  void
3921  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3922  level_dof_access>> & cell,
3923  const typename Triangulation<dim, spacedim>::face_iterator &face,
3924  const typename Triangulation<dim, spacedim>::face_iterator &subface);
3925 
3939  void
3940  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3941  const unsigned int face_no,
3942  const unsigned int subface_no);
3943 
3963  void
3964  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3965  const typename Triangulation<dim, spacedim>::face_iterator &face,
3966  const typename Triangulation<dim, spacedim>::face_iterator &subface);
3967 
3983  get_present_fe_values() const;
3984 
3990  DeclException0(ExcReinitCalledWithBoundaryFace);
3991 
3997  DeclException0(ExcFaceHasNoSubfaces);
3998 
3999 private:
4003  void
4004  initialize(const UpdateFlags update_flags);
4005 
4012  void
4013  do_reinit(const unsigned int face_no, const unsigned int subface_no);
4014 };
4015 
4016 
4017 #ifndef DOXYGEN
4018 
4019 
4020 /*------------------------ Inline functions: namespace FEValuesViews --------*/
4021 
4022 namespace FEValuesViews
4023 {
4024  template <int dim, int spacedim>
4025  inline typename Scalar<dim, spacedim>::value_type
4026  Scalar<dim, spacedim>::value(const unsigned int shape_function,
4027  const unsigned int q_point) const
4028  {
4029  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4030  Assert(
4031  fe_values->update_flags & update_values,
4033  "update_values"))));
4034 
4035  // an adaptation of the FEValuesBase::shape_value_component function
4036  // except that here we know the component as fixed and we have
4037  // pre-computed and cached a bunch of information. See the comments there.
4038  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4039  return fe_values->finite_element_output.shape_values(
4040  shape_function_data[shape_function].row_index, q_point);
4041  else
4042  return 0;
4043  }
4044 
4045 
4046 
4047  template <int dim, int spacedim>
4048  inline typename Scalar<dim, spacedim>::gradient_type
4049  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
4050  const unsigned int q_point) const
4051  {
4052  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4053  Assert(fe_values->update_flags & update_gradients,
4055  "update_gradients")));
4056 
4057  // an adaptation of the FEValuesBase::shape_grad_component
4058  // function except that here we know the component as fixed and we have
4059  // pre-computed and cached a bunch of information. See the comments there.
4060  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4061  return fe_values->finite_element_output
4062  .shape_gradients[shape_function_data[shape_function].row_index]
4063  [q_point];
4064  else
4065  return gradient_type();
4066  }
4067 
4068 
4069 
4070  template <int dim, int spacedim>
4071  inline typename Scalar<dim, spacedim>::hessian_type
4072  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
4073  const unsigned int q_point) const
4074  {
4075  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4076  Assert(fe_values->update_flags & update_hessians,
4078  "update_hessians")));
4079 
4080  // an adaptation of the FEValuesBase::shape_hessian_component
4081  // function except that here we know the component as fixed and we have
4082  // pre-computed and cached a bunch of information. See the comments there.
4083  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4084  return fe_values->finite_element_output
4085  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
4086  else
4087  return hessian_type();
4088  }
4089 
4090 
4091 
4092  template <int dim, int spacedim>
4093  inline typename Scalar<dim, spacedim>::third_derivative_type
4094  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
4095  const unsigned int q_point) const
4096  {
4097  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4098  Assert(fe_values->update_flags & update_3rd_derivatives,
4100  "update_3rd_derivatives")));
4101 
4102  // an adaptation of the FEValuesBase::shape_3rdderivative_component
4103  // function except that here we know the component as fixed and we have
4104  // pre-computed and cached a bunch of information. See the comments there.
4105  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4106  return fe_values->finite_element_output
4107  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
4108  [q_point];
4109  else
4110  return third_derivative_type();
4111  }
4112 
4113 
4114 
4115  template <int dim, int spacedim>
4116  inline typename Vector<dim, spacedim>::value_type
4117  Vector<dim, spacedim>::value(const unsigned int shape_function,
4118  const unsigned int q_point) const
4119  {
4120  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4121  Assert(fe_values->update_flags & update_values,
4123  "update_values")));
4124 
4125  // same as for the scalar case except that we have one more index
4126  const int snc =
4127  shape_function_data[shape_function].single_nonzero_component;
4128  if (snc == -2)
4129  return value_type();
4130  else if (snc != -1)
4131  {
4132  value_type return_value;
4133  return_value[shape_function_data[shape_function]
4134  .single_nonzero_component_index] =
4135  fe_values->finite_element_output.shape_values(snc, q_point);
4136  return return_value;
4137  }
4138  else
4139  {
4140  value_type return_value;
4141  for (unsigned int d = 0; d < dim; ++d)
4142  if (shape_function_data[shape_function]
4143  .is_nonzero_shape_function_component[d])
4144  return_value[d] = fe_values->finite_element_output.shape_values(
4145  shape_function_data[shape_function].row_index[d], q_point);
4146 
4147  return return_value;
4148  }
4149  }
4150 
4151 
4152 
4153  template <int dim, int spacedim>
4154  inline typename Vector<dim, spacedim>::gradient_type
4155  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4156  const unsigned int q_point) const
4157  {
4158  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4159  Assert(fe_values->update_flags & update_gradients,
4161  "update_gradients")));
4162 
4163  // same as for the scalar case except that we have one more index
4164  const int snc =
4165  shape_function_data[shape_function].single_nonzero_component;
4166  if (snc == -2)
4167  return gradient_type();
4168  else if (snc != -1)
4169  {
4170  gradient_type return_value;
4171  return_value[shape_function_data[shape_function]
4172  .single_nonzero_component_index] =
4173  fe_values->finite_element_output.shape_gradients[snc][q_point];
4174  return return_value;
4175  }
4176  else
4177  {
4178  gradient_type return_value;
4179  for (unsigned int d = 0; d < dim; ++d)
4180  if (shape_function_data[shape_function]
4181  .is_nonzero_shape_function_component[d])
4182  return_value[d] =
4183  fe_values->finite_element_output.shape_gradients
4184  [shape_function_data[shape_function].row_index[d]][q_point];
4185 
4186  return return_value;
4187  }
4188  }
4189 
4190 
4191 
4192  template <int dim, int spacedim>
4193  inline typename Vector<dim, spacedim>::divergence_type
4194  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4195  const unsigned int q_point) const
4196  {
4197  // this function works like in the case above
4198  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4199  Assert(fe_values->update_flags & update_gradients,
4201  "update_gradients")));
4202 
4203  // same as for the scalar case except that we have one more index
4204  const int snc =
4205  shape_function_data[shape_function].single_nonzero_component;
4206  if (snc == -2)
4207  return divergence_type();
4208  else if (snc != -1)
4209  return fe_values->finite_element_output
4210  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4211  .single_nonzero_component_index];
4212  else
4213  {
4214  divergence_type return_value = 0;
4215  for (unsigned int d = 0; d < dim; ++d)
4216  if (shape_function_data[shape_function]
4217  .is_nonzero_shape_function_component[d])
4218  return_value +=
4219  fe_values->finite_element_output.shape_gradients
4220  [shape_function_data[shape_function].row_index[d]][q_point][d];
4221 
4222  return return_value;
4223  }
4224  }
4225 
4226 
4227 
4228  template <int dim, int spacedim>
4229  inline typename Vector<dim, spacedim>::curl_type
4230  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4231  const unsigned int q_point) const
4232  {
4233  // this function works like in the case above
4234 
4235  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4236  Assert(fe_values->update_flags & update_gradients,
4238  "update_gradients")));
4239  // same as for the scalar case except that we have one more index
4240  const int snc =
4241  shape_function_data[shape_function].single_nonzero_component;
4242 
4243  if (snc == -2)
4244  return curl_type();
4245 
4246  else
4247  switch (dim)
4248  {
4249  case 1:
4250  {
4251  Assert(false,
4252  ExcMessage(
4253  "Computing the curl in 1d is not a useful operation"));
4254  return curl_type();
4255  }
4256 
4257  case 2:
4258  {
4259  if (snc != -1)
4260  {
4261  curl_type return_value;
4262 
4263  // the single nonzero component can only be zero or one in 2d
4264  if (shape_function_data[shape_function]
4265  .single_nonzero_component_index == 0)
4266  return_value[0] =
4267  -1.0 * fe_values->finite_element_output
4268  .shape_gradients[snc][q_point][1];
4269  else
4270  return_value[0] = fe_values->finite_element_output
4271  .shape_gradients[snc][q_point][0];
4272 
4273  return return_value;
4274  }
4275 
4276  else
4277  {
4278  curl_type return_value;
4279 
4280  return_value[0] = 0.0;
4281 
4282  if (shape_function_data[shape_function]
4283  .is_nonzero_shape_function_component[0])
4284  return_value[0] -=
4285  fe_values->finite_element_output
4286  .shape_gradients[shape_function_data[shape_function]
4287  .row_index[0]][q_point][1];
4288 
4289  if (shape_function_data[shape_function]
4290  .is_nonzero_shape_function_component[1])
4291  return_value[0] +=
4292  fe_values->finite_element_output
4293  .shape_gradients[shape_function_data[shape_function]
4294  .row_index[1]][q_point][0];
4295 
4296  return return_value;
4297  }
4298  }
4299 
4300  case 3:
4301  {
4302  if (snc != -1)
4303  {
4304  curl_type return_value;
4305 
4306  switch (shape_function_data[shape_function]
4307  .single_nonzero_component_index)
4308  {
4309  case 0:
4310  {
4311  return_value[0] = 0;
4312  return_value[1] = fe_values->finite_element_output
4313  .shape_gradients[snc][q_point][2];
4314  return_value[2] =
4315  -1.0 * fe_values->finite_element_output
4316  .shape_gradients[snc][q_point][1];
4317  return return_value;
4318  }
4319 
4320  case 1:
4321  {
4322  return_value[0] =
4323  -1.0 * fe_values->finite_element_output
4324  .shape_gradients[snc][q_point][2];
4325  return_value[1] = 0;
4326  return_value[2] = fe_values->finite_element_output
4327  .shape_gradients[snc][q_point][0];
4328  return return_value;
4329  }
4330 
4331  default:
4332  {
4333  return_value[0] = fe_values->finite_element_output
4334  .shape_gradients[snc][q_point][1];
4335  return_value[1] =
4336  -1.0 * fe_values->finite_element_output
4337  .shape_gradients[snc][q_point][0];
4338  return_value[2] = 0;
4339  return return_value;
4340  }
4341  }
4342  }
4343 
4344  else
4345  {
4346  curl_type return_value;
4347 
4348  for (unsigned int i = 0; i < dim; ++i)
4349  return_value[i] = 0.0;
4350 
4351  if (shape_function_data[shape_function]
4352  .is_nonzero_shape_function_component[0])
4353  {
4354  return_value[1] +=
4355  fe_values->finite_element_output
4356  .shape_gradients[shape_function_data[shape_function]
4357  .row_index[0]][q_point][2];
4358  return_value[2] -=
4359  fe_values->finite_element_output
4360  .shape_gradients[shape_function_data[shape_function]
4361  .row_index[0]][q_point][1];
4362  }
4363 
4364  if (shape_function_data[shape_function]
4365  .is_nonzero_shape_function_component[1])
4366  {
4367  return_value[0] -=
4368  fe_values->finite_element_output
4369  .shape_gradients[shape_function_data[shape_function]
4370  .row_index[1]][q_point][2];
4371  return_value[2] +=
4372  fe_values->finite_element_output
4373  .shape_gradients[shape_function_data[shape_function]
4374  .row_index[1]][q_point][0];
4375  }
4376 
4377  if (shape_function_data[shape_function]
4378  .is_nonzero_shape_function_component[2])
4379  {
4380  return_value[0] +=
4381  fe_values->finite_element_output
4382  .shape_gradients[shape_function_data[shape_function]
4383  .row_index[2]][q_point][1];
4384  return_value[1] -=
4385  fe_values->finite_element_output
4386  .shape_gradients[shape_function_data[shape_function]
4387  .row_index[2]][q_point][0];
4388  }
4389 
4390  return return_value;
4391  }
4392  }
4393  }
4394  // should not end up here
4395  Assert(false, ExcInternalError());
4396  return curl_type();
4397  }
4398 
4399 
4400 
4401  template <int dim, int spacedim>
4402  inline typename Vector<dim, spacedim>::hessian_type
4403  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
4404  const unsigned int q_point) const
4405  {
4406  // this function works like in the case above
4407  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4408  Assert(fe_values->update_flags & update_hessians,
4410  "update_hessians")));
4411 
4412  // same as for the scalar case except that we have one more index
4413  const int snc =
4414  shape_function_data[shape_function].single_nonzero_component;
4415  if (snc == -2)
4416  return hessian_type();
4417  else if (snc != -1)
4418  {
4419  hessian_type return_value;
4420  return_value[shape_function_data[shape_function]
4421  .single_nonzero_component_index] =
4422  fe_values->finite_element_output.shape_hessians[snc][q_point];
4423  return return_value;
4424  }
4425  else
4426  {
4427  hessian_type return_value;
4428  for (unsigned int d = 0; d < dim; ++d)
4429  if (shape_function_data[shape_function]
4430  .is_nonzero_shape_function_component[d])
4431  return_value[d] =
4432  fe_values->finite_element_output.shape_hessians
4433  [shape_function_data[shape_function].row_index[d]][q_point];
4434 
4435  return return_value;
4436  }
4437  }
4438 
4439 
4440 
4441  template <int dim, int spacedim>
4442  inline typename Vector<dim, spacedim>::third_derivative_type
4443  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
4444  const unsigned int q_point) const
4445  {
4446  // this function works like in the case above
4447  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4448  Assert(fe_values->update_flags & update_3rd_derivatives,
4450  "update_3rd_derivatives")));
4451 
4452  // same as for the scalar case except that we have one more index
4453  const int snc =
4454  shape_function_data[shape_function].single_nonzero_component;
4455  if (snc == -2)
4456  return third_derivative_type();
4457  else if (snc != -1)
4458  {
4459  third_derivative_type return_value;
4460  return_value[shape_function_data[shape_function]
4461  .single_nonzero_component_index] =
4462  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
4463  return return_value;
4464  }
4465  else
4466  {
4467  third_derivative_type return_value;
4468  for (unsigned int d = 0; d < dim; ++d)
4469  if (shape_function_data[shape_function]
4470  .is_nonzero_shape_function_component[d])
4471  return_value[d] =
4472  fe_values->finite_element_output.shape_3rd_derivatives
4473  [shape_function_data[shape_function].row_index[d]][q_point];
4474 
4475  return return_value;
4476  }
4477  }
4478 
4479 
4480 
4481  namespace internal
4482  {
4487  inline ::SymmetricTensor<2, 1>
4488  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
4489  {
4490  AssertIndexRange(n, 1);
4491  (void)n;
4492 
4493  return {{t[0]}};
4494  }
4495 
4496 
4497 
4498  inline ::SymmetricTensor<2, 2>
4499  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
4500  {
4501  switch (n)
4502  {
4503  case 0:
4504  {
4505  return {{t[0], 0, t[1] / 2}};
4506  }
4507  case 1:
4508  {
4509  return {{0, t[1], t[0] / 2}};
4510  }
4511  default:
4512  {
4513  AssertIndexRange(n, 2);
4514  return {};
4515  }
4516  }
4517  }
4518 
4519 
4520 
4521  inline ::SymmetricTensor<2, 3>
4522  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
4523  {
4524  switch (n)
4525  {
4526  case 0:
4527  {
4528  return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
4529  }
4530  case 1:
4531  {
4532  return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
4533  }
4534  case 2:
4535  {
4536  return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
4537  }
4538  default:
4539  {
4540  AssertIndexRange(n, 3);
4541  return {};
4542  }
4543  }
4544  }
4545  } // namespace internal
4546 
4547 
4548 
4549  template <int dim, int spacedim>
4550  inline typename Vector<dim, spacedim>::symmetric_gradient_type
4551  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
4552  const unsigned int q_point) const
4553  {
4554  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4555  Assert(fe_values->update_flags & update_gradients,
4557  "update_gradients")));
4558 
4559  // same as for the scalar case except that we have one more index
4560  const int snc =
4561  shape_function_data[shape_function].single_nonzero_component;
4562  if (snc == -2)
4563  return symmetric_gradient_type();
4564  else if (snc != -1)
4565  return internal::symmetrize_single_row(
4566  shape_function_data[shape_function].single_nonzero_component_index,
4567  fe_values->finite_element_output.shape_gradients[snc][q_point]);
4568  else
4569  {
4570  gradient_type return_value;
4571  for (unsigned int d = 0; d < dim; ++d)
4572  if (shape_function_data[shape_function]
4573  .is_nonzero_shape_function_component[d])
4574  return_value[d] =
4575  fe_values->finite_element_output.shape_gradients
4576  [shape_function_data[shape_function].row_index[d]][q_point];
4577 
4578  return symmetrize(return_value);
4579  }
4580  }
4581 
4582 
4583 
4584  template <int dim, int spacedim>
4586  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
4587  const unsigned int q_point) const
4588  {
4589  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4590  Assert(fe_values->update_flags & update_values,
4592  "update_values")));
4593 
4594  // similar to the vector case where we have more then one index and we need
4595  // to convert between unrolled and component indexing for tensors
4596  const int snc =
4597  shape_function_data[shape_function].single_nonzero_component;
4598 
4599  if (snc == -2)
4600  {
4601  // shape function is zero for the selected components
4602  return value_type();
4603  }
4604  else if (snc != -1)
4605  {
4606  value_type return_value;
4607  const unsigned int comp =
4608  shape_function_data[shape_function].single_nonzero_component_index;
4609  return_value[value_type::unrolled_to_component_indices(comp)] =
4610  fe_values->finite_element_output.shape_values(snc, q_point);
4611  return return_value;
4612  }
4613  else
4614  {
4615  value_type return_value;
4616  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
4617  if (shape_function_data[shape_function]
4618  .is_nonzero_shape_function_component[d])
4619  return_value[value_type::unrolled_to_component_indices(d)] =
4620  fe_values->finite_element_output.shape_values(
4621  shape_function_data[shape_function].row_index[d], q_point);
4622  return return_value;
4623  }
4624  }
4625 
4626 
4627 
4628  template <int dim, int spacedim>
4631  const unsigned int shape_function,
4632  const unsigned int q_point) const
4633  {
4634  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4635  Assert(fe_values->update_flags & update_gradients,
4637  "update_gradients")));
4638 
4639  const int snc =
4640  shape_function_data[shape_function].single_nonzero_component;
4641 
4642  if (snc == -2)
4643  {
4644  // shape function is zero for the selected components
4645  return divergence_type();
4646  }
4647  else if (snc != -1)
4648  {
4649  // we have a single non-zero component when the symmetric tensor is
4650  // represented in unrolled form. this implies we potentially have
4651  // two non-zero components when represented in component form! we
4652  // will only have one non-zero entry if the non-zero component lies on
4653  // the diagonal of the tensor.
4654  //
4655  // the divergence of a second-order tensor is a first order tensor.
4656  //
4657  // assume the second-order tensor is A with components A_{ij}. then
4658  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
4659  // entries in the tensorial representation. define the
4660  // divergence as:
4661  // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
4662  // (which is incidentally also
4663  // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
4664  // In both cases, a sum is implied.
4665  //
4666  // Now, we know the nonzero component in unrolled form: it is indicated
4667  // by 'snc'. we can figure out which tensor components belong to this:
4668  const unsigned int comp =
4669  shape_function_data[shape_function].single_nonzero_component_index;
4670  const unsigned int ii =
4671  value_type::unrolled_to_component_indices(comp)[0];
4672  const unsigned int jj =
4673  value_type::unrolled_to_component_indices(comp)[1];
4674 
4675  // given the form of the divergence above, if ii=jj there is only a
4676  // single nonzero component of the full tensor and the gradient
4677  // equals
4678  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
4679  // all other entries of 'b' are zero
4680  //
4681  // on the other hand, if ii!=jj, then there are two nonzero entries in
4682  // the full tensor and
4683  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
4684  // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
4685  // again, all other entries of 'b' are zero
4686  const ::Tensor<1, spacedim> &phi_grad =
4687  fe_values->finite_element_output.shape_gradients[snc][q_point];
4688 
4689  divergence_type return_value;
4690  return_value[ii] = phi_grad[jj];
4691 
4692  if (ii != jj)
4693  return_value[jj] = phi_grad[ii];
4694 
4695  return return_value;
4696  }
4697  else
4698  {
4699  Assert(false, ExcNotImplemented());
4700  divergence_type return_value;
4701  return return_value;
4702  }
4703  }
4704 
4705 
4706 
4707  template <int dim, int spacedim>
4708  inline typename Tensor<2, dim, spacedim>::value_type
4709  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
4710  const unsigned int q_point) const
4711  {
4712  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4713  Assert(fe_values->update_flags & update_values,
4715  "update_values")));
4716 
4717  // similar to the vector case where we have more then one index and we need
4718  // to convert between unrolled and component indexing for tensors
4719  const int snc =
4720  shape_function_data[shape_function].single_nonzero_component;
4721 
4722  if (snc == -2)
4723  {
4724  // shape function is zero for the selected components
4725  return value_type();
4726  }
4727  else if (snc != -1)
4728  {
4729  value_type return_value;
4730  const unsigned int comp =
4731  shape_function_data[shape_function].single_nonzero_component_index;
4732  const TableIndices<2> indices =
4734  return_value[indices] =
4735  fe_values->finite_element_output.shape_values(snc, q_point);
4736  return return_value;
4737  }
4738  else
4739  {
4740  value_type return_value;
4741  for (unsigned int d = 0; d < dim * dim; ++d)
4742  if (shape_function_data[shape_function]
4743  .is_nonzero_shape_function_component[d])
4744  {
4745  const TableIndices<2> indices =
4747  return_value[indices] =
4748  fe_values->finite_element_output.shape_values(
4749  shape_function_data[shape_function].row_index[d], q_point);
4750  }
4751  return return_value;
4752  }
4753  }
4754 
4755 
4756 
4757  template <int dim, int spacedim>
4759  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
4760  const unsigned int q_point) const
4761  {
4762  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4763  Assert(fe_values->update_flags & update_gradients,
4765  "update_gradients")));
4766 
4767  const int snc =
4768  shape_function_data[shape_function].single_nonzero_component;
4769 
4770  if (snc == -2)
4771  {
4772  // shape function is zero for the selected components
4773  return divergence_type();
4774  }
4775  else if (snc != -1)
4776  {
4777  // we have a single non-zero component when the tensor is
4778  // represented in unrolled form.
4779  //
4780  // the divergence of a second-order tensor is a first order tensor.
4781  //
4782  // assume the second-order tensor is A with components A_{ij},
4783  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
4784  //
4785  // Now, we know the nonzero component in unrolled form: it is indicated
4786  // by 'snc'. we can figure out which tensor components belong to this:
4787  const unsigned int comp =
4788  shape_function_data[shape_function].single_nonzero_component_index;
4789  const TableIndices<2> indices =
4791  const unsigned int ii = indices[0];
4792  const unsigned int jj = indices[1];
4793 
4794  const ::Tensor<1, spacedim> &phi_grad =
4795  fe_values->finite_element_output.shape_gradients[snc][q_point];
4796 
4797  divergence_type return_value;
4798  // note that we contract \nabla from the right
4799  return_value[ii] = phi_grad[jj];
4800 
4801  return return_value;
4802  }
4803  else
4804  {
4805  Assert(false, ExcNotImplemented());
4806  divergence_type return_value;
4807  return return_value;
4808  }
4809  }
4810 
4811 
4812 
4813  template <int dim, int spacedim>
4815  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
4816  const unsigned int q_point) const
4817  {
4818  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4819  Assert(fe_values->update_flags & update_gradients,
4821  "update_gradients")));
4822 
4823  const int snc =
4824  shape_function_data[shape_function].single_nonzero_component;
4825 
4826  if (snc == -2)
4827  {
4828  // shape function is zero for the selected components
4829  return gradient_type();
4830  }
4831  else if (snc != -1)
4832  {
4833  // we have a single non-zero component when the tensor is
4834  // represented in unrolled form.
4835  //
4836  // the gradient of a second-order tensor is a third order tensor.
4837  //
4838  // assume the second-order tensor is A with components A_{ij},
4839  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
4840  //
4841  // Now, we know the nonzero component in unrolled form: it is indicated
4842  // by 'snc'. we can figure out which tensor components belong to this:
4843  const unsigned int comp =
4844  shape_function_data[shape_function].single_nonzero_component_index;
4845  const TableIndices<2> indices =
4847  const unsigned int ii = indices[0];
4848  const unsigned int jj = indices[1];
4849 
4850  const ::Tensor<1, spacedim> &phi_grad =
4851  fe_values->finite_element_output.shape_gradients[snc][q_point];
4852 
4853  gradient_type return_value;
4854  return_value[ii][jj] = phi_grad;
4855 
4856  return return_value;
4857  }
4858  else
4859  {
4860  Assert(false, ExcNotImplemented());
4861  gradient_type return_value;
4862  return return_value;
4863  }
4864  }
4865 
4866 } // namespace FEValuesViews
4867 
4868 
4869 
4870 /*---------------------- Inline functions: FEValuesBase ---------------------*/
4871 
4872 
4873 
4874 template <int dim, int spacedim>
4876  operator[](const FEValuesExtractors::Scalar &scalar) const
4877 {
4878  AssertIndexRange(scalar.component, fe_values_views_cache.scalars.size());
4879 
4880  return fe_values_views_cache.scalars[scalar.component];
4881 }
4882 
4883 
4884 
4885 template <int dim, int spacedim>
4887  operator[](const FEValuesExtractors::Vector &vector) const
4888 {
4890  fe_values_views_cache.vectors.size());
4891 
4892  return fe_values_views_cache.vectors[vector.first_vector_component];
4893 }
4894 
4895 
4896 
4897 template <int dim, int spacedim>
4901 {
4902  Assert(
4903  tensor.first_tensor_component <
4904  fe_values_views_cache.symmetric_second_order_tensors.size(),
4906  0,
4907  fe_values_views_cache.symmetric_second_order_tensors.size()));
4908 
4909  return fe_values_views_cache
4910  .symmetric_second_order_tensors[tensor.first_tensor_component];
4911 }
4912 
4913 
4914 
4915 template <int dim, int spacedim>
4918  operator[](const FEValuesExtractors::Tensor<2> &tensor) const
4919 {
4921  fe_values_views_cache.second_order_tensors.size());
4922 
4923  return fe_values_views_cache
4924  .second_order_tensors[tensor.first_tensor_component];
4925 }
4926 
4927 
4928 
4929 template <int dim, int spacedim>
4930 inline const double &
4931 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
4932  const unsigned int j) const
4933 {
4934  AssertIndexRange(i, fe->dofs_per_cell);
4935  Assert(this->update_flags & update_values,
4936  ExcAccessToUninitializedField("update_values"));
4937  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4938  Assert(present_cell.get() != nullptr,
4939  ExcMessage("FEValues object is not reinit'ed to any cell"));
4940  // if the entire FE is primitive,
4941  // then we can take a short-cut:
4942  if (fe->is_primitive())
4943  return this->finite_element_output.shape_values(i, j);
4944  else
4945  {
4946  // otherwise, use the mapping
4947  // between shape function
4948  // numbers and rows. note that
4949  // by the assertions above, we
4950  // know that this particular
4951  // shape function is primitive,
4952  // so we can call
4953  // system_to_component_index
4954  const unsigned int row =
4955  this->finite_element_output
4956  .shape_function_to_row_table[i * fe->n_components() +
4957  fe->system_to_component_index(i).first];
4958  return this->finite_element_output.shape_values(row, j);
4959  }
4960 }
4961 
4962 
4963 
4964 template <int dim, int spacedim>
4965 inline double
4967  const unsigned int i,
4968  const unsigned int j,
4969  const unsigned int component) const
4970 {
4971  AssertIndexRange(i, fe->dofs_per_cell);
4972  Assert(this->update_flags & update_values,
4973  ExcAccessToUninitializedField("update_values"));
4974  AssertIndexRange(component, fe->n_components());
4975  Assert(present_cell.get() != nullptr,
4976  ExcMessage("FEValues object is not reinit'ed to any cell"));
4977 
4978  // check whether the shape function
4979  // is non-zero at all within
4980  // this component:
4981  if (fe->get_nonzero_components(i)[component] == false)
4982  return 0;
4983 
4984  // look up the right row in the
4985  // table and take the data from
4986  // there
4987  const unsigned int row =
4988  this->finite_element_output
4989  .shape_function_to_row_table[i * fe->n_components() + component];
4990  return this->finite_element_output.shape_values(row, j);
4991 }
4992 
4993 
4994 
4995 template <int dim, int spacedim>
4996 inline const Tensor<1, spacedim> &
4997 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
4998  const unsigned int j) const
4999 {
5000  AssertIndexRange(i, fe->dofs_per_cell);
5001  Assert(this->update_flags & update_gradients,
5002  ExcAccessToUninitializedField("update_gradients"));
5003  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5004  Assert(present_cell.get() != nullptr,
5005  ExcMessage("FEValues object is not reinit'ed to any cell"));
5006  // if the entire FE is primitive,
5007  // then we can take a short-cut:
5008  if (fe->is_primitive())
5009  return this->finite_element_output.shape_gradients[i][j];
5010  else
5011  {
5012  // otherwise, use the mapping
5013  // between shape function
5014  // numbers and rows. note that
5015  // by the assertions above, we
5016  // know that this particular
5017  // shape function is primitive,
5018  // so we can call
5019  // system_to_component_index
5020  const unsigned int row =
5021  this->finite_element_output
5022  .shape_function_to_row_table[i * fe->n_components() +
5023  fe->system_to_component_index(i).first];
5024  return this->finite_element_output.shape_gradients[row][j];
5025  }
5026 }
5027 
5028 
5029 
5030 template <int dim, int spacedim>
5031 inline Tensor<1, spacedim>
5033  const unsigned int i,
5034  const unsigned int j,
5035  const unsigned int component) const
5036 {
5037  AssertIndexRange(i, fe->dofs_per_cell);
5038  Assert(this->update_flags & update_gradients,
5039  ExcAccessToUninitializedField("update_gradients"));
5040  AssertIndexRange(component, fe->n_components());
5041  Assert(present_cell.get() != nullptr,
5042  ExcMessage("FEValues object is not reinit'ed to any cell"));
5043  // check whether the shape function
5044  // is non-zero at all within
5045  // this component:
5046  if (fe->get_nonzero_components(i)[component] == false)
5047  return Tensor<1, spacedim>();
5048 
5049  // look up the right row in the
5050  // table and take the data from
5051  // there
5052  const unsigned int row =
5053  this->finite_element_output
5054  .shape_function_to_row_table[i * fe->n_components() + component];
5055  return this->finite_element_output.shape_gradients[row][j];
5056 }
5057 
5058 
5059 
5060 template <int dim, int spacedim>
5061 inline const Tensor<2, spacedim> &
5062 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
5063  const unsigned int j) const
5064 {
5065  AssertIndexRange(i, fe->dofs_per_cell);
5066  Assert(this->update_flags & update_hessians,
5067  ExcAccessToUninitializedField("update_hessians"));
5068  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5069  Assert(present_cell.get() != nullptr,
5070  ExcMessage("FEValues object is not reinit'ed to any cell"));
5071  // if the entire FE is primitive,
5072  // then we can take a short-cut:
5073  if (fe->is_primitive())
5074  return this->finite_element_output.shape_hessians[i][j];
5075  else
5076  {
5077  // otherwise, use the mapping
5078  // between shape function
5079  // numbers and rows. note that
5080  // by the assertions above, we
5081  // know that this particular
5082  // shape function is primitive,
5083  // so we can call
5084  // system_to_component_index
5085  const unsigned int row =
5086  this->finite_element_output
5087  .shape_function_to_row_table[i * fe->n_components() +
5088  fe->system_to_component_index(i).first];
5089  return this->finite_element_output.shape_hessians[row][j];
5090  }
5091 }
5092 
5093 
5094 
5095 template <int dim, int spacedim>
5096 inline Tensor<2, spacedim>
5098  const unsigned int i,
5099  const unsigned int j,
5100  const unsigned int component) const
5101 {
5102  AssertIndexRange(i, fe->dofs_per_cell);
5103  Assert(this->update_flags & update_hessians,
5104  ExcAccessToUninitializedField("update_hessians"));
5105  AssertIndexRange(component, fe->n_components());
5106  Assert(present_cell.get() != nullptr,
5107  ExcMessage("FEValues object is not reinit'ed to any cell"));
5108  // check whether the shape function
5109  // is non-zero at all within
5110  // this component:
5111  if (fe->get_nonzero_components(i)[component] == false)
5112  return Tensor<2, spacedim>();
5113 
5114  // look up the right row in the
5115  // table and take the data from
5116  // there
5117  const unsigned int row =
5118  this->finite_element_output
5119  .shape_function_to_row_table[i * fe->n_components() + component];
5120  return this->finite_element_output.shape_hessians[row][j];
5121 }
5122 
5123 
5124 
5125 template <int dim, int spacedim>
5126 inline const Tensor<3, spacedim> &
5128  const unsigned int j) const
5129 {
5130  AssertIndexRange(i, fe->dofs_per_cell);
5131  Assert(this->update_flags & update_3rd_derivatives,
5132  ExcAccessToUninitializedField("update_3rd_derivatives"));
5133  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5134  Assert(present_cell.get() != nullptr,
5135  ExcMessage("FEValues object is not reinit'ed to any cell"));
5136  // if the entire FE is primitive,
5137  // then we can take a short-cut:
5138  if (fe->is_primitive())
5139  return this->finite_element_output.shape_3rd_derivatives[i][j];
5140  else
5141  {
5142  // otherwise, use the mapping
5143  // between shape function
5144  // numbers and rows. note that
5145  // by the assertions above, we
5146  // know that this particular
5147  // shape function is primitive,
5148  // so we can call
5149  // system_to_component_index
5150  const unsigned int row =
5151  this->finite_element_output
5152  .shape_function_to_row_table[i * fe->n_components() +
5153  fe->system_to_component_index(i).first];
5154  return this->finite_element_output.shape_3rd_derivatives[row][j];
5155  }
5156 }
5157 
5158 
5159 
5160 template <int dim, int spacedim>
5161 inline Tensor<3, spacedim>
5163  const unsigned int i,
5164  const unsigned int j,
5165  const unsigned int component) const
5166 {
5167  AssertIndexRange(i, fe->dofs_per_cell);
5168  Assert(this->update_flags & update_3rd_derivatives,
5169  ExcAccessToUninitializedField("update_3rd_derivatives"));
5170  AssertIndexRange(component, fe->n_components());
5171  Assert(present_cell.get() != nullptr,
5172  ExcMessage("FEValues object is not reinit'ed to any cell"));
5173  // check whether the shape function
5174  // is non-zero at all within
5175  // this component:
5176  if (fe->get_nonzero_components(i)[component] == false)
5177  return Tensor<3, spacedim>();
5178 
5179  // look up the right row in the
5180  // table and take the data from
5181  // there
5182  const unsigned int row =
5183  this->finite_element_output
5184  .shape_function_to_row_table[i * fe->n_components() + component];
5185  return this->finite_element_output.shape_3rd_derivatives[row][j];
5186 }
5187 
5188 
5189 
5190 template <int dim, int spacedim>
5191 inline const FiniteElement<dim, spacedim> &
5193 {
5194  return *fe;
5195 }
5196 
5197 
5198 
5199 template <int dim, int spacedim>
5200 inline const Mapping<dim, spacedim> &
5202 {
5203  return *mapping;
5204 }
5205 
5206 
5207 
5208 template <int dim, int spacedim>
5209 inline UpdateFlags
5211 {
5212  return this->update_flags;
5213 }
5214 
5215 
5216 
5217 template <int dim, int spacedim>
5218 inline const std::vector<Point<spacedim>> &
5220 {
5221  Assert(this->update_flags & update_quadrature_points,
5222  ExcAccessToUninitializedField("update_quadrature_points"));
5223  Assert(present_cell.get() != nullptr,
5224  ExcMessage("FEValues object is not reinit'ed to any cell"));
5225  return this->mapping_output.quadrature_points;
5226 }
5227 
5228 
5229 
5230 template <int dim, int spacedim>
5231 inline const std::vector<double> &
5233 {
5234  Assert(this->update_flags & update_JxW_values,
5235  ExcAccessToUninitializedField("update_JxW_values"));
5236  Assert(present_cell.get() != nullptr,
5237  ExcMessage("FEValues object is not reinit'ed to any cell"));
5238  return this->mapping_output.JxW_values;
5239 }
5240 
5241 
5242 
5243 template <int dim, int spacedim>
5244 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5246 {
5247  Assert(this->update_flags & update_jacobians,
5248  ExcAccessToUninitializedField("update_jacobians"));
5249  Assert(present_cell.get() != nullptr,
5250  ExcMessage("FEValues object is not reinit'ed to any cell"));
5251  return this->mapping_output.jacobians;
5252 }
5253 
5254 
5255 
5256 template <int dim, int spacedim>
5257 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5259 {
5260  Assert(this->update_flags & update_jacobian_grads,
5261  ExcAccessToUninitializedField("update_jacobians_grads"));
5262  Assert(present_cell.get() != nullptr,
5263  ExcMessage("FEValues object is not reinit'ed to any cell"));
5264  return this->mapping_output.jacobian_grads;
5265 }
5266 
5267 
5268 
5269 template <int dim, int spacedim>
5270 inline const Tensor<3, spacedim> &
5272  const unsigned int i) const
5273 {
5274  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5275  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5276  Assert(present_cell.get() != nullptr,
5277  ExcMessage("FEValues object is not reinit'ed to any cell"));
5278  return this->mapping_output.jacobian_pushed_forward_grads[i];
5279 }
5280 
5281 
5282 
5283 template <int dim, int spacedim>
5284 inline const std::vector<Tensor<3, spacedim>> &
5286 {
5287  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5288  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5289  Assert(present_cell.get() != nullptr,
5290  ExcMessage("FEValues object is not reinit'ed to any cell"));
5291  return this->mapping_output.jacobian_pushed_forward_grads;
5292 }
5293 
5294 
5295 
5296 template <int dim, int spacedim>
5297 inline const DerivativeForm<3, dim, spacedim> &
5298 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5299 {
5300  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5301  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5302  Assert(present_cell.get() != nullptr,
5303  ExcMessage("FEValues object is not reinit'ed to any cell"));
5304  return this->mapping_output.jacobian_2nd_derivatives[i];
5305 }
5306 
5307 
5308 
5309 template <int dim, int spacedim>
5310 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5312 {
5313  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5314  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5315  Assert(present_cell.get() != nullptr,
5316  ExcMessage("FEValues object is not reinit'ed to any cell"));
5317  return this->mapping_output.jacobian_2nd_derivatives;
5318 }
5319 
5320 
5321 
5322 template <int dim, int spacedim>
5323 inline const Tensor<4, spacedim> &
5325  const unsigned int i) const
5326 {
5329  "update_jacobian_pushed_forward_2nd_derivatives"));
5330  Assert(present_cell.get() != nullptr,
5331  ExcMessage("FEValues object is not reinit'ed to any cell"));
5332  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5333 }
5334 
5335 
5336 
5337 template <int dim, int spacedim>
5338 inline const std::vector<Tensor<4, spacedim>> &
5340 {
5341  Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
5343  "update_jacobian_pushed_forward_2nd_derivatives"));
5344  Assert(present_cell.get() != nullptr,
5345  ExcMessage("FEValues object is not reinit'ed to any cell"));
5346  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5347 }
5348 
5349 
5350 
5351 template <int dim, int spacedim>
5352 inline const DerivativeForm<4, dim, spacedim> &
5353 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5354 {
5355  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5356  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5357  Assert(present_cell.get() != nullptr,
5358  ExcMessage("FEValues object is not reinit'ed to any cell"));
5359  return this->mapping_output.jacobian_3rd_derivatives[i];
5360 }
5361 
5362 
5363 
5364 template <int dim, int spacedim>
5365 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5367 {
5368  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5369  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5370  Assert(present_cell.get() != nullptr,
5371  ExcMessage("FEValues object is not reinit'ed to any cell"));
5372  return this->mapping_output.jacobian_3rd_derivatives;
5373 }
5374 
5375 
5376 
5377 template <int dim, int spacedim>
5378 inline const Tensor<5, spacedim> &
5380  const unsigned int i) const
5381 {
5384  "update_jacobian_pushed_forward_3rd_derivatives"));
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_3rd_derivatives[i];
5388 }
5389 
5390 
5391 
5392 template <int dim, int spacedim>
5393 inline const std::vector<Tensor<5, spacedim>> &
5395 {
5396  Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
5398  "update_jacobian_pushed_forward_3rd_derivatives"));
5399  Assert(present_cell.get() != nullptr,
5400  ExcMessage("FEValues object is not reinit'ed to any cell"));
5401  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
5402 }
5403 
5404 
5405 
5406 template <int dim, int spacedim>
5407 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
5409 {
5410  Assert(this->update_flags & update_inverse_jacobians,
5411  ExcAccessToUninitializedField("update_inverse_jacobians"));
5412  Assert(present_cell.get() != nullptr,
5413  ExcMessage("FEValues object is not reinit'ed to any cell"));
5414  return this->mapping_output.inverse_jacobians;
5415 }
5416 
5417 
5418 
5419 template <int dim, int spacedim>
5420 inline const Point<spacedim> &
5421 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
5422 {
5423  Assert(this->update_flags & update_quadrature_points,
5424  ExcAccessToUninitializedField("update_quadrature_points"));
5425  AssertIndexRange(i, this->mapping_output.quadrature_points.size());
5426  Assert(present_cell.get() != nullptr,
5427  ExcMessage("FEValues object is not reinit'ed to any cell"));
5428 
5429  return this->mapping_output.quadrature_points[i];
5430 }
5431 
5432 
5433 
5434 template <int dim, int spacedim>
5435 inline double
5436 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
5437 {
5438  Assert(this->update_flags & update_JxW_values,
5439  ExcAccessToUninitializedField("update_JxW_values"));
5440  AssertIndexRange(i, this->mapping_output.JxW_values.size());
5441  Assert(present_cell.get() != nullptr,
5442  ExcMessage("FEValues object is not reinit'ed to any cell"));
5443 
5444  return this->mapping_output.JxW_values[i];
5445 }
5446 
5447 
5448 
5449 template <int dim, int spacedim>
5450 inline const DerivativeForm<1, dim, spacedim> &
5451 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
5452 {
5453  Assert(this->update_flags & update_jacobians,
5454  ExcAccessToUninitializedField("update_jacobians"));
5455  AssertIndexRange(i, this->mapping_output.jacobians.size());
5456  Assert(present_cell.get() != nullptr,
5457  ExcMessage("FEValues object is not reinit'ed to any cell"));
5458 
5459  return this->mapping_output.jacobians[i];
5460 }
5461 
5462 
5463 
5464 template <int dim, int spacedim>
5465 inline const DerivativeForm<2, dim, spacedim> &
5466 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
5467 {
5468  Assert(this->update_flags & update_jacobian_grads,
5469  ExcAccessToUninitializedField("update_jacobians_grads"));
5470  AssertIndexRange(i, this->mapping_output.jacobian_grads.size());
5471  Assert(present_cell.get() != nullptr,
5472  ExcMessage("FEValues object is not reinit'ed to any cell"));
5473 
5474  return this->mapping_output.jacobian_grads[i];
5475 }
5476 
5477 
5478 
5479 template <int dim, int spacedim>
5480 inline const DerivativeForm<1, spacedim, dim> &
5481 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
5482 {
5483  Assert(this->update_flags & update_inverse_jacobians,
5484  ExcAccessToUninitializedField("update_inverse_jacobians"));
5485  AssertIndexRange(i, this->mapping_output.inverse_jacobians.size());
5486  Assert(present_cell.get() != nullptr,
5487  ExcMessage("FEValues object is not reinit'ed to any cell"));
5488 
5489  return this->mapping_output.inverse_jacobians[i];
5490 }
5491 
5492 
5493 
5494 template <int dim, int spacedim>
5495 inline const Tensor<1, spacedim> &
5496 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
5497 {
5498  Assert(this->update_flags & update_normal_vectors,
5500  "update_normal_vectors")));
5501  AssertIndexRange(i, this->mapping_output.normal_vectors.size());
5502  Assert(present_cell.get() != nullptr,
5503  ExcMessage("FEValues object is not reinit'ed to any cell"));
5504 
5505  return this->mapping_output.normal_vectors[i];
5506 }
5507 
5508 
5509 
5510 /*--------------------- Inline functions: FEValues --------------------------*/
5511 
5512 
5513 template <int dim, int spacedim>
5514 inline const Quadrature<dim> &
5516 {
5517  return quadrature;
5518 }
5519 
5520 
5521 
5522 template <int dim, int spacedim>
5523 inline const FEValues<dim, spacedim> &
5525 {
5526  return *this;
5527 }
5528 
5529 
5530 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
5531 
5532 
5533 template <int dim, int spacedim>
5534 inline unsigned int
5536 {
5537  return present_face_index;
5538 }
5539 
5540 
5541 /*----------------------- Inline functions: FE*FaceValues -------------------*/
5542 
5543 template <int dim, int spacedim>
5544 inline const Quadrature<dim - 1> &
5546 {
5547  return quadrature;
5548 }
5549 
5550 
5551 
5552 template <int dim, int spacedim>
5553 inline const FEFaceValues<dim, spacedim> &
5555 {
5556  return *this;
5557 }
5558 
5559 
5560 
5561 template <int dim, int spacedim>
5562 inline const FESubfaceValues<dim, spacedim> &
5564 {
5565  return *this;
5566 }
5567 
5568 
5569 
5570 template <int dim, int spacedim>
5571 inline const Tensor<1, spacedim> &
5572 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
5573 {
5574  AssertIndexRange(i, this->mapping_output.boundary_forms.size());
5575  Assert(this->update_flags & update_boundary_forms,
5577  "update_boundary_forms")));
5578 
5579  return this->mapping_output.boundary_forms[i];
5580 }
5581 
5582 #endif // DOXYGEN
5583 
5584 DEAL_II_NAMESPACE_CLOSE
5585 
5586 #endif
Transformed quadrature weights.
constexpr Tensor()=default
Shape function values.
typename ProductType< Number, typename Vector< dim, spacedim >::curl_type >::type curl_type
Definition: fe_values.h:696
const FEFaceValues< dim, spacedim > & get_present_fe_values() 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:3429
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:3461
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:1597
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1294
unsigned int present_face_index
Definition: fe_values.h:3694
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1520
static ::ExceptionBase & ExcAccessToUninitializedField()
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:536
const unsigned int dofs_per_cell
Definition: fe_values.h:2107
typename ::internal::CurlType< spacedim >::type curl_type
Definition: fe_values.h:627
const unsigned int component
Definition: fe_values.h:542
const Quadrature< dim - 1 > quadrature
Definition: fe_values.h:3699
Volume element.
#define AssertIndexRange(index, range)
Definition: exceptions.h:1641
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:213
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:3345
UpdateFlags get_update_flags() const
Transformed quadrature points.
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:189
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:3397
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:3476
LinearAlgebra::distributed::Vector< Number > Vector
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
typename ProductType< Number, typename Vector< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:664
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:197
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1581
static ::ExceptionBase & ExcMessage(std::string arg1)
constexpr SymmetricTensor()=default
Number value_type
Definition: vector.h:125
typename ProductType< Number, typename Vector< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:712
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:518
const double & shape_value(const unsigned int function_no, const unsigned int point_no) 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:1867
unsigned int get_face_index() const
#define Assert(cond, exc)
Definition: exceptions.h:1411
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:302
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1227
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
const Quadrature< dim > quadrature
Definition: fe_values.h:3592
const unsigned int first_vector_component
Definition: fe_values.h:1222
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
#define DeclException0(Exception0)
Definition: exceptions.h:473
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition: fe_values.h:3405
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:452
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:1302
std::vector<::FEValuesViews::Scalar< dim, spacedim > > scalars
Definition: fe_values.h:1954
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:2100
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:3370
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:656
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
Definition: tensor.h:422
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:672
typename ProductType< Number, typename Vector< dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:680
typename ProductType< Number, typename Vector< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:704
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:1589
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:1856
Definition: fe.h:38
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1509
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
static DEAL_II_CONSTEXPR TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
Definition: tensor.h:1528
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1216
DEAL_II_CONSTEXPR SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
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:3361
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3437
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:688
::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:3412
const std::vector< Tensor< 5, spacedim > > & get_jacobian_pushed_forward_3rd_derivatives() const
const FESubfaceValues< dim, spacedim > & get_present_fe_values() const
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:205
typename ProductType< Number, typename Scalar< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:221
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:547
UpdateFlags update_flags
Definition: fe_values.h:3443
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:3421
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