Reference documentation for deal.II version Git f81eda9982 2020-03-28 21:30:57 -0400
\(\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 
3142  const std::vector<Tensor<1, spacedim>> &
3143  get_normal_vectors() const;
3144 
3146 
3148 
3149 
3159  operator[](const FEValuesExtractors::Scalar &scalar) const;
3160 
3170  operator[](const FEValuesExtractors::Vector &vector) const;
3171 
3182  operator[](const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
3183 
3184 
3194  operator[](const FEValuesExtractors::Tensor<2> &tensor) const;
3195 
3197 
3199 
3200 
3204  const Mapping<dim, spacedim> &
3205  get_mapping() const;
3206 
3211  get_fe() const;
3212 
3216  UpdateFlags
3217  get_update_flags() const;
3218 
3223  get_cell() const;
3224 
3231  get_cell_similarity() const;
3232 
3237  std::size_t
3238  memory_consumption() const;
3240 
3241 
3250  std::string,
3251  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3252  << "object for which this kind of information has not been computed. What "
3253  << "information these objects compute is determined by the update_* flags you "
3254  << "pass to the constructor. Here, the operation you are attempting requires "
3255  << "the <" << arg1
3256  << "> flag to be set, but it was apparently not specified "
3257  << "upon construction.");
3258 
3266  ExcFEDontMatch,
3267  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3268  "to the DoFHandler that provided the cell iterator do not match.");
3274  DeclException1(ExcShapeFunctionNotPrimitive,
3275  int,
3276  << "The shape function with index " << arg1
3277  << " is not primitive, i.e. it is vector-valued and "
3278  << "has more than one non-zero vector component. This "
3279  << "function cannot be called for these shape functions. "
3280  << "Maybe you want to use the same function with the "
3281  << "_component suffix?");
3282 
3291  "The given FiniteElement is not a primitive element but the requested operation "
3292  "only works for those. See FiniteElement::is_primitive() for more information.");
3293 
3294 protected:
3327  class CellIteratorBase;
3328 
3333  template <typename CI>
3335  class TriaCellIterator;
3336 
3342  std::unique_ptr<const CellIteratorBase> present_cell;
3343 
3351  boost::signals2::connection tria_listener_refinement;
3352 
3360  boost::signals2::connection tria_listener_mesh_transform;
3361 
3367  void
3368  invalidate_present_cell();
3369 
3379  void
3380  maybe_invalidate_previous_present_cell(
3381  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3382 
3388 
3394  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3396 
3403 
3404 
3412 
3418  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3420 
3426  spacedim>
3428 
3429 
3434 
3443  UpdateFlags
3444  compute_update_flags(const UpdateFlags update_flags) const;
3445 
3452 
3458  void
3459  check_cell_similarity(
3460  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3461 
3462 private:
3467 
3468  // Make the view classes friends of this class, since they access internal
3469  // data.
3470  template <int, int>
3471  friend class FEValuesViews::Scalar;
3472  template <int, int>
3473  friend class FEValuesViews::Vector;
3474  template <int, int, int>
3475  friend class FEValuesViews::SymmetricTensor;
3476  template <int, int, int>
3477  friend class FEValuesViews::Tensor;
3478 };
3479 
3480 
3481 
3492 template <int dim, int spacedim = dim>
3493 class FEValues : public FEValuesBase<dim, spacedim>
3494 {
3495 public:
3500  static const unsigned int integral_dimension = dim;
3501 
3506  FEValues(const Mapping<dim, spacedim> & mapping,
3507  const FiniteElement<dim, spacedim> &fe,
3508  const Quadrature<dim> & quadrature,
3509  const UpdateFlags update_flags);
3510 
3517  const Quadrature<dim> & quadrature,
3518  const UpdateFlags update_flags);
3519 
3526  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3527  void
3528  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3529  level_dof_access>> &cell);
3530 
3544  void
3545  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3546 
3551  const Quadrature<dim> &
3552  get_quadrature() const;
3553 
3558  std::size_t
3559  memory_consumption() const;
3560 
3575  const FEValues<dim, spacedim> &
3576  get_present_fe_values() const;
3577 
3578 private:
3583 
3587  void
3588  initialize(const UpdateFlags update_flags);
3589 
3596  void
3597  do_reinit();
3598 };
3599 
3600 
3611 template <int dim, int spacedim = dim>
3612 class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
3613 {
3614 public:
3619  static const unsigned int integral_dimension = dim - 1;
3620 
3632  FEFaceValuesBase(const unsigned int n_q_points,
3633  const unsigned int dofs_per_cell,
3634  const UpdateFlags update_flags,
3635  const Mapping<dim, spacedim> & mapping,
3636  const FiniteElement<dim, spacedim> &fe,
3637  const Quadrature<dim - 1> & quadrature);
3638 
3646  const Tensor<1, spacedim> &
3647  boundary_form(const unsigned int i) const;
3648 
3655  const std::vector<Tensor<1, spacedim>> &
3656  get_boundary_forms() const;
3657 
3662  unsigned int
3663  get_face_index() const;
3664 
3669  const Quadrature<dim - 1> &
3670  get_quadrature() const;
3671 
3676  std::size_t
3677  memory_consumption() const;
3678 
3679 protected:
3684  unsigned int present_face_index;
3685 
3689  const Quadrature<dim - 1> quadrature;
3690 };
3691 
3692 
3693 
3708 template <int dim, int spacedim = dim>
3709 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
3710 {
3711 public:
3716  static const unsigned int dimension = dim;
3717 
3718  static const unsigned int space_dimension = spacedim;
3719 
3724  static const unsigned int integral_dimension = dim - 1;
3725 
3730  FEFaceValues(const Mapping<dim, spacedim> & mapping,
3731  const FiniteElement<dim, spacedim> &fe,
3732  const Quadrature<dim - 1> & quadrature,
3733  const UpdateFlags update_flags);
3734 
3741  const Quadrature<dim - 1> & quadrature,
3742  const UpdateFlags update_flags);
3743 
3748  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3749  void
3750  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3751  level_dof_access>> &cell,
3752  const unsigned int face_no);
3753 
3760  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3761  void
3762  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3763  level_dof_access>> & cell,
3764  const typename Triangulation<dim, spacedim>::face_iterator &face);
3765 
3779  void
3780  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3781  const unsigned int face_no);
3782 
3783  /*
3784  * Reinitialize the gradients, Jacobi determinants, etc for the given face
3785  * on a given cell of type "iterator into a Triangulation object", and the
3786  * given finite element. Since iterators into a triangulation alone only
3787  * convey information about the geometry of a cell, but not about degrees of
3788  * freedom possibly associated with this cell, you will not be able to call
3789  * some functions of this class if they need information about degrees of
3790  * freedom. These functions are, above all, the
3791  * <tt>get_function_value/gradients/hessians/third_derivatives</tt>
3792  * functions. If you want to call these functions, you have to call the @p
3793  * reinit variants that take iterators into DoFHandler or other DoF handler
3794  * type objects.
3795  *
3796  * @note @p face must be one of @p cell's face iterators.
3797  */
3798  void
3799  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3800  const typename Triangulation<dim, spacedim>::face_iterator &face);
3801 
3817  get_present_fe_values() const;
3818 
3819 private:
3823  void
3824  initialize(const UpdateFlags update_flags);
3825 
3832  void
3833  do_reinit(const unsigned int face_no);
3834 };
3835 
3836 
3854 template <int dim, int spacedim = dim>
3855 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
3856 {
3857 public:
3861  static const unsigned int dimension = dim;
3862 
3866  static const unsigned int space_dimension = spacedim;
3867 
3872  static const unsigned int integral_dimension = dim - 1;
3873 
3878  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
3879  const FiniteElement<dim, spacedim> &fe,
3880  const Quadrature<dim - 1> & face_quadrature,
3881  const UpdateFlags update_flags);
3882 
3889  const Quadrature<dim - 1> & face_quadrature,
3890  const UpdateFlags update_flags);
3891 
3898  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3899  void
3900  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3901  level_dof_access>> &cell,
3902  const unsigned int face_no,
3903  const unsigned int subface_no);
3904 
3909  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3910  void
3911  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3912  level_dof_access>> & cell,
3913  const typename Triangulation<dim, spacedim>::face_iterator &face,
3914  const typename Triangulation<dim, spacedim>::face_iterator &subface);
3915 
3929  void
3930  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3931  const unsigned int face_no,
3932  const unsigned int subface_no);
3933 
3953  void
3954  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3955  const typename Triangulation<dim, spacedim>::face_iterator &face,
3956  const typename Triangulation<dim, spacedim>::face_iterator &subface);
3957 
3973  get_present_fe_values() const;
3974 
3980  DeclException0(ExcReinitCalledWithBoundaryFace);
3981 
3987  DeclException0(ExcFaceHasNoSubfaces);
3988 
3989 private:
3993  void
3994  initialize(const UpdateFlags update_flags);
3995 
4002  void
4003  do_reinit(const unsigned int face_no, const unsigned int subface_no);
4004 };
4005 
4006 
4007 #ifndef DOXYGEN
4008 
4009 
4010 /*------------------------ Inline functions: namespace FEValuesViews --------*/
4011 
4012 namespace FEValuesViews
4013 {
4014  template <int dim, int spacedim>
4015  inline typename Scalar<dim, spacedim>::value_type
4016  Scalar<dim, spacedim>::value(const unsigned int shape_function,
4017  const unsigned int q_point) const
4018  {
4019  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4020  Assert(
4021  fe_values->update_flags & update_values,
4023  "update_values"))));
4024 
4025  // an adaptation of the FEValuesBase::shape_value_component function
4026  // except that here we know the component as fixed and we have
4027  // pre-computed and cached a bunch of information. See the comments there.
4028  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4029  return fe_values->finite_element_output.shape_values(
4030  shape_function_data[shape_function].row_index, q_point);
4031  else
4032  return 0;
4033  }
4034 
4035 
4036 
4037  template <int dim, int spacedim>
4038  inline typename Scalar<dim, spacedim>::gradient_type
4039  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
4040  const unsigned int q_point) const
4041  {
4042  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4043  Assert(fe_values->update_flags & update_gradients,
4045  "update_gradients")));
4046 
4047  // an adaptation of the FEValuesBase::shape_grad_component
4048  // function except that here we know the component as fixed and we have
4049  // pre-computed and cached a bunch of information. See the comments there.
4050  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4051  return fe_values->finite_element_output
4052  .shape_gradients[shape_function_data[shape_function].row_index]
4053  [q_point];
4054  else
4055  return gradient_type();
4056  }
4057 
4058 
4059 
4060  template <int dim, int spacedim>
4061  inline typename Scalar<dim, spacedim>::hessian_type
4062  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
4063  const unsigned int q_point) const
4064  {
4065  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4066  Assert(fe_values->update_flags & update_hessians,
4068  "update_hessians")));
4069 
4070  // an adaptation of the FEValuesBase::shape_hessian_component
4071  // function except that here we know the component as fixed and we have
4072  // pre-computed and cached a bunch of information. See the comments there.
4073  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4074  return fe_values->finite_element_output
4075  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
4076  else
4077  return hessian_type();
4078  }
4079 
4080 
4081 
4082  template <int dim, int spacedim>
4083  inline typename Scalar<dim, spacedim>::third_derivative_type
4084  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
4085  const unsigned int q_point) const
4086  {
4087  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4088  Assert(fe_values->update_flags & update_3rd_derivatives,
4090  "update_3rd_derivatives")));
4091 
4092  // an adaptation of the FEValuesBase::shape_3rdderivative_component
4093  // function except that here we know the component as fixed and we have
4094  // pre-computed and cached a bunch of information. See the comments there.
4095  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4096  return fe_values->finite_element_output
4097  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
4098  [q_point];
4099  else
4100  return third_derivative_type();
4101  }
4102 
4103 
4104 
4105  template <int dim, int spacedim>
4106  inline typename Vector<dim, spacedim>::value_type
4107  Vector<dim, spacedim>::value(const unsigned int shape_function,
4108  const unsigned int q_point) const
4109  {
4110  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4111  Assert(fe_values->update_flags & update_values,
4113  "update_values")));
4114 
4115  // same as for the scalar case except that we have one more index
4116  const int snc =
4117  shape_function_data[shape_function].single_nonzero_component;
4118  if (snc == -2)
4119  return value_type();
4120  else if (snc != -1)
4121  {
4122  value_type return_value;
4123  return_value[shape_function_data[shape_function]
4124  .single_nonzero_component_index] =
4125  fe_values->finite_element_output.shape_values(snc, q_point);
4126  return return_value;
4127  }
4128  else
4129  {
4130  value_type return_value;
4131  for (unsigned int d = 0; d < dim; ++d)
4132  if (shape_function_data[shape_function]
4133  .is_nonzero_shape_function_component[d])
4134  return_value[d] = fe_values->finite_element_output.shape_values(
4135  shape_function_data[shape_function].row_index[d], q_point);
4136 
4137  return return_value;
4138  }
4139  }
4140 
4141 
4142 
4143  template <int dim, int spacedim>
4144  inline typename Vector<dim, spacedim>::gradient_type
4145  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4146  const unsigned int q_point) const
4147  {
4148  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4149  Assert(fe_values->update_flags & update_gradients,
4151  "update_gradients")));
4152 
4153  // same as for the scalar case except that we have one more index
4154  const int snc =
4155  shape_function_data[shape_function].single_nonzero_component;
4156  if (snc == -2)
4157  return gradient_type();
4158  else if (snc != -1)
4159  {
4160  gradient_type return_value;
4161  return_value[shape_function_data[shape_function]
4162  .single_nonzero_component_index] =
4163  fe_values->finite_element_output.shape_gradients[snc][q_point];
4164  return return_value;
4165  }
4166  else
4167  {
4168  gradient_type return_value;
4169  for (unsigned int d = 0; d < dim; ++d)
4170  if (shape_function_data[shape_function]
4171  .is_nonzero_shape_function_component[d])
4172  return_value[d] =
4173  fe_values->finite_element_output.shape_gradients
4174  [shape_function_data[shape_function].row_index[d]][q_point];
4175 
4176  return return_value;
4177  }
4178  }
4179 
4180 
4181 
4182  template <int dim, int spacedim>
4183  inline typename Vector<dim, spacedim>::divergence_type
4184  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4185  const unsigned int q_point) const
4186  {
4187  // this function works like in the case above
4188  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4189  Assert(fe_values->update_flags & update_gradients,
4191  "update_gradients")));
4192 
4193  // same as for the scalar case except that we have one more index
4194  const int snc =
4195  shape_function_data[shape_function].single_nonzero_component;
4196  if (snc == -2)
4197  return divergence_type();
4198  else if (snc != -1)
4199  return fe_values->finite_element_output
4200  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4201  .single_nonzero_component_index];
4202  else
4203  {
4204  divergence_type return_value = 0;
4205  for (unsigned int d = 0; d < dim; ++d)
4206  if (shape_function_data[shape_function]
4207  .is_nonzero_shape_function_component[d])
4208  return_value +=
4209  fe_values->finite_element_output.shape_gradients
4210  [shape_function_data[shape_function].row_index[d]][q_point][d];
4211 
4212  return return_value;
4213  }
4214  }
4215 
4216 
4217 
4218  template <int dim, int spacedim>
4219  inline typename Vector<dim, spacedim>::curl_type
4220  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4221  const unsigned int q_point) const
4222  {
4223  // this function works like in the case above
4224 
4225  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4226  Assert(fe_values->update_flags & update_gradients,
4228  "update_gradients")));
4229  // same as for the scalar case except that we have one more index
4230  const int snc =
4231  shape_function_data[shape_function].single_nonzero_component;
4232 
4233  if (snc == -2)
4234  return curl_type();
4235 
4236  else
4237  switch (dim)
4238  {
4239  case 1:
4240  {
4241  Assert(false,
4242  ExcMessage(
4243  "Computing the curl in 1d is not a useful operation"));
4244  return curl_type();
4245  }
4246 
4247  case 2:
4248  {
4249  if (snc != -1)
4250  {
4251  curl_type return_value;
4252 
4253  // the single nonzero component can only be zero or one in 2d
4254  if (shape_function_data[shape_function]
4255  .single_nonzero_component_index == 0)
4256  return_value[0] =
4257  -1.0 * fe_values->finite_element_output
4258  .shape_gradients[snc][q_point][1];
4259  else
4260  return_value[0] = fe_values->finite_element_output
4261  .shape_gradients[snc][q_point][0];
4262 
4263  return return_value;
4264  }
4265 
4266  else
4267  {
4268  curl_type return_value;
4269 
4270  return_value[0] = 0.0;
4271 
4272  if (shape_function_data[shape_function]
4273  .is_nonzero_shape_function_component[0])
4274  return_value[0] -=
4275  fe_values->finite_element_output
4276  .shape_gradients[shape_function_data[shape_function]
4277  .row_index[0]][q_point][1];
4278 
4279  if (shape_function_data[shape_function]
4280  .is_nonzero_shape_function_component[1])
4281  return_value[0] +=
4282  fe_values->finite_element_output
4283  .shape_gradients[shape_function_data[shape_function]
4284  .row_index[1]][q_point][0];
4285 
4286  return return_value;
4287  }
4288  }
4289 
4290  case 3:
4291  {
4292  if (snc != -1)
4293  {
4294  curl_type return_value;
4295 
4296  switch (shape_function_data[shape_function]
4297  .single_nonzero_component_index)
4298  {
4299  case 0:
4300  {
4301  return_value[0] = 0;
4302  return_value[1] = fe_values->finite_element_output
4303  .shape_gradients[snc][q_point][2];
4304  return_value[2] =
4305  -1.0 * fe_values->finite_element_output
4306  .shape_gradients[snc][q_point][1];
4307  return return_value;
4308  }
4309 
4310  case 1:
4311  {
4312  return_value[0] =
4313  -1.0 * fe_values->finite_element_output
4314  .shape_gradients[snc][q_point][2];
4315  return_value[1] = 0;
4316  return_value[2] = fe_values->finite_element_output
4317  .shape_gradients[snc][q_point][0];
4318  return return_value;
4319  }
4320 
4321  default:
4322  {
4323  return_value[0] = fe_values->finite_element_output
4324  .shape_gradients[snc][q_point][1];
4325  return_value[1] =
4326  -1.0 * fe_values->finite_element_output
4327  .shape_gradients[snc][q_point][0];
4328  return_value[2] = 0;
4329  return return_value;
4330  }
4331  }
4332  }
4333 
4334  else
4335  {
4336  curl_type return_value;
4337 
4338  for (unsigned int i = 0; i < dim; ++i)
4339  return_value[i] = 0.0;
4340 
4341  if (shape_function_data[shape_function]
4342  .is_nonzero_shape_function_component[0])
4343  {
4344  return_value[1] +=
4345  fe_values->finite_element_output
4346  .shape_gradients[shape_function_data[shape_function]
4347  .row_index[0]][q_point][2];
4348  return_value[2] -=
4349  fe_values->finite_element_output
4350  .shape_gradients[shape_function_data[shape_function]
4351  .row_index[0]][q_point][1];
4352  }
4353 
4354  if (shape_function_data[shape_function]
4355  .is_nonzero_shape_function_component[1])
4356  {
4357  return_value[0] -=
4358  fe_values->finite_element_output
4359  .shape_gradients[shape_function_data[shape_function]
4360  .row_index[1]][q_point][2];
4361  return_value[2] +=
4362  fe_values->finite_element_output
4363  .shape_gradients[shape_function_data[shape_function]
4364  .row_index[1]][q_point][0];
4365  }
4366 
4367  if (shape_function_data[shape_function]
4368  .is_nonzero_shape_function_component[2])
4369  {
4370  return_value[0] +=
4371  fe_values->finite_element_output
4372  .shape_gradients[shape_function_data[shape_function]
4373  .row_index[2]][q_point][1];
4374  return_value[1] -=
4375  fe_values->finite_element_output
4376  .shape_gradients[shape_function_data[shape_function]
4377  .row_index[2]][q_point][0];
4378  }
4379 
4380  return return_value;
4381  }
4382  }
4383  }
4384  // should not end up here
4385  Assert(false, ExcInternalError());
4386  return curl_type();
4387  }
4388 
4389 
4390 
4391  template <int dim, int spacedim>
4392  inline typename Vector<dim, spacedim>::hessian_type
4393  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
4394  const unsigned int q_point) const
4395  {
4396  // this function works like in the case above
4397  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4398  Assert(fe_values->update_flags & update_hessians,
4400  "update_hessians")));
4401 
4402  // same as for the scalar case except that we have one more index
4403  const int snc =
4404  shape_function_data[shape_function].single_nonzero_component;
4405  if (snc == -2)
4406  return hessian_type();
4407  else if (snc != -1)
4408  {
4409  hessian_type return_value;
4410  return_value[shape_function_data[shape_function]
4411  .single_nonzero_component_index] =
4412  fe_values->finite_element_output.shape_hessians[snc][q_point];
4413  return return_value;
4414  }
4415  else
4416  {
4417  hessian_type return_value;
4418  for (unsigned int d = 0; d < dim; ++d)
4419  if (shape_function_data[shape_function]
4420  .is_nonzero_shape_function_component[d])
4421  return_value[d] =
4422  fe_values->finite_element_output.shape_hessians
4423  [shape_function_data[shape_function].row_index[d]][q_point];
4424 
4425  return return_value;
4426  }
4427  }
4428 
4429 
4430 
4431  template <int dim, int spacedim>
4432  inline typename Vector<dim, spacedim>::third_derivative_type
4433  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
4434  const unsigned int q_point) const
4435  {
4436  // this function works like in the case above
4437  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4438  Assert(fe_values->update_flags & update_3rd_derivatives,
4440  "update_3rd_derivatives")));
4441 
4442  // same as for the scalar case except that we have one more index
4443  const int snc =
4444  shape_function_data[shape_function].single_nonzero_component;
4445  if (snc == -2)
4446  return third_derivative_type();
4447  else if (snc != -1)
4448  {
4449  third_derivative_type return_value;
4450  return_value[shape_function_data[shape_function]
4451  .single_nonzero_component_index] =
4452  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
4453  return return_value;
4454  }
4455  else
4456  {
4457  third_derivative_type return_value;
4458  for (unsigned int d = 0; d < dim; ++d)
4459  if (shape_function_data[shape_function]
4460  .is_nonzero_shape_function_component[d])
4461  return_value[d] =
4462  fe_values->finite_element_output.shape_3rd_derivatives
4463  [shape_function_data[shape_function].row_index[d]][q_point];
4464 
4465  return return_value;
4466  }
4467  }
4468 
4469 
4470 
4471  namespace internal
4472  {
4477  inline ::SymmetricTensor<2, 1>
4478  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
4479  {
4480  AssertIndexRange(n, 1);
4481  (void)n;
4482 
4483  return {{t[0]}};
4484  }
4485 
4486 
4487 
4488  inline ::SymmetricTensor<2, 2>
4489  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
4490  {
4491  switch (n)
4492  {
4493  case 0:
4494  {
4495  return {{t[0], 0, t[1] / 2}};
4496  }
4497  case 1:
4498  {
4499  return {{0, t[1], t[0] / 2}};
4500  }
4501  default:
4502  {
4503  AssertIndexRange(n, 2);
4504  return {};
4505  }
4506  }
4507  }
4508 
4509 
4510 
4511  inline ::SymmetricTensor<2, 3>
4512  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
4513  {
4514  switch (n)
4515  {
4516  case 0:
4517  {
4518  return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
4519  }
4520  case 1:
4521  {
4522  return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
4523  }
4524  case 2:
4525  {
4526  return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
4527  }
4528  default:
4529  {
4530  AssertIndexRange(n, 3);
4531  return {};
4532  }
4533  }
4534  }
4535  } // namespace internal
4536 
4537 
4538 
4539  template <int dim, int spacedim>
4540  inline typename Vector<dim, spacedim>::symmetric_gradient_type
4541  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
4542  const unsigned int q_point) const
4543  {
4544  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4545  Assert(fe_values->update_flags & update_gradients,
4547  "update_gradients")));
4548 
4549  // same as for the scalar case except that we have one more index
4550  const int snc =
4551  shape_function_data[shape_function].single_nonzero_component;
4552  if (snc == -2)
4553  return symmetric_gradient_type();
4554  else if (snc != -1)
4555  return internal::symmetrize_single_row(
4556  shape_function_data[shape_function].single_nonzero_component_index,
4557  fe_values->finite_element_output.shape_gradients[snc][q_point]);
4558  else
4559  {
4560  gradient_type return_value;
4561  for (unsigned int d = 0; d < dim; ++d)
4562  if (shape_function_data[shape_function]
4563  .is_nonzero_shape_function_component[d])
4564  return_value[d] =
4565  fe_values->finite_element_output.shape_gradients
4566  [shape_function_data[shape_function].row_index[d]][q_point];
4567 
4568  return symmetrize(return_value);
4569  }
4570  }
4571 
4572 
4573 
4574  template <int dim, int spacedim>
4576  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
4577  const unsigned int q_point) const
4578  {
4579  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4580  Assert(fe_values->update_flags & update_values,
4582  "update_values")));
4583 
4584  // similar to the vector case where we have more then one index and we need
4585  // to convert between unrolled and component indexing for tensors
4586  const int snc =
4587  shape_function_data[shape_function].single_nonzero_component;
4588 
4589  if (snc == -2)
4590  {
4591  // shape function is zero for the selected components
4592  return value_type();
4593  }
4594  else if (snc != -1)
4595  {
4596  value_type return_value;
4597  const unsigned int comp =
4598  shape_function_data[shape_function].single_nonzero_component_index;
4599  return_value[value_type::unrolled_to_component_indices(comp)] =
4600  fe_values->finite_element_output.shape_values(snc, q_point);
4601  return return_value;
4602  }
4603  else
4604  {
4605  value_type return_value;
4606  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
4607  if (shape_function_data[shape_function]
4608  .is_nonzero_shape_function_component[d])
4609  return_value[value_type::unrolled_to_component_indices(d)] =
4610  fe_values->finite_element_output.shape_values(
4611  shape_function_data[shape_function].row_index[d], q_point);
4612  return return_value;
4613  }
4614  }
4615 
4616 
4617 
4618  template <int dim, int spacedim>
4621  const unsigned int shape_function,
4622  const unsigned int q_point) const
4623  {
4624  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4625  Assert(fe_values->update_flags & update_gradients,
4627  "update_gradients")));
4628 
4629  const int snc =
4630  shape_function_data[shape_function].single_nonzero_component;
4631 
4632  if (snc == -2)
4633  {
4634  // shape function is zero for the selected components
4635  return divergence_type();
4636  }
4637  else if (snc != -1)
4638  {
4639  // we have a single non-zero component when the symmetric tensor is
4640  // represented in unrolled form. this implies we potentially have
4641  // two non-zero components when represented in component form! we
4642  // will only have one non-zero entry if the non-zero component lies on
4643  // the diagonal of the tensor.
4644  //
4645  // the divergence of a second-order tensor is a first order tensor.
4646  //
4647  // assume the second-order tensor is A with components A_{ij}. then
4648  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
4649  // entries in the tensorial representation. define the
4650  // divergence as:
4651  // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
4652  // (which is incidentally also
4653  // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
4654  // In both cases, a sum is implied.
4655  //
4656  // Now, we know the nonzero component in unrolled form: it is indicated
4657  // by 'snc'. we can figure out which tensor components belong to this:
4658  const unsigned int comp =
4659  shape_function_data[shape_function].single_nonzero_component_index;
4660  const unsigned int ii =
4661  value_type::unrolled_to_component_indices(comp)[0];
4662  const unsigned int jj =
4663  value_type::unrolled_to_component_indices(comp)[1];
4664 
4665  // given the form of the divergence above, if ii=jj there is only a
4666  // single nonzero component of the full tensor and the gradient
4667  // equals
4668  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
4669  // all other entries of 'b' are zero
4670  //
4671  // on the other hand, if ii!=jj, then there are two nonzero entries in
4672  // the full tensor and
4673  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
4674  // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
4675  // again, all other entries of 'b' are zero
4676  const ::Tensor<1, spacedim> &phi_grad =
4677  fe_values->finite_element_output.shape_gradients[snc][q_point];
4678 
4679  divergence_type return_value;
4680  return_value[ii] = phi_grad[jj];
4681 
4682  if (ii != jj)
4683  return_value[jj] = phi_grad[ii];
4684 
4685  return return_value;
4686  }
4687  else
4688  {
4689  Assert(false, ExcNotImplemented());
4690  divergence_type return_value;
4691  return return_value;
4692  }
4693  }
4694 
4695 
4696 
4697  template <int dim, int spacedim>
4698  inline typename Tensor<2, dim, spacedim>::value_type
4699  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
4700  const unsigned int q_point) const
4701  {
4702  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4703  Assert(fe_values->update_flags & update_values,
4705  "update_values")));
4706 
4707  // similar to the vector case where we have more then one index and we need
4708  // to convert between unrolled and component indexing for tensors
4709  const int snc =
4710  shape_function_data[shape_function].single_nonzero_component;
4711 
4712  if (snc == -2)
4713  {
4714  // shape function is zero for the selected components
4715  return value_type();
4716  }
4717  else if (snc != -1)
4718  {
4719  value_type return_value;
4720  const unsigned int comp =
4721  shape_function_data[shape_function].single_nonzero_component_index;
4722  const TableIndices<2> indices =
4724  return_value[indices] =
4725  fe_values->finite_element_output.shape_values(snc, q_point);
4726  return return_value;
4727  }
4728  else
4729  {
4730  value_type return_value;
4731  for (unsigned int d = 0; d < dim * dim; ++d)
4732  if (shape_function_data[shape_function]
4733  .is_nonzero_shape_function_component[d])
4734  {
4735  const TableIndices<2> indices =
4737  return_value[indices] =
4738  fe_values->finite_element_output.shape_values(
4739  shape_function_data[shape_function].row_index[d], q_point);
4740  }
4741  return return_value;
4742  }
4743  }
4744 
4745 
4746 
4747  template <int dim, int spacedim>
4749  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
4750  const unsigned int q_point) const
4751  {
4752  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4753  Assert(fe_values->update_flags & update_gradients,
4755  "update_gradients")));
4756 
4757  const int snc =
4758  shape_function_data[shape_function].single_nonzero_component;
4759 
4760  if (snc == -2)
4761  {
4762  // shape function is zero for the selected components
4763  return divergence_type();
4764  }
4765  else if (snc != -1)
4766  {
4767  // we have a single non-zero component when the tensor is
4768  // represented in unrolled form.
4769  //
4770  // the divergence of a second-order tensor is a first order tensor.
4771  //
4772  // assume the second-order tensor is A with components A_{ij},
4773  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
4774  //
4775  // Now, we know the nonzero component in unrolled form: it is indicated
4776  // by 'snc'. we can figure out which tensor components belong to this:
4777  const unsigned int comp =
4778  shape_function_data[shape_function].single_nonzero_component_index;
4779  const TableIndices<2> indices =
4781  const unsigned int ii = indices[0];
4782  const unsigned int jj = indices[1];
4783 
4784  const ::Tensor<1, spacedim> &phi_grad =
4785  fe_values->finite_element_output.shape_gradients[snc][q_point];
4786 
4787  divergence_type return_value;
4788  // note that we contract \nabla from the right
4789  return_value[ii] = phi_grad[jj];
4790 
4791  return return_value;
4792  }
4793  else
4794  {
4795  Assert(false, ExcNotImplemented());
4796  divergence_type return_value;
4797  return return_value;
4798  }
4799  }
4800 
4801 
4802 
4803  template <int dim, int spacedim>
4805  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
4806  const unsigned int q_point) const
4807  {
4808  AssertIndexRange(shape_function, fe_values->fe->dofs_per_cell);
4809  Assert(fe_values->update_flags & update_gradients,
4811  "update_gradients")));
4812 
4813  const int snc =
4814  shape_function_data[shape_function].single_nonzero_component;
4815 
4816  if (snc == -2)
4817  {
4818  // shape function is zero for the selected components
4819  return gradient_type();
4820  }
4821  else if (snc != -1)
4822  {
4823  // we have a single non-zero component when the tensor is
4824  // represented in unrolled form.
4825  //
4826  // the gradient of a second-order tensor is a third order tensor.
4827  //
4828  // assume the second-order tensor is A with components A_{ij},
4829  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
4830  //
4831  // Now, we know the nonzero component in unrolled form: it is indicated
4832  // by 'snc'. we can figure out which tensor components belong to this:
4833  const unsigned int comp =
4834  shape_function_data[shape_function].single_nonzero_component_index;
4835  const TableIndices<2> indices =
4837  const unsigned int ii = indices[0];
4838  const unsigned int jj = indices[1];
4839 
4840  const ::Tensor<1, spacedim> &phi_grad =
4841  fe_values->finite_element_output.shape_gradients[snc][q_point];
4842 
4843  gradient_type return_value;
4844  return_value[ii][jj] = phi_grad;
4845 
4846  return return_value;
4847  }
4848  else
4849  {
4850  Assert(false, ExcNotImplemented());
4851  gradient_type return_value;
4852  return return_value;
4853  }
4854  }
4855 
4856 } // namespace FEValuesViews
4857 
4858 
4859 
4860 /*---------------------- Inline functions: FEValuesBase ---------------------*/
4861 
4862 
4863 
4864 template <int dim, int spacedim>
4866  operator[](const FEValuesExtractors::Scalar &scalar) const
4867 {
4868  AssertIndexRange(scalar.component, fe_values_views_cache.scalars.size());
4869 
4870  return fe_values_views_cache.scalars[scalar.component];
4871 }
4872 
4873 
4874 
4875 template <int dim, int spacedim>
4877  operator[](const FEValuesExtractors::Vector &vector) const
4878 {
4880  fe_values_views_cache.vectors.size());
4881 
4882  return fe_values_views_cache.vectors[vector.first_vector_component];
4883 }
4884 
4885 
4886 
4887 template <int dim, int spacedim>
4891 {
4892  Assert(
4893  tensor.first_tensor_component <
4894  fe_values_views_cache.symmetric_second_order_tensors.size(),
4896  0,
4897  fe_values_views_cache.symmetric_second_order_tensors.size()));
4898 
4899  return fe_values_views_cache
4900  .symmetric_second_order_tensors[tensor.first_tensor_component];
4901 }
4902 
4903 
4904 
4905 template <int dim, int spacedim>
4908  operator[](const FEValuesExtractors::Tensor<2> &tensor) const
4909 {
4911  fe_values_views_cache.second_order_tensors.size());
4912 
4913  return fe_values_views_cache
4914  .second_order_tensors[tensor.first_tensor_component];
4915 }
4916 
4917 
4918 
4919 template <int dim, int spacedim>
4920 inline const double &
4921 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
4922  const unsigned int j) const
4923 {
4924  AssertIndexRange(i, fe->dofs_per_cell);
4925  Assert(this->update_flags & update_values,
4926  ExcAccessToUninitializedField("update_values"));
4927  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4928  Assert(present_cell.get() != nullptr,
4929  ExcMessage("FEValues object is not reinit'ed to any cell"));
4930  // if the entire FE is primitive,
4931  // then we can take a short-cut:
4932  if (fe->is_primitive())
4933  return this->finite_element_output.shape_values(i, j);
4934  else
4935  {
4936  // otherwise, use the mapping
4937  // between shape function
4938  // numbers and rows. note that
4939  // by the assertions above, we
4940  // know that this particular
4941  // shape function is primitive,
4942  // so we can call
4943  // system_to_component_index
4944  const unsigned int row =
4945  this->finite_element_output
4946  .shape_function_to_row_table[i * fe->n_components() +
4947  fe->system_to_component_index(i).first];
4948  return this->finite_element_output.shape_values(row, j);
4949  }
4950 }
4951 
4952 
4953 
4954 template <int dim, int spacedim>
4955 inline double
4957  const unsigned int i,
4958  const unsigned int j,
4959  const unsigned int component) const
4960 {
4961  AssertIndexRange(i, fe->dofs_per_cell);
4962  Assert(this->update_flags & update_values,
4963  ExcAccessToUninitializedField("update_values"));
4964  AssertIndexRange(component, fe->n_components());
4965  Assert(present_cell.get() != nullptr,
4966  ExcMessage("FEValues object is not reinit'ed to any cell"));
4967 
4968  // check whether the shape function
4969  // is non-zero at all within
4970  // this component:
4971  if (fe->get_nonzero_components(i)[component] == false)
4972  return 0;
4973 
4974  // look up the right row in the
4975  // table and take the data from
4976  // there
4977  const unsigned int row =
4978  this->finite_element_output
4979  .shape_function_to_row_table[i * fe->n_components() + component];
4980  return this->finite_element_output.shape_values(row, j);
4981 }
4982 
4983 
4984 
4985 template <int dim, int spacedim>
4986 inline const Tensor<1, spacedim> &
4987 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
4988  const unsigned int j) const
4989 {
4990  AssertIndexRange(i, fe->dofs_per_cell);
4991  Assert(this->update_flags & update_gradients,
4992  ExcAccessToUninitializedField("update_gradients"));
4993  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4994  Assert(present_cell.get() != nullptr,
4995  ExcMessage("FEValues object is not reinit'ed to any cell"));
4996  // if the entire FE is primitive,
4997  // then we can take a short-cut:
4998  if (fe->is_primitive())
4999  return this->finite_element_output.shape_gradients[i][j];
5000  else
5001  {
5002  // otherwise, use the mapping
5003  // between shape function
5004  // numbers and rows. note that
5005  // by the assertions above, we
5006  // know that this particular
5007  // shape function is primitive,
5008  // so we can call
5009  // system_to_component_index
5010  const unsigned int row =
5011  this->finite_element_output
5012  .shape_function_to_row_table[i * fe->n_components() +
5013  fe->system_to_component_index(i).first];
5014  return this->finite_element_output.shape_gradients[row][j];
5015  }
5016 }
5017 
5018 
5019 
5020 template <int dim, int spacedim>
5021 inline Tensor<1, spacedim>
5023  const unsigned int i,
5024  const unsigned int j,
5025  const unsigned int component) const
5026 {
5027  AssertIndexRange(i, fe->dofs_per_cell);
5028  Assert(this->update_flags & update_gradients,
5029  ExcAccessToUninitializedField("update_gradients"));
5030  AssertIndexRange(component, fe->n_components());
5031  Assert(present_cell.get() != nullptr,
5032  ExcMessage("FEValues object is not reinit'ed to any cell"));
5033  // check whether the shape function
5034  // is non-zero at all within
5035  // this component:
5036  if (fe->get_nonzero_components(i)[component] == false)
5037  return Tensor<1, spacedim>();
5038 
5039  // look up the right row in the
5040  // table and take the data from
5041  // there
5042  const unsigned int row =
5043  this->finite_element_output
5044  .shape_function_to_row_table[i * fe->n_components() + component];
5045  return this->finite_element_output.shape_gradients[row][j];
5046 }
5047 
5048 
5049 
5050 template <int dim, int spacedim>
5051 inline const Tensor<2, spacedim> &
5052 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
5053  const unsigned int j) const
5054 {
5055  AssertIndexRange(i, fe->dofs_per_cell);
5056  Assert(this->update_flags & update_hessians,
5057  ExcAccessToUninitializedField("update_hessians"));
5058  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5059  Assert(present_cell.get() != nullptr,
5060  ExcMessage("FEValues object is not reinit'ed to any cell"));
5061  // if the entire FE is primitive,
5062  // then we can take a short-cut:
5063  if (fe->is_primitive())
5064  return this->finite_element_output.shape_hessians[i][j];
5065  else
5066  {
5067  // otherwise, use the mapping
5068  // between shape function
5069  // numbers and rows. note that
5070  // by the assertions above, we
5071  // know that this particular
5072  // shape function is primitive,
5073  // so we can call
5074  // system_to_component_index
5075  const unsigned int row =
5076  this->finite_element_output
5077  .shape_function_to_row_table[i * fe->n_components() +
5078  fe->system_to_component_index(i).first];
5079  return this->finite_element_output.shape_hessians[row][j];
5080  }
5081 }
5082 
5083 
5084 
5085 template <int dim, int spacedim>
5086 inline Tensor<2, spacedim>
5088  const unsigned int i,
5089  const unsigned int j,
5090  const unsigned int component) const
5091 {
5092  AssertIndexRange(i, fe->dofs_per_cell);
5093  Assert(this->update_flags & update_hessians,
5094  ExcAccessToUninitializedField("update_hessians"));
5095  AssertIndexRange(component, fe->n_components());
5096  Assert(present_cell.get() != nullptr,
5097  ExcMessage("FEValues object is not reinit'ed to any cell"));
5098  // check whether the shape function
5099  // is non-zero at all within
5100  // this component:
5101  if (fe->get_nonzero_components(i)[component] == false)
5102  return Tensor<2, spacedim>();
5103 
5104  // look up the right row in the
5105  // table and take the data from
5106  // there
5107  const unsigned int row =
5108  this->finite_element_output
5109  .shape_function_to_row_table[i * fe->n_components() + component];
5110  return this->finite_element_output.shape_hessians[row][j];
5111 }
5112 
5113 
5114 
5115 template <int dim, int spacedim>
5116 inline const Tensor<3, spacedim> &
5118  const unsigned int j) const
5119 {
5120  AssertIndexRange(i, fe->dofs_per_cell);
5121  Assert(this->update_flags & update_3rd_derivatives,
5122  ExcAccessToUninitializedField("update_3rd_derivatives"));
5123  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5124  Assert(present_cell.get() != nullptr,
5125  ExcMessage("FEValues object is not reinit'ed to any cell"));
5126  // if the entire FE is primitive,
5127  // then we can take a short-cut:
5128  if (fe->is_primitive())
5129  return this->finite_element_output.shape_3rd_derivatives[i][j];
5130  else
5131  {
5132  // otherwise, use the mapping
5133  // between shape function
5134  // numbers and rows. note that
5135  // by the assertions above, we
5136  // know that this particular
5137  // shape function is primitive,
5138  // so we can call
5139  // system_to_component_index
5140  const unsigned int row =
5141  this->finite_element_output
5142  .shape_function_to_row_table[i * fe->n_components() +
5143  fe->system_to_component_index(i).first];
5144  return this->finite_element_output.shape_3rd_derivatives[row][j];
5145  }
5146 }
5147 
5148 
5149 
5150 template <int dim, int spacedim>
5151 inline Tensor<3, spacedim>
5153  const unsigned int i,
5154  const unsigned int j,
5155  const unsigned int component) const
5156 {
5157  AssertIndexRange(i, fe->dofs_per_cell);
5158  Assert(this->update_flags & update_3rd_derivatives,
5159  ExcAccessToUninitializedField("update_3rd_derivatives"));
5160  AssertIndexRange(component, fe->n_components());
5161  Assert(present_cell.get() != nullptr,
5162  ExcMessage("FEValues object is not reinit'ed to any cell"));
5163  // check whether the shape function
5164  // is non-zero at all within
5165  // this component:
5166  if (fe->get_nonzero_components(i)[component] == false)
5167  return Tensor<3, spacedim>();
5168 
5169  // look up the right row in the
5170  // table and take the data from
5171  // there
5172  const unsigned int row =
5173  this->finite_element_output
5174  .shape_function_to_row_table[i * fe->n_components() + component];
5175  return this->finite_element_output.shape_3rd_derivatives[row][j];
5176 }
5177 
5178 
5179 
5180 template <int dim, int spacedim>
5181 inline const FiniteElement<dim, spacedim> &
5183 {
5184  return *fe;
5185 }
5186 
5187 
5188 
5189 template <int dim, int spacedim>
5190 inline const Mapping<dim, spacedim> &
5192 {
5193  return *mapping;
5194 }
5195 
5196 
5197 
5198 template <int dim, int spacedim>
5199 inline UpdateFlags
5201 {
5202  return this->update_flags;
5203 }
5204 
5205 
5206 
5207 template <int dim, int spacedim>
5208 inline const std::vector<Point<spacedim>> &
5210 {
5211  Assert(this->update_flags & update_quadrature_points,
5212  ExcAccessToUninitializedField("update_quadrature_points"));
5213  Assert(present_cell.get() != nullptr,
5214  ExcMessage("FEValues object is not reinit'ed to any cell"));
5215  return this->mapping_output.quadrature_points;
5216 }
5217 
5218 
5219 
5220 template <int dim, int spacedim>
5221 inline const std::vector<double> &
5223 {
5224  Assert(this->update_flags & update_JxW_values,
5225  ExcAccessToUninitializedField("update_JxW_values"));
5226  Assert(present_cell.get() != nullptr,
5227  ExcMessage("FEValues object is not reinit'ed to any cell"));
5228  return this->mapping_output.JxW_values;
5229 }
5230 
5231 
5232 
5233 template <int dim, int spacedim>
5234 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5236 {
5237  Assert(this->update_flags & update_jacobians,
5238  ExcAccessToUninitializedField("update_jacobians"));
5239  Assert(present_cell.get() != nullptr,
5240  ExcMessage("FEValues object is not reinit'ed to any cell"));
5241  return this->mapping_output.jacobians;
5242 }
5243 
5244 
5245 
5246 template <int dim, int spacedim>
5247 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5249 {
5250  Assert(this->update_flags & update_jacobian_grads,
5251  ExcAccessToUninitializedField("update_jacobians_grads"));
5252  Assert(present_cell.get() != nullptr,
5253  ExcMessage("FEValues object is not reinit'ed to any cell"));
5254  return this->mapping_output.jacobian_grads;
5255 }
5256 
5257 
5258 
5259 template <int dim, int spacedim>
5260 inline const Tensor<3, spacedim> &
5262  const unsigned int i) const
5263 {
5264  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5265  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5266  Assert(present_cell.get() != nullptr,
5267  ExcMessage("FEValues object is not reinit'ed to any cell"));
5268  return this->mapping_output.jacobian_pushed_forward_grads[i];
5269 }
5270 
5271 
5272 
5273 template <int dim, int spacedim>
5274 inline const std::vector<Tensor<3, spacedim>> &
5276 {
5277  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5278  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5279  Assert(present_cell.get() != nullptr,
5280  ExcMessage("FEValues object is not reinit'ed to any cell"));
5281  return this->mapping_output.jacobian_pushed_forward_grads;
5282 }
5283 
5284 
5285 
5286 template <int dim, int spacedim>
5287 inline const DerivativeForm<3, dim, spacedim> &
5288 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5289 {
5290  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5291  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5292  Assert(present_cell.get() != nullptr,
5293  ExcMessage("FEValues object is not reinit'ed to any cell"));
5294  return this->mapping_output.jacobian_2nd_derivatives[i];
5295 }
5296 
5297 
5298 
5299 template <int dim, int spacedim>
5300 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5302 {
5303  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5304  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5305  Assert(present_cell.get() != nullptr,
5306  ExcMessage("FEValues object is not reinit'ed to any cell"));
5307  return this->mapping_output.jacobian_2nd_derivatives;
5308 }
5309 
5310 
5311 
5312 template <int dim, int spacedim>
5313 inline const Tensor<4, spacedim> &
5315  const unsigned int i) const
5316 {
5319  "update_jacobian_pushed_forward_2nd_derivatives"));
5320  Assert(present_cell.get() != nullptr,
5321  ExcMessage("FEValues object is not reinit'ed to any cell"));
5322  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5323 }
5324 
5325 
5326 
5327 template <int dim, int spacedim>
5328 inline const std::vector<Tensor<4, spacedim>> &
5330 {
5331  Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
5333  "update_jacobian_pushed_forward_2nd_derivatives"));
5334  Assert(present_cell.get() != nullptr,
5335  ExcMessage("FEValues object is not reinit'ed to any cell"));
5336  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5337 }
5338 
5339 
5340 
5341 template <int dim, int spacedim>
5342 inline const DerivativeForm<4, dim, spacedim> &
5343 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5344 {
5345  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5346  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5347  Assert(present_cell.get() != nullptr,
5348  ExcMessage("FEValues object is not reinit'ed to any cell"));
5349  return this->mapping_output.jacobian_3rd_derivatives[i];
5350 }
5351 
5352 
5353 
5354 template <int dim, int spacedim>
5355 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5357 {
5358  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5359  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5360  Assert(present_cell.get() != nullptr,
5361  ExcMessage("FEValues object is not reinit'ed to any cell"));
5362  return this->mapping_output.jacobian_3rd_derivatives;
5363 }
5364 
5365 
5366 
5367 template <int dim, int spacedim>
5368 inline const Tensor<5, spacedim> &
5370  const unsigned int i) const
5371 {
5374  "update_jacobian_pushed_forward_3rd_derivatives"));
5375  Assert(present_cell.get() != nullptr,
5376  ExcMessage("FEValues object is not reinit'ed to any cell"));
5377  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
5378 }
5379 
5380 
5381 
5382 template <int dim, int spacedim>
5383 inline const std::vector<Tensor<5, spacedim>> &
5385 {
5386  Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
5388  "update_jacobian_pushed_forward_3rd_derivatives"));
5389  Assert(present_cell.get() != nullptr,
5390  ExcMessage("FEValues object is not reinit'ed to any cell"));
5391  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
5392 }
5393 
5394 
5395 
5396 template <int dim, int spacedim>
5397 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
5399 {
5400  Assert(this->update_flags & update_inverse_jacobians,
5401  ExcAccessToUninitializedField("update_inverse_jacobians"));
5402  Assert(present_cell.get() != nullptr,
5403  ExcMessage("FEValues object is not reinit'ed to any cell"));
5404  return this->mapping_output.inverse_jacobians;
5405 }
5406 
5407 
5408 
5409 template <int dim, int spacedim>
5410 inline const Point<spacedim> &
5411 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
5412 {
5413  Assert(this->update_flags & update_quadrature_points,
5414  ExcAccessToUninitializedField("update_quadrature_points"));
5415  AssertIndexRange(i, this->mapping_output.quadrature_points.size());
5416  Assert(present_cell.get() != nullptr,
5417  ExcMessage("FEValues object is not reinit'ed to any cell"));
5418 
5419  return this->mapping_output.quadrature_points[i];
5420 }
5421 
5422 
5423 
5424 template <int dim, int spacedim>
5425 inline double
5426 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
5427 {
5428  Assert(this->update_flags & update_JxW_values,
5429  ExcAccessToUninitializedField("update_JxW_values"));
5430  AssertIndexRange(i, this->mapping_output.JxW_values.size());
5431  Assert(present_cell.get() != nullptr,
5432  ExcMessage("FEValues object is not reinit'ed to any cell"));
5433 
5434  return this->mapping_output.JxW_values[i];
5435 }
5436 
5437 
5438 
5439 template <int dim, int spacedim>
5440 inline const DerivativeForm<1, dim, spacedim> &
5441 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
5442 {
5443  Assert(this->update_flags & update_jacobians,
5444  ExcAccessToUninitializedField("update_jacobians"));
5445  AssertIndexRange(i, this->mapping_output.jacobians.size());
5446  Assert(present_cell.get() != nullptr,
5447  ExcMessage("FEValues object is not reinit'ed to any cell"));
5448 
5449  return this->mapping_output.jacobians[i];
5450 }
5451 
5452 
5453 
5454 template <int dim, int spacedim>
5455 inline const DerivativeForm<2, dim, spacedim> &
5456 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
5457 {
5458  Assert(this->update_flags & update_jacobian_grads,
5459  ExcAccessToUninitializedField("update_jacobians_grads"));
5460  AssertIndexRange(i, this->mapping_output.jacobian_grads.size());
5461  Assert(present_cell.get() != nullptr,
5462  ExcMessage("FEValues object is not reinit'ed to any cell"));
5463 
5464  return this->mapping_output.jacobian_grads[i];
5465 }
5466 
5467 
5468 
5469 template <int dim, int spacedim>
5470 inline const DerivativeForm<1, spacedim, dim> &
5471 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
5472 {
5473  Assert(this->update_flags & update_inverse_jacobians,
5474  ExcAccessToUninitializedField("update_inverse_jacobians"));
5475  AssertIndexRange(i, this->mapping_output.inverse_jacobians.size());
5476  Assert(present_cell.get() != nullptr,
5477  ExcMessage("FEValues object is not reinit'ed to any cell"));
5478 
5479  return this->mapping_output.inverse_jacobians[i];
5480 }
5481 
5482 
5483 
5484 template <int dim, int spacedim>
5485 inline const Tensor<1, spacedim> &
5486 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
5487 {
5488  Assert(this->update_flags & update_normal_vectors,
5490  "update_normal_vectors")));
5491  AssertIndexRange(i, this->mapping_output.normal_vectors.size());
5492  Assert(present_cell.get() != nullptr,
5493  ExcMessage("FEValues object is not reinit'ed to any cell"));
5494 
5495  return this->mapping_output.normal_vectors[i];
5496 }
5497 
5498 
5499 
5500 /*--------------------- Inline functions: FEValues --------------------------*/
5501 
5502 
5503 template <int dim, int spacedim>
5504 inline const Quadrature<dim> &
5506 {
5507  return quadrature;
5508 }
5509 
5510 
5511 
5512 template <int dim, int spacedim>
5513 inline const FEValues<dim, spacedim> &
5515 {
5516  return *this;
5517 }
5518 
5519 
5520 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
5521 
5522 
5523 template <int dim, int spacedim>
5524 inline unsigned int
5526 {
5527  return present_face_index;
5528 }
5529 
5530 
5531 /*----------------------- Inline functions: FE*FaceValues -------------------*/
5532 
5533 template <int dim, int spacedim>
5534 inline const Quadrature<dim - 1> &
5536 {
5537  return quadrature;
5538 }
5539 
5540 
5541 
5542 template <int dim, int spacedim>
5543 inline const FEFaceValues<dim, spacedim> &
5545 {
5546  return *this;
5547 }
5548 
5549 
5550 
5551 template <int dim, int spacedim>
5552 inline const FESubfaceValues<dim, spacedim> &
5554 {
5555  return *this;
5556 }
5557 
5558 
5559 
5560 template <int dim, int spacedim>
5561 inline const Tensor<1, spacedim> &
5562 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
5563 {
5564  AssertIndexRange(i, this->mapping_output.boundary_forms.size());
5565  Assert(this->update_flags & update_boundary_forms,
5567  "update_boundary_forms")));
5568 
5569  return this->mapping_output.boundary_forms[i];
5570 }
5571 
5572 #endif // DOXYGEN
5573 
5574 DEAL_II_NAMESPACE_CLOSE
5575 
5576 #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:3419
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:3451
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:3684
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1520
constexpr SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
static ::ExceptionBase & ExcAccessToUninitializedField()
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h: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:3689
Volume element.
#define AssertIndexRange(index, range)
Definition: exceptions.h:1649
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:3335
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:3387
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:3466
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:1419
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:3582
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:3395
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:451
const std::vector< DerivativeForm< 4, dim, spacedim > > & get_jacobian_3rd_derivatives() const
static constexpr TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
Definition: tensor.h:1541
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:3360
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:417
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
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1216
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:3351
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3427
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:3402
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:3433
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:3411
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