Reference documentation for deal.II version Git 7026f387cc 2019-10-15 14:19:01 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
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 
2085 template <int dim, int spacedim>
2086 class FEValuesBase : public Subscriptor
2087 {
2088 public:
2092  static const unsigned int dimension = dim;
2093 
2097  static const unsigned int space_dimension = spacedim;
2098 
2102  const unsigned int n_quadrature_points;
2103 
2109  const unsigned int dofs_per_cell;
2110 
2111 
2119  FEValuesBase(const unsigned int n_q_points,
2120  const unsigned int dofs_per_cell,
2121  const UpdateFlags update_flags,
2122  const Mapping<dim, spacedim> & mapping,
2123  const FiniteElement<dim, spacedim> &fe);
2124 
2129  FEValuesBase &
2130  operator=(const FEValuesBase &) = delete;
2131 
2136  FEValuesBase(const FEValuesBase &) = delete;
2137 
2141  virtual ~FEValuesBase() override;
2142 
2143 
2147 
2148 
2169  const double &
2170  shape_value(const unsigned int function_no,
2171  const unsigned int point_no) const;
2172 
2193  double
2194  shape_value_component(const unsigned int function_no,
2195  const unsigned int point_no,
2196  const unsigned int component) const;
2197 
2223  const Tensor<1, spacedim> &
2224  shape_grad(const unsigned int function_no,
2225  const unsigned int quadrature_point) const;
2226 
2244  shape_grad_component(const unsigned int function_no,
2245  const unsigned int point_no,
2246  const unsigned int component) const;
2247 
2267  const Tensor<2, spacedim> &
2268  shape_hessian(const unsigned int function_no,
2269  const unsigned int point_no) const;
2270 
2288  shape_hessian_component(const unsigned int function_no,
2289  const unsigned int point_no,
2290  const unsigned int component) const;
2291 
2311  const Tensor<3, spacedim> &
2312  shape_3rd_derivative(const unsigned int function_no,
2313  const unsigned int point_no) const;
2314 
2332  shape_3rd_derivative_component(const unsigned int function_no,
2333  const unsigned int point_no,
2334  const unsigned int component) const;
2335 
2337 
2339 
2376  template <class InputVector>
2377  void
2378  get_function_values(
2379  const InputVector & fe_function,
2380  std::vector<typename InputVector::value_type> &values) const;
2381 
2395  template <class InputVector>
2396  void
2397  get_function_values(
2398  const InputVector & fe_function,
2399  std::vector<Vector<typename InputVector::value_type>> &values) const;
2400 
2419  template <class InputVector>
2420  void
2421  get_function_values(
2422  const InputVector & fe_function,
2424  std::vector<typename InputVector::value_type> & values) const;
2425 
2447  template <class InputVector>
2448  void
2449  get_function_values(
2450  const InputVector & fe_function,
2452  std::vector<Vector<typename InputVector::value_type>> &values) const;
2453 
2454 
2485  template <class InputVector>
2486  void
2487  get_function_values(
2488  const InputVector & fe_function,
2490  ArrayView<std::vector<typename InputVector::value_type>> values,
2491  const bool quadrature_points_fastest) const;
2492 
2494 
2496 
2533  template <class InputVector>
2534  void
2535  get_function_gradients(
2536  const InputVector &fe_function,
2538  &gradients) const;
2539 
2556  template <class InputVector>
2557  void
2558  get_function_gradients(
2559  const InputVector &fe_function,
2560  std::vector<
2562  &gradients) const;
2563 
2570  template <class InputVector>
2571  void
2572  get_function_gradients(
2573  const InputVector & fe_function,
2576  &gradients) const;
2577 
2584  template <class InputVector>
2585  void
2586  get_function_gradients(
2587  const InputVector & fe_function,
2589  ArrayView<
2591  gradients,
2592  bool quadrature_points_fastest = false) const;
2593 
2595 
2599 
2637  template <class InputVector>
2638  void
2639  get_function_hessians(
2640  const InputVector &fe_function,
2642  &hessians) const;
2643 
2661  template <class InputVector>
2662  void
2663  get_function_hessians(
2664  const InputVector &fe_function,
2665  std::vector<
2667  & hessians,
2668  bool quadrature_points_fastest = false) const;
2669 
2674  template <class InputVector>
2675  void
2676  get_function_hessians(
2677  const InputVector & fe_function,
2680  &hessians) const;
2681 
2688  template <class InputVector>
2689  void
2690  get_function_hessians(
2691  const InputVector & fe_function,
2693  ArrayView<
2695  hessians,
2696  bool quadrature_points_fastest = false) const;
2697 
2738  template <class InputVector>
2739  void
2740  get_function_laplacians(
2741  const InputVector & fe_function,
2742  std::vector<typename InputVector::value_type> &laplacians) const;
2743 
2763  template <class InputVector>
2764  void
2765  get_function_laplacians(
2766  const InputVector & fe_function,
2767  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2768 
2775  template <class InputVector>
2776  void
2777  get_function_laplacians(
2778  const InputVector & fe_function,
2780  std::vector<typename InputVector::value_type> & laplacians) const;
2781 
2788  template <class InputVector>
2789  void
2790  get_function_laplacians(
2791  const InputVector & fe_function,
2793  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2794 
2801  template <class InputVector>
2802  void
2803  get_function_laplacians(
2804  const InputVector & fe_function,
2806  std::vector<std::vector<typename InputVector::value_type>> &laplacians,
2807  bool quadrature_points_fastest = false) const;
2808 
2810 
2812 
2851  template <class InputVector>
2852  void
2853  get_function_third_derivatives(
2854  const InputVector &fe_function,
2856  &third_derivatives) const;
2857 
2876  template <class InputVector>
2877  void
2878  get_function_third_derivatives(
2879  const InputVector &fe_function,
2880  std::vector<
2882  & third_derivatives,
2883  bool quadrature_points_fastest = false) const;
2884 
2889  template <class InputVector>
2890  void
2891  get_function_third_derivatives(
2892  const InputVector & fe_function,
2895  &third_derivatives) const;
2896 
2903  template <class InputVector>
2904  void
2905  get_function_third_derivatives(
2906  const InputVector & fe_function,
2908  ArrayView<
2910  third_derivatives,
2911  bool quadrature_points_fastest = false) const;
2913 
2915 
2916 
2922  const Point<spacedim> &
2923  quadrature_point(const unsigned int q) const;
2924 
2930  const std::vector<Point<spacedim>> &
2931  get_quadrature_points() const;
2932 
2948  double
2949  JxW(const unsigned int quadrature_point) const;
2950 
2954  const std::vector<double> &
2955  get_JxW_values() const;
2956 
2964  jacobian(const unsigned int quadrature_point) const;
2965 
2972  const std::vector<DerivativeForm<1, dim, spacedim>> &
2973  get_jacobians() const;
2974 
2983  jacobian_grad(const unsigned int quadrature_point) const;
2984 
2991  const std::vector<DerivativeForm<2, dim, spacedim>> &
2992  get_jacobian_grads() const;
2993 
3002  const Tensor<3, spacedim> &
3003  jacobian_pushed_forward_grad(const unsigned int quadrature_point) const;
3004 
3011  const std::vector<Tensor<3, spacedim>> &
3012  get_jacobian_pushed_forward_grads() const;
3013 
3022  jacobian_2nd_derivative(const unsigned int quadrature_point) const;
3023 
3030  const std::vector<DerivativeForm<3, dim, spacedim>> &
3031  get_jacobian_2nd_derivatives() const;
3032 
3042  const Tensor<4, spacedim> &
3043  jacobian_pushed_forward_2nd_derivative(
3044  const unsigned int quadrature_point) const;
3045 
3052  const std::vector<Tensor<4, spacedim>> &
3053  get_jacobian_pushed_forward_2nd_derivatives() const;
3054 
3064  jacobian_3rd_derivative(const unsigned int quadrature_point) const;
3065 
3072  const std::vector<DerivativeForm<4, dim, spacedim>> &
3073  get_jacobian_3rd_derivatives() const;
3074 
3084  const Tensor<5, spacedim> &
3085  jacobian_pushed_forward_3rd_derivative(
3086  const unsigned int quadrature_point) const;
3087 
3094  const std::vector<Tensor<5, spacedim>> &
3095  get_jacobian_pushed_forward_3rd_derivatives() const;
3096 
3104  inverse_jacobian(const unsigned int quadrature_point) const;
3105 
3112  const std::vector<DerivativeForm<1, spacedim, dim>> &
3113  get_inverse_jacobians() const;
3114 
3128  const Tensor<1, spacedim> &
3129  normal_vector(const unsigned int i) const;
3130 
3141  DEAL_II_DEPRECATED
3142  const std::vector<Tensor<1, spacedim>> &
3143  get_all_normal_vectors() const;
3144 
3152  const std::vector<Tensor<1, spacedim>> &
3153  get_normal_vectors() const;
3154 
3156 
3158 
3159 
3169  operator[](const FEValuesExtractors::Scalar &scalar) const;
3170 
3180  operator[](const FEValuesExtractors::Vector &vector) const;
3181 
3192  operator[](const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
3193 
3194 
3204  operator[](const FEValuesExtractors::Tensor<2> &tensor) const;
3205 
3207 
3209 
3210 
3214  const Mapping<dim, spacedim> &
3215  get_mapping() const;
3216 
3221  get_fe() const;
3222 
3226  UpdateFlags
3227  get_update_flags() const;
3228 
3233  get_cell() const;
3234 
3241  get_cell_similarity() const;
3242 
3247  std::size_t
3248  memory_consumption() const;
3250 
3251 
3260  std::string,
3261  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3262  << "object for which this kind of information has not been computed. What "
3263  << "information these objects compute is determined by the update_* flags you "
3264  << "pass to the constructor. Here, the operation you are attempting requires "
3265  << "the <" << arg1
3266  << "> flag to be set, but it was apparently not specified "
3267  << "upon construction.");
3268 
3276  ExcFEDontMatch,
3277  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3278  "to the DoFHandler that provided the cell iterator do not match.");
3284  DeclException1(ExcShapeFunctionNotPrimitive,
3285  int,
3286  << "The shape function with index " << arg1
3287  << " is not primitive, i.e. it is vector-valued and "
3288  << "has more than one non-zero vector component. This "
3289  << "function cannot be called for these shape functions. "
3290  << "Maybe you want to use the same function with the "
3291  << "_component suffix?");
3292 
3301  "The given FiniteElement is not a primitive element but the requested operation "
3302  "only works for those. See FiniteElement::is_primitive() for more information.");
3303 
3304 protected:
3333  class CellIteratorBase;
3334 
3339  template <typename CI>
3341  class TriaCellIterator;
3342 
3348  std::unique_ptr<const CellIteratorBase> present_cell;
3349 
3357  boost::signals2::connection tria_listener_refinement;
3358 
3366  boost::signals2::connection tria_listener_mesh_transform;
3367 
3373  void
3374  invalidate_present_cell();
3375 
3385  void
3386  maybe_invalidate_previous_present_cell(
3387  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3388 
3394 
3400  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3402 
3409 
3410 
3418 
3424  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3426 
3432  spacedim>
3434 
3435 
3440 
3449  UpdateFlags
3450  compute_update_flags(const UpdateFlags update_flags) const;
3451 
3458 
3464  void
3465  check_cell_similarity(
3466  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3467 
3468 private:
3473 
3474  // Make the view classes friends of this class, since they access internal
3475  // data.
3476  template <int, int>
3477  friend class FEValuesViews::Scalar;
3478  template <int, int>
3479  friend class FEValuesViews::Vector;
3480  template <int, int, int>
3481  friend class FEValuesViews::SymmetricTensor;
3482  template <int, int, int>
3483  friend class FEValuesViews::Tensor;
3484 };
3485 
3486 
3487 
3498 template <int dim, int spacedim = dim>
3499 class FEValues : public FEValuesBase<dim, spacedim>
3500 {
3501 public:
3506  static const unsigned int integral_dimension = dim;
3507 
3512  FEValues(const Mapping<dim, spacedim> & mapping,
3513  const FiniteElement<dim, spacedim> &fe,
3514  const Quadrature<dim> & quadrature,
3515  const UpdateFlags update_flags);
3516 
3523  const Quadrature<dim> & quadrature,
3524  const UpdateFlags update_flags);
3525 
3532  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3533  void
3534  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3535  level_dof_access>> &cell);
3536 
3550  void
3551  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3552 
3557  const Quadrature<dim> &
3558  get_quadrature() const;
3559 
3564  std::size_t
3565  memory_consumption() const;
3566 
3581  const FEValues<dim, spacedim> &
3582  get_present_fe_values() const;
3583 
3584 private:
3589 
3593  void
3594  initialize(const UpdateFlags update_flags);
3595 
3602  void
3603  do_reinit();
3604 };
3605 
3606 
3617 template <int dim, int spacedim = dim>
3618 class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
3619 {
3620 public:
3625  static const unsigned int integral_dimension = dim - 1;
3626 
3638  FEFaceValuesBase(const unsigned int n_q_points,
3639  const unsigned int dofs_per_cell,
3640  const UpdateFlags update_flags,
3641  const Mapping<dim, spacedim> & mapping,
3642  const FiniteElement<dim, spacedim> &fe,
3643  const Quadrature<dim - 1> & quadrature);
3644 
3652  const Tensor<1, spacedim> &
3653  boundary_form(const unsigned int i) const;
3654 
3661  const std::vector<Tensor<1, spacedim>> &
3662  get_boundary_forms() const;
3663 
3668  unsigned int
3669  get_face_index() const;
3670 
3675  const Quadrature<dim - 1> &
3676  get_quadrature() const;
3677 
3682  std::size_t
3683  memory_consumption() const;
3684 
3685 protected:
3690  unsigned int present_face_index;
3691 
3695  const Quadrature<dim - 1> quadrature;
3696 };
3697 
3698 
3699 
3714 template <int dim, int spacedim = dim>
3715 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
3716 {
3717 public:
3722  static const unsigned int dimension = dim;
3723 
3724  static const unsigned int space_dimension = spacedim;
3725 
3730  static const unsigned int integral_dimension = dim - 1;
3731 
3736  FEFaceValues(const Mapping<dim, spacedim> & mapping,
3737  const FiniteElement<dim, spacedim> &fe,
3738  const Quadrature<dim - 1> & quadrature,
3739  const UpdateFlags update_flags);
3740 
3747  const Quadrature<dim - 1> & quadrature,
3748  const UpdateFlags update_flags);
3749 
3754  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3755  void
3756  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3757  level_dof_access>> &cell,
3758  const unsigned int face_no);
3759 
3773  void
3774  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3775  const unsigned int face_no);
3776 
3792  get_present_fe_values() const;
3793 
3794 private:
3798  void
3799  initialize(const UpdateFlags update_flags);
3800 
3807  void
3808  do_reinit(const unsigned int face_no);
3809 };
3810 
3811 
3829 template <int dim, int spacedim = dim>
3830 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
3831 {
3832 public:
3836  static const unsigned int dimension = dim;
3837 
3841  static const unsigned int space_dimension = spacedim;
3842 
3847  static const unsigned int integral_dimension = dim - 1;
3848 
3853  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
3854  const FiniteElement<dim, spacedim> &fe,
3855  const Quadrature<dim - 1> & face_quadrature,
3856  const UpdateFlags update_flags);
3857 
3864  const Quadrature<dim - 1> & face_quadrature,
3865  const UpdateFlags update_flags);
3866 
3873  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3874  void
3875  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3876  level_dof_access>> &cell,
3877  const unsigned int face_no,
3878  const unsigned int subface_no);
3879 
3893  void
3894  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3895  const unsigned int face_no,
3896  const unsigned int subface_no);
3897 
3913  get_present_fe_values() const;
3914 
3920  DeclException0(ExcReinitCalledWithBoundaryFace);
3921 
3927  DeclException0(ExcFaceHasNoSubfaces);
3928 
3929 private:
3933  void
3934  initialize(const UpdateFlags update_flags);
3935 
3942  void
3943  do_reinit(const unsigned int face_no, const unsigned int subface_no);
3944 };
3945 
3946 
3947 #ifndef DOXYGEN
3948 
3949 
3950 /*------------------------ Inline functions: namespace FEValuesViews --------*/
3951 
3952 namespace FEValuesViews
3953 {
3954  template <int dim, int spacedim>
3955  inline typename Scalar<dim, spacedim>::value_type
3956  Scalar<dim, spacedim>::value(const unsigned int shape_function,
3957  const unsigned int q_point) const
3958  {
3959  Assert(shape_function < fe_values->fe->dofs_per_cell,
3960  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3961  Assert(
3962  fe_values->update_flags & update_values,
3964  "update_values"))));
3965 
3966  // an adaptation of the FEValuesBase::shape_value_component function
3967  // except that here we know the component as fixed and we have
3968  // pre-computed and cached a bunch of information. See the comments there.
3969  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3970  return fe_values->finite_element_output.shape_values(
3971  shape_function_data[shape_function].row_index, q_point);
3972  else
3973  return 0;
3974  }
3975 
3976 
3977 
3978  template <int dim, int spacedim>
3979  inline typename Scalar<dim, spacedim>::gradient_type
3980  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
3981  const unsigned int q_point) const
3982  {
3983  Assert(shape_function < fe_values->fe->dofs_per_cell,
3984  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3985  Assert(fe_values->update_flags & update_gradients,
3987  "update_gradients")));
3988 
3989  // an adaptation of the FEValuesBase::shape_grad_component
3990  // function except that here we know the component as fixed and we have
3991  // pre-computed and cached a bunch of information. See the comments there.
3992  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3993  return fe_values->finite_element_output
3994  .shape_gradients[shape_function_data[shape_function].row_index]
3995  [q_point];
3996  else
3997  return gradient_type();
3998  }
3999 
4000 
4001 
4002  template <int dim, int spacedim>
4003  inline typename Scalar<dim, spacedim>::hessian_type
4004  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
4005  const unsigned int q_point) const
4006  {
4007  Assert(shape_function < fe_values->fe->dofs_per_cell,
4008  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4009  Assert(fe_values->update_flags & update_hessians,
4011  "update_hessians")));
4012 
4013  // an adaptation of the FEValuesBase::shape_hessian_component
4014  // function except that here we know the component as fixed and we have
4015  // pre-computed and cached a bunch of information. See the comments there.
4016  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4017  return fe_values->finite_element_output
4018  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
4019  else
4020  return hessian_type();
4021  }
4022 
4023 
4024 
4025  template <int dim, int spacedim>
4026  inline typename Scalar<dim, spacedim>::third_derivative_type
4027  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
4028  const unsigned int q_point) const
4029  {
4030  Assert(shape_function < fe_values->fe->dofs_per_cell,
4031  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4032  Assert(fe_values->update_flags & update_3rd_derivatives,
4034  "update_3rd_derivatives")));
4035 
4036  // an adaptation of the FEValuesBase::shape_3rdderivative_component
4037  // function except that here we know the component as fixed and we have
4038  // pre-computed and cached a bunch of information. See the comments there.
4039  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4040  return fe_values->finite_element_output
4041  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
4042  [q_point];
4043  else
4044  return third_derivative_type();
4045  }
4046 
4047 
4048 
4049  template <int dim, int spacedim>
4050  inline typename Vector<dim, spacedim>::value_type
4051  Vector<dim, spacedim>::value(const unsigned int shape_function,
4052  const unsigned int q_point) const
4053  {
4054  Assert(shape_function < fe_values->fe->dofs_per_cell,
4055  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4056  Assert(fe_values->update_flags & update_values,
4058  "update_values")));
4059 
4060  // same as for the scalar case except that we have one more index
4061  const int snc =
4062  shape_function_data[shape_function].single_nonzero_component;
4063  if (snc == -2)
4064  return value_type();
4065  else if (snc != -1)
4066  {
4067  value_type return_value;
4068  return_value[shape_function_data[shape_function]
4069  .single_nonzero_component_index] =
4070  fe_values->finite_element_output.shape_values(snc, q_point);
4071  return return_value;
4072  }
4073  else
4074  {
4075  value_type return_value;
4076  for (unsigned int d = 0; d < dim; ++d)
4077  if (shape_function_data[shape_function]
4078  .is_nonzero_shape_function_component[d])
4079  return_value[d] = fe_values->finite_element_output.shape_values(
4080  shape_function_data[shape_function].row_index[d], q_point);
4081 
4082  return return_value;
4083  }
4084  }
4085 
4086 
4087 
4088  template <int dim, int spacedim>
4089  inline typename Vector<dim, spacedim>::gradient_type
4090  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4091  const unsigned int q_point) const
4092  {
4093  Assert(shape_function < fe_values->fe->dofs_per_cell,
4094  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4095  Assert(fe_values->update_flags & update_gradients,
4097  "update_gradients")));
4098 
4099  // same as for the scalar case except that we have one more index
4100  const int snc =
4101  shape_function_data[shape_function].single_nonzero_component;
4102  if (snc == -2)
4103  return gradient_type();
4104  else if (snc != -1)
4105  {
4106  gradient_type return_value;
4107  return_value[shape_function_data[shape_function]
4108  .single_nonzero_component_index] =
4109  fe_values->finite_element_output.shape_gradients[snc][q_point];
4110  return return_value;
4111  }
4112  else
4113  {
4114  gradient_type return_value;
4115  for (unsigned int d = 0; d < dim; ++d)
4116  if (shape_function_data[shape_function]
4117  .is_nonzero_shape_function_component[d])
4118  return_value[d] =
4119  fe_values->finite_element_output.shape_gradients
4120  [shape_function_data[shape_function].row_index[d]][q_point];
4121 
4122  return return_value;
4123  }
4124  }
4125 
4126 
4127 
4128  template <int dim, int spacedim>
4129  inline typename Vector<dim, spacedim>::divergence_type
4130  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4131  const unsigned int q_point) const
4132  {
4133  // this function works like in the case above
4134  Assert(shape_function < fe_values->fe->dofs_per_cell,
4135  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4136  Assert(fe_values->update_flags & update_gradients,
4138  "update_gradients")));
4139 
4140  // same as for the scalar case except that we have one more index
4141  const int snc =
4142  shape_function_data[shape_function].single_nonzero_component;
4143  if (snc == -2)
4144  return divergence_type();
4145  else if (snc != -1)
4146  return fe_values->finite_element_output
4147  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4148  .single_nonzero_component_index];
4149  else
4150  {
4151  divergence_type return_value = 0;
4152  for (unsigned int d = 0; d < dim; ++d)
4153  if (shape_function_data[shape_function]
4154  .is_nonzero_shape_function_component[d])
4155  return_value +=
4156  fe_values->finite_element_output.shape_gradients
4157  [shape_function_data[shape_function].row_index[d]][q_point][d];
4158 
4159  return return_value;
4160  }
4161  }
4162 
4163 
4164 
4165  template <int dim, int spacedim>
4166  inline typename Vector<dim, spacedim>::curl_type
4167  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4168  const unsigned int q_point) const
4169  {
4170  // this function works like in the case above
4171 
4172  Assert(shape_function < fe_values->fe->dofs_per_cell,
4173  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4174  Assert(fe_values->update_flags & update_gradients,
4176  "update_gradients")));
4177  // same as for the scalar case except that we have one more index
4178  const int snc =
4179  shape_function_data[shape_function].single_nonzero_component;
4180 
4181  if (snc == -2)
4182  return curl_type();
4183 
4184  else
4185  switch (dim)
4186  {
4187  case 1:
4188  {
4189  Assert(false,
4190  ExcMessage(
4191  "Computing the curl in 1d is not a useful operation"));
4192  return curl_type();
4193  }
4194 
4195  case 2:
4196  {
4197  if (snc != -1)
4198  {
4199  curl_type return_value;
4200 
4201  // the single nonzero component can only be zero or one in 2d
4202  if (shape_function_data[shape_function]
4203  .single_nonzero_component_index == 0)
4204  return_value[0] =
4205  -1.0 * fe_values->finite_element_output
4206  .shape_gradients[snc][q_point][1];
4207  else
4208  return_value[0] = fe_values->finite_element_output
4209  .shape_gradients[snc][q_point][0];
4210 
4211  return return_value;
4212  }
4213 
4214  else
4215  {
4216  curl_type return_value;
4217 
4218  return_value[0] = 0.0;
4219 
4220  if (shape_function_data[shape_function]
4221  .is_nonzero_shape_function_component[0])
4222  return_value[0] -=
4223  fe_values->finite_element_output
4224  .shape_gradients[shape_function_data[shape_function]
4225  .row_index[0]][q_point][1];
4226 
4227  if (shape_function_data[shape_function]
4228  .is_nonzero_shape_function_component[1])
4229  return_value[0] +=
4230  fe_values->finite_element_output
4231  .shape_gradients[shape_function_data[shape_function]
4232  .row_index[1]][q_point][0];
4233 
4234  return return_value;
4235  }
4236  }
4237 
4238  case 3:
4239  {
4240  if (snc != -1)
4241  {
4242  curl_type return_value;
4243 
4244  switch (shape_function_data[shape_function]
4245  .single_nonzero_component_index)
4246  {
4247  case 0:
4248  {
4249  return_value[0] = 0;
4250  return_value[1] = fe_values->finite_element_output
4251  .shape_gradients[snc][q_point][2];
4252  return_value[2] =
4253  -1.0 * fe_values->finite_element_output
4254  .shape_gradients[snc][q_point][1];
4255  return return_value;
4256  }
4257 
4258  case 1:
4259  {
4260  return_value[0] =
4261  -1.0 * fe_values->finite_element_output
4262  .shape_gradients[snc][q_point][2];
4263  return_value[1] = 0;
4264  return_value[2] = fe_values->finite_element_output
4265  .shape_gradients[snc][q_point][0];
4266  return return_value;
4267  }
4268 
4269  default:
4270  {
4271  return_value[0] = fe_values->finite_element_output
4272  .shape_gradients[snc][q_point][1];
4273  return_value[1] =
4274  -1.0 * fe_values->finite_element_output
4275  .shape_gradients[snc][q_point][0];
4276  return_value[2] = 0;
4277  return return_value;
4278  }
4279  }
4280  }
4281 
4282  else
4283  {
4284  curl_type return_value;
4285 
4286  for (unsigned int i = 0; i < dim; ++i)
4287  return_value[i] = 0.0;
4288 
4289  if (shape_function_data[shape_function]
4290  .is_nonzero_shape_function_component[0])
4291  {
4292  return_value[1] +=
4293  fe_values->finite_element_output
4294  .shape_gradients[shape_function_data[shape_function]
4295  .row_index[0]][q_point][2];
4296  return_value[2] -=
4297  fe_values->finite_element_output
4298  .shape_gradients[shape_function_data[shape_function]
4299  .row_index[0]][q_point][1];
4300  }
4301 
4302  if (shape_function_data[shape_function]
4303  .is_nonzero_shape_function_component[1])
4304  {
4305  return_value[0] -=
4306  fe_values->finite_element_output
4307  .shape_gradients[shape_function_data[shape_function]
4308  .row_index[1]][q_point][2];
4309  return_value[2] +=
4310  fe_values->finite_element_output
4311  .shape_gradients[shape_function_data[shape_function]
4312  .row_index[1]][q_point][0];
4313  }
4314 
4315  if (shape_function_data[shape_function]
4316  .is_nonzero_shape_function_component[2])
4317  {
4318  return_value[0] +=
4319  fe_values->finite_element_output
4320  .shape_gradients[shape_function_data[shape_function]
4321  .row_index[2]][q_point][1];
4322  return_value[1] -=
4323  fe_values->finite_element_output
4324  .shape_gradients[shape_function_data[shape_function]
4325  .row_index[2]][q_point][0];
4326  }
4327 
4328  return return_value;
4329  }
4330  }
4331  }
4332  // should not end up here
4333  Assert(false, ExcInternalError());
4334  return curl_type();
4335  }
4336 
4337 
4338 
4339  template <int dim, int spacedim>
4340  inline typename Vector<dim, spacedim>::hessian_type
4341  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
4342  const unsigned int q_point) const
4343  {
4344  // this function works like in the case above
4345  Assert(shape_function < fe_values->fe->dofs_per_cell,
4346  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4347  Assert(fe_values->update_flags & update_hessians,
4349  "update_hessians")));
4350 
4351  // same as for the scalar case except that we have one more index
4352  const int snc =
4353  shape_function_data[shape_function].single_nonzero_component;
4354  if (snc == -2)
4355  return hessian_type();
4356  else if (snc != -1)
4357  {
4358  hessian_type return_value;
4359  return_value[shape_function_data[shape_function]
4360  .single_nonzero_component_index] =
4361  fe_values->finite_element_output.shape_hessians[snc][q_point];
4362  return return_value;
4363  }
4364  else
4365  {
4366  hessian_type return_value;
4367  for (unsigned int d = 0; d < dim; ++d)
4368  if (shape_function_data[shape_function]
4369  .is_nonzero_shape_function_component[d])
4370  return_value[d] =
4371  fe_values->finite_element_output.shape_hessians
4372  [shape_function_data[shape_function].row_index[d]][q_point];
4373 
4374  return return_value;
4375  }
4376  }
4377 
4378 
4379 
4380  template <int dim, int spacedim>
4381  inline typename Vector<dim, spacedim>::third_derivative_type
4382  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
4383  const unsigned int q_point) const
4384  {
4385  // this function works like in the case above
4386  Assert(shape_function < fe_values->fe->dofs_per_cell,
4387  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4388  Assert(fe_values->update_flags & update_3rd_derivatives,
4390  "update_3rd_derivatives")));
4391 
4392  // same as for the scalar case except that we have one more index
4393  const int snc =
4394  shape_function_data[shape_function].single_nonzero_component;
4395  if (snc == -2)
4396  return third_derivative_type();
4397  else if (snc != -1)
4398  {
4399  third_derivative_type return_value;
4400  return_value[shape_function_data[shape_function]
4401  .single_nonzero_component_index] =
4402  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
4403  return return_value;
4404  }
4405  else
4406  {
4407  third_derivative_type return_value;
4408  for (unsigned int d = 0; d < dim; ++d)
4409  if (shape_function_data[shape_function]
4410  .is_nonzero_shape_function_component[d])
4411  return_value[d] =
4412  fe_values->finite_element_output.shape_3rd_derivatives
4413  [shape_function_data[shape_function].row_index[d]][q_point];
4414 
4415  return return_value;
4416  }
4417  }
4418 
4419 
4420 
4421  namespace internal
4422  {
4427  inline ::SymmetricTensor<2, 1>
4428  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
4429  {
4430  Assert(n < 1, ExcIndexRange(n, 0, 1));
4431  (void)n;
4432 
4433  return {{t[0]}};
4434  }
4435 
4436 
4437 
4438  inline ::SymmetricTensor<2, 2>
4439  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
4440  {
4441  switch (n)
4442  {
4443  case 0:
4444  {
4445  return {{t[0], 0, t[1] / 2}};
4446  }
4447  case 1:
4448  {
4449  return {{0, t[1], t[0] / 2}};
4450  }
4451  default:
4452  {
4453  Assert(false, ExcIndexRange(n, 0, 2));
4454  return {};
4455  }
4456  }
4457  }
4458 
4459 
4460 
4461  inline ::SymmetricTensor<2, 3>
4462  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
4463  {
4464  switch (n)
4465  {
4466  case 0:
4467  {
4468  return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
4469  }
4470  case 1:
4471  {
4472  return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
4473  }
4474  case 2:
4475  {
4476  return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
4477  }
4478  default:
4479  {
4480  Assert(false, ExcIndexRange(n, 0, 3));
4481  return {};
4482  }
4483  }
4484  }
4485  } // namespace internal
4486 
4487 
4488 
4489  template <int dim, int spacedim>
4490  inline typename Vector<dim, spacedim>::symmetric_gradient_type
4491  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
4492  const unsigned int q_point) const
4493  {
4494  Assert(shape_function < fe_values->fe->dofs_per_cell,
4495  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4496  Assert(fe_values->update_flags & update_gradients,
4498  "update_gradients")));
4499 
4500  // same as for the scalar case except that we have one more index
4501  const int snc =
4502  shape_function_data[shape_function].single_nonzero_component;
4503  if (snc == -2)
4504  return symmetric_gradient_type();
4505  else if (snc != -1)
4506  return internal::symmetrize_single_row(
4507  shape_function_data[shape_function].single_nonzero_component_index,
4508  fe_values->finite_element_output.shape_gradients[snc][q_point]);
4509  else
4510  {
4511  gradient_type return_value;
4512  for (unsigned int d = 0; d < dim; ++d)
4513  if (shape_function_data[shape_function]
4514  .is_nonzero_shape_function_component[d])
4515  return_value[d] =
4516  fe_values->finite_element_output.shape_gradients
4517  [shape_function_data[shape_function].row_index[d]][q_point];
4518 
4519  return symmetrize(return_value);
4520  }
4521  }
4522 
4523 
4524 
4525  template <int dim, int spacedim>
4527  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
4528  const unsigned int q_point) const
4529  {
4530  Assert(shape_function < fe_values->fe->dofs_per_cell,
4531  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4532  Assert(fe_values->update_flags & update_values,
4534  "update_values")));
4535 
4536  // similar to the vector case where we have more then one index and we need
4537  // to convert between unrolled and component indexing for tensors
4538  const int snc =
4539  shape_function_data[shape_function].single_nonzero_component;
4540 
4541  if (snc == -2)
4542  {
4543  // shape function is zero for the selected components
4544  return value_type();
4545  }
4546  else if (snc != -1)
4547  {
4548  value_type return_value;
4549  const unsigned int comp =
4550  shape_function_data[shape_function].single_nonzero_component_index;
4551  return_value[value_type::unrolled_to_component_indices(comp)] =
4552  fe_values->finite_element_output.shape_values(snc, q_point);
4553  return return_value;
4554  }
4555  else
4556  {
4557  value_type return_value;
4558  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
4559  if (shape_function_data[shape_function]
4560  .is_nonzero_shape_function_component[d])
4561  return_value[value_type::unrolled_to_component_indices(d)] =
4562  fe_values->finite_element_output.shape_values(
4563  shape_function_data[shape_function].row_index[d], q_point);
4564  return return_value;
4565  }
4566  }
4567 
4568 
4569 
4570  template <int dim, int spacedim>
4573  const unsigned int shape_function,
4574  const unsigned int q_point) const
4575  {
4576  Assert(shape_function < fe_values->fe->dofs_per_cell,
4577  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4578  Assert(fe_values->update_flags & update_gradients,
4580  "update_gradients")));
4581 
4582  const int snc =
4583  shape_function_data[shape_function].single_nonzero_component;
4584 
4585  if (snc == -2)
4586  {
4587  // shape function is zero for the selected components
4588  return divergence_type();
4589  }
4590  else if (snc != -1)
4591  {
4592  // we have a single non-zero component when the symmetric tensor is
4593  // represented in unrolled form. this implies we potentially have
4594  // two non-zero components when represented in component form! we
4595  // will only have one non-zero entry if the non-zero component lies on
4596  // the diagonal of the tensor.
4597  //
4598  // the divergence of a second-order tensor is a first order tensor.
4599  //
4600  // assume the second-order tensor is A with components A_{ij}. then
4601  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
4602  // entries in the tensorial representation. define the
4603  // divergence as:
4604  // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
4605  // (which is incidentally also
4606  // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
4607  // In both cases, a sum is implied.
4608  //
4609  // Now, we know the nonzero component in unrolled form: it is indicated
4610  // by 'snc'. we can figure out which tensor components belong to this:
4611  const unsigned int comp =
4612  shape_function_data[shape_function].single_nonzero_component_index;
4613  const unsigned int ii =
4614  value_type::unrolled_to_component_indices(comp)[0];
4615  const unsigned int jj =
4616  value_type::unrolled_to_component_indices(comp)[1];
4617 
4618  // given the form of the divergence above, if ii=jj there is only a
4619  // single nonzero component of the full tensor and the gradient
4620  // equals
4621  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
4622  // all other entries of 'b' are zero
4623  //
4624  // on the other hand, if ii!=jj, then there are two nonzero entries in
4625  // the full tensor and
4626  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
4627  // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
4628  // again, all other entries of 'b' are zero
4629  const ::Tensor<1, spacedim> &phi_grad =
4630  fe_values->finite_element_output.shape_gradients[snc][q_point];
4631 
4632  divergence_type return_value;
4633  return_value[ii] = phi_grad[jj];
4634 
4635  if (ii != jj)
4636  return_value[jj] = phi_grad[ii];
4637 
4638  return return_value;
4639  }
4640  else
4641  {
4642  Assert(false, ExcNotImplemented());
4643  divergence_type return_value;
4644  return return_value;
4645  }
4646  }
4647 
4648 
4649 
4650  template <int dim, int spacedim>
4651  inline typename Tensor<2, dim, spacedim>::value_type
4652  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
4653  const unsigned int q_point) const
4654  {
4655  Assert(shape_function < fe_values->fe->dofs_per_cell,
4656  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4657  Assert(fe_values->update_flags & update_values,
4659  "update_values")));
4660 
4661  // similar to the vector case where we have more then one index and we need
4662  // to convert between unrolled and component indexing for tensors
4663  const int snc =
4664  shape_function_data[shape_function].single_nonzero_component;
4665 
4666  if (snc == -2)
4667  {
4668  // shape function is zero for the selected components
4669  return value_type();
4670  }
4671  else if (snc != -1)
4672  {
4673  value_type return_value;
4674  const unsigned int comp =
4675  shape_function_data[shape_function].single_nonzero_component_index;
4676  const TableIndices<2> indices =
4678  return_value[indices] =
4679  fe_values->finite_element_output.shape_values(snc, q_point);
4680  return return_value;
4681  }
4682  else
4683  {
4684  value_type return_value;
4685  for (unsigned int d = 0; d < dim * dim; ++d)
4686  if (shape_function_data[shape_function]
4687  .is_nonzero_shape_function_component[d])
4688  {
4689  const TableIndices<2> indices =
4691  return_value[indices] =
4692  fe_values->finite_element_output.shape_values(
4693  shape_function_data[shape_function].row_index[d], q_point);
4694  }
4695  return return_value;
4696  }
4697  }
4698 
4699 
4700 
4701  template <int dim, int spacedim>
4703  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
4704  const unsigned int q_point) const
4705  {
4706  Assert(shape_function < fe_values->fe->dofs_per_cell,
4707  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4708  Assert(fe_values->update_flags & update_gradients,
4710  "update_gradients")));
4711 
4712  const int snc =
4713  shape_function_data[shape_function].single_nonzero_component;
4714 
4715  if (snc == -2)
4716  {
4717  // shape function is zero for the selected components
4718  return divergence_type();
4719  }
4720  else if (snc != -1)
4721  {
4722  // we have a single non-zero component when the tensor is
4723  // represented in unrolled form.
4724  //
4725  // the divergence of a second-order tensor is a first order tensor.
4726  //
4727  // assume the second-order tensor is A with components A_{ij},
4728  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
4729  //
4730  // Now, we know the nonzero component in unrolled form: it is indicated
4731  // by 'snc'. we can figure out which tensor components belong to this:
4732  const unsigned int comp =
4733  shape_function_data[shape_function].single_nonzero_component_index;
4734  const TableIndices<2> indices =
4736  const unsigned int ii = indices[0];
4737  const unsigned int jj = indices[1];
4738 
4739  const ::Tensor<1, spacedim> &phi_grad =
4740  fe_values->finite_element_output.shape_gradients[snc][q_point];
4741 
4742  divergence_type return_value;
4743  // note that we contract \nabla from the right
4744  return_value[ii] = phi_grad[jj];
4745 
4746  return return_value;
4747  }
4748  else
4749  {
4750  Assert(false, ExcNotImplemented());
4751  divergence_type return_value;
4752  return return_value;
4753  }
4754  }
4755 
4756 
4757 
4758  template <int dim, int spacedim>
4760  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
4761  const unsigned int q_point) const
4762  {
4763  Assert(shape_function < fe_values->fe->dofs_per_cell,
4764  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4765  Assert(fe_values->update_flags & update_gradients,
4767  "update_gradients")));
4768 
4769  const int snc =
4770  shape_function_data[shape_function].single_nonzero_component;
4771 
4772  if (snc == -2)
4773  {
4774  // shape function is zero for the selected components
4775  return gradient_type();
4776  }
4777  else if (snc != -1)
4778  {
4779  // we have a single non-zero component when the tensor is
4780  // represented in unrolled form.
4781  //
4782  // the gradient of a second-order tensor is a third order tensor.
4783  //
4784  // assume the second-order tensor is A with components A_{ij},
4785  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
4786  //
4787  // Now, we know the nonzero component in unrolled form: it is indicated
4788  // by 'snc'. we can figure out which tensor components belong to this:
4789  const unsigned int comp =
4790  shape_function_data[shape_function].single_nonzero_component_index;
4791  const TableIndices<2> indices =
4793  const unsigned int ii = indices[0];
4794  const unsigned int jj = indices[1];
4795 
4796  const ::Tensor<1, spacedim> &phi_grad =
4797  fe_values->finite_element_output.shape_gradients[snc][q_point];
4798 
4799  gradient_type return_value;
4800  return_value[ii][jj] = phi_grad;
4801 
4802  return return_value;
4803  }
4804  else
4805  {
4806  Assert(false, ExcNotImplemented());
4807  gradient_type return_value;
4808  return return_value;
4809  }
4810  }
4811 
4812 } // namespace FEValuesViews
4813 
4814 
4815 
4816 /*---------------------- Inline functions: FEValuesBase ---------------------*/
4817 
4818 
4819 
4820 template <int dim, int spacedim>
4822  operator[](const FEValuesExtractors::Scalar &scalar) const
4823 {
4824  Assert(scalar.component < fe_values_views_cache.scalars.size(),
4825  ExcIndexRange(scalar.component,
4826  0,
4827  fe_values_views_cache.scalars.size()));
4828 
4829  return fe_values_views_cache.scalars[scalar.component];
4830 }
4831 
4832 
4833 
4834 template <int dim, int spacedim>
4836  operator[](const FEValuesExtractors::Vector &vector) const
4837 {
4838  Assert(vector.first_vector_component < fe_values_views_cache.vectors.size(),
4840  0,
4841  fe_values_views_cache.vectors.size()));
4842 
4843  return fe_values_views_cache.vectors[vector.first_vector_component];
4844 }
4845 
4846 
4847 
4848 template <int dim, int spacedim>
4852 {
4853  Assert(
4854  tensor.first_tensor_component <
4855  fe_values_views_cache.symmetric_second_order_tensors.size(),
4857  0,
4858  fe_values_views_cache.symmetric_second_order_tensors.size()));
4859 
4860  return fe_values_views_cache
4861  .symmetric_second_order_tensors[tensor.first_tensor_component];
4862 }
4863 
4864 
4865 
4866 template <int dim, int spacedim>
4869  operator[](const FEValuesExtractors::Tensor<2> &tensor) const
4870 {
4872  fe_values_views_cache.second_order_tensors.size(),
4874  0,
4875  fe_values_views_cache.second_order_tensors.size()));
4876 
4877  return fe_values_views_cache
4878  .second_order_tensors[tensor.first_tensor_component];
4879 }
4880 
4881 
4882 
4883 template <int dim, int spacedim>
4884 inline const double &
4885 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
4886  const unsigned int j) const
4887 {
4888  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4889  Assert(this->update_flags & update_values,
4890  ExcAccessToUninitializedField("update_values"));
4891  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4892  Assert(present_cell.get() != nullptr,
4893  ExcMessage("FEValues object is not reinit'ed to any cell"));
4894  // if the entire FE is primitive,
4895  // then we can take a short-cut:
4896  if (fe->is_primitive())
4897  return this->finite_element_output.shape_values(i, j);
4898  else
4899  {
4900  // otherwise, use the mapping
4901  // between shape function
4902  // numbers and rows. note that
4903  // by the assertions above, we
4904  // know that this particular
4905  // shape function is primitive,
4906  // so we can call
4907  // system_to_component_index
4908  const unsigned int row =
4909  this->finite_element_output
4910  .shape_function_to_row_table[i * fe->n_components() +
4911  fe->system_to_component_index(i).first];
4912  return this->finite_element_output.shape_values(row, j);
4913  }
4914 }
4915 
4916 
4917 
4918 template <int dim, int spacedim>
4919 inline double
4921  const unsigned int i,
4922  const unsigned int j,
4923  const unsigned int component) const
4924 {
4925  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4926  Assert(this->update_flags & update_values,
4927  ExcAccessToUninitializedField("update_values"));
4928  Assert(component < fe->n_components(),
4929  ExcIndexRange(component, 0, fe->n_components()));
4930  Assert(present_cell.get() != nullptr,
4931  ExcMessage("FEValues object is not reinit'ed to any cell"));
4932 
4933  // check whether the shape function
4934  // is non-zero at all within
4935  // this component:
4936  if (fe->get_nonzero_components(i)[component] == false)
4937  return 0;
4938 
4939  // look up the right row in the
4940  // table and take the data from
4941  // there
4942  const unsigned int row =
4943  this->finite_element_output
4944  .shape_function_to_row_table[i * fe->n_components() + component];
4945  return this->finite_element_output.shape_values(row, j);
4946 }
4947 
4948 
4949 
4950 template <int dim, int spacedim>
4951 inline const Tensor<1, spacedim> &
4952 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
4953  const unsigned int j) const
4954 {
4955  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4956  Assert(this->update_flags & update_gradients,
4957  ExcAccessToUninitializedField("update_gradients"));
4958  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4959  Assert(present_cell.get() != nullptr,
4960  ExcMessage("FEValues object is not reinit'ed to any cell"));
4961  // if the entire FE is primitive,
4962  // then we can take a short-cut:
4963  if (fe->is_primitive())
4964  return this->finite_element_output.shape_gradients[i][j];
4965  else
4966  {
4967  // otherwise, use the mapping
4968  // between shape function
4969  // numbers and rows. note that
4970  // by the assertions above, we
4971  // know that this particular
4972  // shape function is primitive,
4973  // so we can call
4974  // system_to_component_index
4975  const unsigned int row =
4976  this->finite_element_output
4977  .shape_function_to_row_table[i * fe->n_components() +
4978  fe->system_to_component_index(i).first];
4979  return this->finite_element_output.shape_gradients[row][j];
4980  }
4981 }
4982 
4983 
4984 
4985 template <int dim, int spacedim>
4986 inline Tensor<1, spacedim>
4988  const unsigned int i,
4989  const unsigned int j,
4990  const unsigned int component) const
4991 {
4992  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4993  Assert(this->update_flags & update_gradients,
4994  ExcAccessToUninitializedField("update_gradients"));
4995  Assert(component < fe->n_components(),
4996  ExcIndexRange(component, 0, fe->n_components()));
4997  Assert(present_cell.get() != nullptr,
4998  ExcMessage("FEValues object is not reinit'ed to any cell"));
4999  // check whether the shape function
5000  // is non-zero at all within
5001  // this component:
5002  if (fe->get_nonzero_components(i)[component] == false)
5003  return Tensor<1, spacedim>();
5004 
5005  // look up the right row in the
5006  // table and take the data from
5007  // there
5008  const unsigned int row =
5009  this->finite_element_output
5010  .shape_function_to_row_table[i * fe->n_components() + component];
5011  return this->finite_element_output.shape_gradients[row][j];
5012 }
5013 
5014 
5015 
5016 template <int dim, int spacedim>
5017 inline const Tensor<2, spacedim> &
5018 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
5019  const unsigned int j) const
5020 {
5021  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5022  Assert(this->update_flags & update_hessians,
5023  ExcAccessToUninitializedField("update_hessians"));
5024  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5025  Assert(present_cell.get() != nullptr,
5026  ExcMessage("FEValues object is not reinit'ed to any cell"));
5027  // if the entire FE is primitive,
5028  // then we can take a short-cut:
5029  if (fe->is_primitive())
5030  return this->finite_element_output.shape_hessians[i][j];
5031  else
5032  {
5033  // otherwise, use the mapping
5034  // between shape function
5035  // numbers and rows. note that
5036  // by the assertions above, we
5037  // know that this particular
5038  // shape function is primitive,
5039  // so we can call
5040  // system_to_component_index
5041  const unsigned int row =
5042  this->finite_element_output
5043  .shape_function_to_row_table[i * fe->n_components() +
5044  fe->system_to_component_index(i).first];
5045  return this->finite_element_output.shape_hessians[row][j];
5046  }
5047 }
5048 
5049 
5050 
5051 template <int dim, int spacedim>
5052 inline Tensor<2, spacedim>
5054  const unsigned int i,
5055  const unsigned int j,
5056  const unsigned int component) const
5057 {
5058  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5059  Assert(this->update_flags & update_hessians,
5060  ExcAccessToUninitializedField("update_hessians"));
5061  Assert(component < fe->n_components(),
5062  ExcIndexRange(component, 0, fe->n_components()));
5063  Assert(present_cell.get() != nullptr,
5064  ExcMessage("FEValues object is not reinit'ed to any cell"));
5065  // check whether the shape function
5066  // is non-zero at all within
5067  // this component:
5068  if (fe->get_nonzero_components(i)[component] == false)
5069  return Tensor<2, spacedim>();
5070 
5071  // look up the right row in the
5072  // table and take the data from
5073  // there
5074  const unsigned int row =
5075  this->finite_element_output
5076  .shape_function_to_row_table[i * fe->n_components() + component];
5077  return this->finite_element_output.shape_hessians[row][j];
5078 }
5079 
5080 
5081 
5082 template <int dim, int spacedim>
5083 inline const Tensor<3, spacedim> &
5085  const unsigned int j) const
5086 {
5087  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5088  Assert(this->update_flags & update_hessians,
5089  ExcAccessToUninitializedField("update_3rd_derivatives"));
5090  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5091  Assert(present_cell.get() != nullptr,
5092  ExcMessage("FEValues object is not reinit'ed to any cell"));
5093  // if the entire FE is primitive,
5094  // then we can take a short-cut:
5095  if (fe->is_primitive())
5096  return this->finite_element_output.shape_3rd_derivatives[i][j];
5097  else
5098  {
5099  // otherwise, use the mapping
5100  // between shape function
5101  // numbers and rows. note that
5102  // by the assertions above, we
5103  // know that this particular
5104  // shape function is primitive,
5105  // so we can call
5106  // system_to_component_index
5107  const unsigned int row =
5108  this->finite_element_output
5109  .shape_function_to_row_table[i * fe->n_components() +
5110  fe->system_to_component_index(i).first];
5111  return this->finite_element_output.shape_3rd_derivatives[row][j];
5112  }
5113 }
5114 
5115 
5116 
5117 template <int dim, int spacedim>
5118 inline Tensor<3, spacedim>
5120  const unsigned int i,
5121  const unsigned int j,
5122  const unsigned int component) const
5123 {
5124  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5125  Assert(this->update_flags & update_hessians,
5126  ExcAccessToUninitializedField("update_3rd_derivatives"));
5127  Assert(component < fe->n_components(),
5128  ExcIndexRange(component, 0, fe->n_components()));
5129  Assert(present_cell.get() != nullptr,
5130  ExcMessage("FEValues object is not reinit'ed to any cell"));
5131  // check whether the shape function
5132  // is non-zero at all within
5133  // this component:
5134  if (fe->get_nonzero_components(i)[component] == false)
5135  return Tensor<3, spacedim>();
5136 
5137  // look up the right row in the
5138  // table and take the data from
5139  // there
5140  const unsigned int row =
5141  this->finite_element_output
5142  .shape_function_to_row_table[i * fe->n_components() + component];
5143  return this->finite_element_output.shape_3rd_derivatives[row][j];
5144 }
5145 
5146 
5147 
5148 template <int dim, int spacedim>
5149 inline const FiniteElement<dim, spacedim> &
5151 {
5152  return *fe;
5153 }
5154 
5155 
5156 
5157 template <int dim, int spacedim>
5158 inline const Mapping<dim, spacedim> &
5160 {
5161  return *mapping;
5162 }
5163 
5164 
5165 
5166 template <int dim, int spacedim>
5167 inline UpdateFlags
5169 {
5170  return this->update_flags;
5171 }
5172 
5173 
5174 
5175 template <int dim, int spacedim>
5176 inline const std::vector<Point<spacedim>> &
5178 {
5179  Assert(this->update_flags & update_quadrature_points,
5180  ExcAccessToUninitializedField("update_quadrature_points"));
5181  Assert(present_cell.get() != nullptr,
5182  ExcMessage("FEValues object is not reinit'ed to any cell"));
5183  return this->mapping_output.quadrature_points;
5184 }
5185 
5186 
5187 
5188 template <int dim, int spacedim>
5189 inline const std::vector<double> &
5191 {
5192  Assert(this->update_flags & update_JxW_values,
5193  ExcAccessToUninitializedField("update_JxW_values"));
5194  Assert(present_cell.get() != nullptr,
5195  ExcMessage("FEValues object is not reinit'ed to any cell"));
5196  return this->mapping_output.JxW_values;
5197 }
5198 
5199 
5200 
5201 template <int dim, int spacedim>
5202 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5204 {
5205  Assert(this->update_flags & update_jacobians,
5206  ExcAccessToUninitializedField("update_jacobians"));
5207  Assert(present_cell.get() != nullptr,
5208  ExcMessage("FEValues object is not reinit'ed to any cell"));
5209  return this->mapping_output.jacobians;
5210 }
5211 
5212 
5213 
5214 template <int dim, int spacedim>
5215 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5217 {
5218  Assert(this->update_flags & update_jacobian_grads,
5219  ExcAccessToUninitializedField("update_jacobians_grads"));
5220  Assert(present_cell.get() != nullptr,
5221  ExcMessage("FEValues object is not reinit'ed to any cell"));
5222  return this->mapping_output.jacobian_grads;
5223 }
5224 
5225 
5226 
5227 template <int dim, int spacedim>
5228 inline const Tensor<3, spacedim> &
5230  const unsigned int i) const
5231 {
5232  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5233  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5234  Assert(present_cell.get() != nullptr,
5235  ExcMessage("FEValues object is not reinit'ed to any cell"));
5236  return this->mapping_output.jacobian_pushed_forward_grads[i];
5237 }
5238 
5239 
5240 
5241 template <int dim, int spacedim>
5242 inline const std::vector<Tensor<3, spacedim>> &
5244 {
5245  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5246  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5247  Assert(present_cell.get() != nullptr,
5248  ExcMessage("FEValues object is not reinit'ed to any cell"));
5249  return this->mapping_output.jacobian_pushed_forward_grads;
5250 }
5251 
5252 
5253 
5254 template <int dim, int spacedim>
5255 inline const DerivativeForm<3, dim, spacedim> &
5256 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5257 {
5258  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5259  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5260  Assert(present_cell.get() != nullptr,
5261  ExcMessage("FEValues object is not reinit'ed to any cell"));
5262  return this->mapping_output.jacobian_2nd_derivatives[i];
5263 }
5264 
5265 
5266 
5267 template <int dim, int spacedim>
5268 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5270 {
5271  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5272  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5273  Assert(present_cell.get() != nullptr,
5274  ExcMessage("FEValues object is not reinit'ed to any cell"));
5275  return this->mapping_output.jacobian_2nd_derivatives;
5276 }
5277 
5278 
5279 
5280 template <int dim, int spacedim>
5281 inline const Tensor<4, spacedim> &
5283  const unsigned int i) const
5284 {
5287  "update_jacobian_pushed_forward_2nd_derivatives"));
5288  Assert(present_cell.get() != nullptr,
5289  ExcMessage("FEValues object is not reinit'ed to any cell"));
5290  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5291 }
5292 
5293 
5294 
5295 template <int dim, int spacedim>
5296 inline const std::vector<Tensor<4, spacedim>> &
5298 {
5299  Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
5301  "update_jacobian_pushed_forward_2nd_derivatives"));
5302  Assert(present_cell.get() != nullptr,
5303  ExcMessage("FEValues object is not reinit'ed to any cell"));
5304  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5305 }
5306 
5307 
5308 
5309 template <int dim, int spacedim>
5310 inline const DerivativeForm<4, dim, spacedim> &
5311 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5312 {
5313  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5314  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5315  Assert(present_cell.get() != nullptr,
5316  ExcMessage("FEValues object is not reinit'ed to any cell"));
5317  return this->mapping_output.jacobian_3rd_derivatives[i];
5318 }
5319 
5320 
5321 
5322 template <int dim, int spacedim>
5323 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5325 {
5326  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5327  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5328  Assert(present_cell.get() != nullptr,
5329  ExcMessage("FEValues object is not reinit'ed to any cell"));
5330  return this->mapping_output.jacobian_3rd_derivatives;
5331 }
5332 
5333 
5334 
5335 template <int dim, int spacedim>
5336 inline const Tensor<5, spacedim> &
5338  const unsigned int i) const
5339 {
5342  "update_jacobian_pushed_forward_3rd_derivatives"));
5343  Assert(present_cell.get() != nullptr,
5344  ExcMessage("FEValues object is not reinit'ed to any cell"));
5345  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
5346 }
5347 
5348 
5349 
5350 template <int dim, int spacedim>
5351 inline const std::vector<Tensor<5, spacedim>> &
5353 {
5354  Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
5356  "update_jacobian_pushed_forward_3rd_derivatives"));
5357  Assert(present_cell.get() != nullptr,
5358  ExcMessage("FEValues object is not reinit'ed to any cell"));
5359  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
5360 }
5361 
5362 
5363 
5364 template <int dim, int spacedim>
5365 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
5367 {
5368  Assert(this->update_flags & update_inverse_jacobians,
5369  ExcAccessToUninitializedField("update_inverse_jacobians"));
5370  Assert(present_cell.get() != nullptr,
5371  ExcMessage("FEValues object is not reinit'ed to any cell"));
5372  return this->mapping_output.inverse_jacobians;
5373 }
5374 
5375 
5376 
5377 template <int dim, int spacedim>
5378 inline const Point<spacedim> &
5379 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
5380 {
5381  Assert(this->update_flags & update_quadrature_points,
5382  ExcAccessToUninitializedField("update_quadrature_points"));
5383  Assert(i < this->mapping_output.quadrature_points.size(),
5384  ExcIndexRange(i, 0, this->mapping_output.quadrature_points.size()));
5385  Assert(present_cell.get() != nullptr,
5386  ExcMessage("FEValues object is not reinit'ed to any cell"));
5387 
5388  return this->mapping_output.quadrature_points[i];
5389 }
5390 
5391 
5392 
5393 template <int dim, int spacedim>
5394 inline double
5395 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
5396 {
5397  Assert(this->update_flags & update_JxW_values,
5398  ExcAccessToUninitializedField("update_JxW_values"));
5399  Assert(i < this->mapping_output.JxW_values.size(),
5400  ExcIndexRange(i, 0, this->mapping_output.JxW_values.size()));
5401  Assert(present_cell.get() != nullptr,
5402  ExcMessage("FEValues object is not reinit'ed to any cell"));
5403 
5404  return this->mapping_output.JxW_values[i];
5405 }
5406 
5407 
5408 
5409 template <int dim, int spacedim>
5410 inline const DerivativeForm<1, dim, spacedim> &
5411 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
5412 {
5413  Assert(this->update_flags & update_jacobians,
5414  ExcAccessToUninitializedField("update_jacobians"));
5415  Assert(i < this->mapping_output.jacobians.size(),
5416  ExcIndexRange(i, 0, this->mapping_output.jacobians.size()));
5417  Assert(present_cell.get() != nullptr,
5418  ExcMessage("FEValues object is not reinit'ed to any cell"));
5419 
5420  return this->mapping_output.jacobians[i];
5421 }
5422 
5423 
5424 
5425 template <int dim, int spacedim>
5426 inline const DerivativeForm<2, dim, spacedim> &
5427 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
5428 {
5429  Assert(this->update_flags & update_jacobian_grads,
5430  ExcAccessToUninitializedField("update_jacobians_grads"));
5431  Assert(i < this->mapping_output.jacobian_grads.size(),
5432  ExcIndexRange(i, 0, this->mapping_output.jacobian_grads.size()));
5433  Assert(present_cell.get() != nullptr,
5434  ExcMessage("FEValues object is not reinit'ed to any cell"));
5435 
5436  return this->mapping_output.jacobian_grads[i];
5437 }
5438 
5439 
5440 
5441 template <int dim, int spacedim>
5442 inline const DerivativeForm<1, spacedim, dim> &
5443 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
5444 {
5445  Assert(this->update_flags & update_inverse_jacobians,
5446  ExcAccessToUninitializedField("update_inverse_jacobians"));
5447  Assert(i < this->mapping_output.inverse_jacobians.size(),
5448  ExcIndexRange(i, 0, this->mapping_output.inverse_jacobians.size()));
5449  Assert(present_cell.get() != nullptr,
5450  ExcMessage("FEValues object is not reinit'ed to any cell"));
5451 
5452  return this->mapping_output.inverse_jacobians[i];
5453 }
5454 
5455 
5456 
5457 template <int dim, int spacedim>
5458 inline const Tensor<1, spacedim> &
5459 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
5460 {
5461  Assert(this->update_flags & update_normal_vectors,
5463  "update_normal_vectors")));
5464  Assert(i < this->mapping_output.normal_vectors.size(),
5465  ExcIndexRange(i, 0, this->mapping_output.normal_vectors.size()));
5466  Assert(present_cell.get() != nullptr,
5467  ExcMessage("FEValues object is not reinit'ed to any cell"));
5468 
5469  return this->mapping_output.normal_vectors[i];
5470 }
5471 
5472 
5473 
5474 /*--------------------- Inline functions: FEValues --------------------------*/
5475 
5476 
5477 template <int dim, int spacedim>
5478 inline const Quadrature<dim> &
5480 {
5481  return quadrature;
5482 }
5483 
5484 
5485 
5486 template <int dim, int spacedim>
5487 inline const FEValues<dim, spacedim> &
5489 {
5490  return *this;
5491 }
5492 
5493 
5494 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
5495 
5496 
5497 template <int dim, int spacedim>
5498 inline unsigned int
5500 {
5501  return present_face_index;
5502 }
5503 
5504 
5505 /*----------------------- Inline functions: FE*FaceValues -------------------*/
5506 
5507 template <int dim, int spacedim>
5508 inline const Quadrature<dim - 1> &
5510 {
5511  return quadrature;
5512 }
5513 
5514 
5515 
5516 template <int dim, int spacedim>
5517 inline const FEFaceValues<dim, spacedim> &
5519 {
5520  return *this;
5521 }
5522 
5523 
5524 
5525 template <int dim, int spacedim>
5526 inline const FESubfaceValues<dim, spacedim> &
5528 {
5529  return *this;
5530 }
5531 
5532 
5533 
5534 template <int dim, int spacedim>
5535 inline const Tensor<1, spacedim> &
5536 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
5537 {
5538  Assert(i < this->mapping_output.boundary_forms.size(),
5539  ExcIndexRange(i, 0, this->mapping_output.boundary_forms.size()));
5540  Assert(this->update_flags & update_boundary_forms,
5542  "update_boundary_forms")));
5543 
5544  return this->mapping_output.boundary_forms[i];
5545 }
5546 
5547 #endif // DOXYGEN
5548 
5549 DEAL_II_NAMESPACE_CLOSE
5550 
5551 #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:3425
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:3457
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:3690
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1520
static ::ExceptionBase & ExcAccessToUninitializedField()
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:536
const unsigned int dofs_per_cell
Definition: fe_values.h:2109
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:3695
Volume element.
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:3341
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:3393
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:3472
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:1407
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:3588
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:3401
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int quadrature_point) const
typename Tensor< rank_ - 1, dim, Number >::tensor_type value_type
Definition: tensor.h:452
const std::vector< DerivativeForm< 4, dim, spacedim > > & get_jacobian_3rd_derivatives() const
Second derivatives of shape functions.
Gradient of volume element.
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1302
std::vector<::FEValuesViews::Scalar< dim, spacedim > > scalars
Definition: fe_values.h:1954
const Tensor< 1, spacedim > & boundary_form(const unsigned int i) const
const Quadrature< dim > & get_quadrature() const
const unsigned int n_quadrature_points
Definition: fe_values.h:2102
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:3366
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:656
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
Definition: tensor.h:422
const std::vector< double > & get_JxW_values() const
double JxW(const unsigned int quadrature_point) const
typename ProductType< Number, typename Vector< dim, spacedim >::symmetric_gradient_type >::type symmetric_gradient_type
Definition: fe_values.h:672
typename ProductType< Number, typename Vector< dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:680
typename ProductType< Number, typename Vector< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:704
Shape function gradients.
Normal vectors.
const std::vector< DerivativeForm< 2, dim, spacedim > > & get_jacobian_grads() const
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1589
const Tensor< 4, spacedim > & jacobian_pushed_forward_2nd_derivative(const unsigned int quadrature_point) const
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1856
Definition: fe.h:38
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1509
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int quadrature_point) const
static ::ExceptionBase & ExcNotImplemented()
const std::vector< DerivativeForm< 3, dim, spacedim > > & get_jacobian_2nd_derivatives() const
static DEAL_II_CONSTEXPR TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
Definition: tensor.h:1515
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1216
DEAL_II_CONSTEXPR SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
const Tensor< 1, spacedim > & shape_grad(const unsigned int function_no, const unsigned int quadrature_point) const
boost::signals2::connection tria_listener_refinement
Definition: fe_values.h:3357
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3433
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:3408
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:3439
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:3417
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