Reference documentation for deal.II version Git cd5f69ed86 2019-12-07 22:13:35 -0500
\(\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 
3134  const Tensor<1, spacedim> &
3135  normal_vector(const unsigned int i) const;
3136 
3147  DEAL_II_DEPRECATED
3148  const std::vector<Tensor<1, spacedim>> &
3149  get_all_normal_vectors() const;
3150 
3158  const std::vector<Tensor<1, spacedim>> &
3159  get_normal_vectors() const;
3160 
3162 
3164 
3165 
3175  operator[](const FEValuesExtractors::Scalar &scalar) const;
3176 
3186  operator[](const FEValuesExtractors::Vector &vector) const;
3187 
3198  operator[](const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
3199 
3200 
3210  operator[](const FEValuesExtractors::Tensor<2> &tensor) const;
3211 
3213 
3215 
3216 
3220  const Mapping<dim, spacedim> &
3221  get_mapping() const;
3222 
3227  get_fe() const;
3228 
3232  UpdateFlags
3233  get_update_flags() const;
3234 
3239  get_cell() const;
3240 
3247  get_cell_similarity() const;
3248 
3253  std::size_t
3254  memory_consumption() const;
3256 
3257 
3266  std::string,
3267  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3268  << "object for which this kind of information has not been computed. What "
3269  << "information these objects compute is determined by the update_* flags you "
3270  << "pass to the constructor. Here, the operation you are attempting requires "
3271  << "the <" << arg1
3272  << "> flag to be set, but it was apparently not specified "
3273  << "upon construction.");
3274 
3282  ExcFEDontMatch,
3283  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3284  "to the DoFHandler that provided the cell iterator do not match.");
3290  DeclException1(ExcShapeFunctionNotPrimitive,
3291  int,
3292  << "The shape function with index " << arg1
3293  << " is not primitive, i.e. it is vector-valued and "
3294  << "has more than one non-zero vector component. This "
3295  << "function cannot be called for these shape functions. "
3296  << "Maybe you want to use the same function with the "
3297  << "_component suffix?");
3298 
3307  "The given FiniteElement is not a primitive element but the requested operation "
3308  "only works for those. See FiniteElement::is_primitive() for more information.");
3309 
3310 protected:
3339  class CellIteratorBase;
3340 
3345  template <typename CI>
3347  class TriaCellIterator;
3348 
3354  std::unique_ptr<const CellIteratorBase> present_cell;
3355 
3363  boost::signals2::connection tria_listener_refinement;
3364 
3372  boost::signals2::connection tria_listener_mesh_transform;
3373 
3379  void
3380  invalidate_present_cell();
3381 
3391  void
3392  maybe_invalidate_previous_present_cell(
3393  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3394 
3400 
3406  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3408 
3415 
3416 
3424 
3430  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3432 
3438  spacedim>
3440 
3441 
3446 
3455  UpdateFlags
3456  compute_update_flags(const UpdateFlags update_flags) const;
3457 
3464 
3470  void
3471  check_cell_similarity(
3472  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3473 
3474 private:
3479 
3480  // Make the view classes friends of this class, since they access internal
3481  // data.
3482  template <int, int>
3483  friend class FEValuesViews::Scalar;
3484  template <int, int>
3485  friend class FEValuesViews::Vector;
3486  template <int, int, int>
3487  friend class FEValuesViews::SymmetricTensor;
3488  template <int, int, int>
3489  friend class FEValuesViews::Tensor;
3490 };
3491 
3492 
3493 
3504 template <int dim, int spacedim = dim>
3505 class FEValues : public FEValuesBase<dim, spacedim>
3506 {
3507 public:
3512  static const unsigned int integral_dimension = dim;
3513 
3518  FEValues(const Mapping<dim, spacedim> & mapping,
3519  const FiniteElement<dim, spacedim> &fe,
3520  const Quadrature<dim> & quadrature,
3521  const UpdateFlags update_flags);
3522 
3529  const Quadrature<dim> & quadrature,
3530  const UpdateFlags update_flags);
3531 
3538  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3539  void
3540  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3541  level_dof_access>> &cell);
3542 
3556  void
3557  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3558 
3563  const Quadrature<dim> &
3564  get_quadrature() const;
3565 
3570  std::size_t
3571  memory_consumption() const;
3572 
3587  const FEValues<dim, spacedim> &
3588  get_present_fe_values() const;
3589 
3590 private:
3595 
3599  void
3600  initialize(const UpdateFlags update_flags);
3601 
3608  void
3609  do_reinit();
3610 };
3611 
3612 
3623 template <int dim, int spacedim = dim>
3624 class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
3625 {
3626 public:
3631  static const unsigned int integral_dimension = dim - 1;
3632 
3644  FEFaceValuesBase(const unsigned int n_q_points,
3645  const unsigned int dofs_per_cell,
3646  const UpdateFlags update_flags,
3647  const Mapping<dim, spacedim> & mapping,
3648  const FiniteElement<dim, spacedim> &fe,
3649  const Quadrature<dim - 1> & quadrature);
3650 
3658  const Tensor<1, spacedim> &
3659  boundary_form(const unsigned int i) const;
3660 
3667  const std::vector<Tensor<1, spacedim>> &
3668  get_boundary_forms() const;
3669 
3674  unsigned int
3675  get_face_index() const;
3676 
3681  const Quadrature<dim - 1> &
3682  get_quadrature() const;
3683 
3688  std::size_t
3689  memory_consumption() const;
3690 
3691 protected:
3696  unsigned int present_face_index;
3697 
3701  const Quadrature<dim - 1> quadrature;
3702 };
3703 
3704 
3705 
3720 template <int dim, int spacedim = dim>
3721 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
3722 {
3723 public:
3728  static const unsigned int dimension = dim;
3729 
3730  static const unsigned int space_dimension = spacedim;
3731 
3736  static const unsigned int integral_dimension = dim - 1;
3737 
3742  FEFaceValues(const Mapping<dim, spacedim> & mapping,
3743  const FiniteElement<dim, spacedim> &fe,
3744  const Quadrature<dim - 1> & quadrature,
3745  const UpdateFlags update_flags);
3746 
3753  const Quadrature<dim - 1> & quadrature,
3754  const UpdateFlags update_flags);
3755 
3760  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3761  void
3762  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3763  level_dof_access>> &cell,
3764  const unsigned int face_no);
3765 
3772  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3773  void
3774  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3775  level_dof_access>> & cell,
3776  const typename Triangulation<dim, spacedim>::face_iterator &face);
3777 
3791  void
3792  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3793  const unsigned int face_no);
3794 
3795  /*
3796  * Reinitialize the gradients, Jacobi determinants, etc for the given face
3797  * on a given cell of type "iterator into a Triangulation object", and the
3798  * given finite element. Since iterators into a triangulation alone only
3799  * convey information about the geometry of a cell, but not about degrees of
3800  * freedom possibly associated with this cell, you will not be able to call
3801  * some functions of this class if they need information about degrees of
3802  * freedom. These functions are, above all, the
3803  * <tt>get_function_value/gradients/hessians/third_derivatives</tt>
3804  * functions. If you want to call these functions, you have to call the @p
3805  * reinit variants that take iterators into DoFHandler or other DoF handler
3806  * type objects.
3807  *
3808  * @note @p face must be one of @p cell's face iterators.
3809  */
3810  void
3811  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3812  const typename Triangulation<dim, spacedim>::face_iterator &face);
3813 
3829  get_present_fe_values() const;
3830 
3831 private:
3835  void
3836  initialize(const UpdateFlags update_flags);
3837 
3844  void
3845  do_reinit(const unsigned int face_no);
3846 };
3847 
3848 
3866 template <int dim, int spacedim = dim>
3867 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
3868 {
3869 public:
3873  static const unsigned int dimension = dim;
3874 
3878  static const unsigned int space_dimension = spacedim;
3879 
3884  static const unsigned int integral_dimension = dim - 1;
3885 
3890  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
3891  const FiniteElement<dim, spacedim> &fe,
3892  const Quadrature<dim - 1> & face_quadrature,
3893  const UpdateFlags update_flags);
3894 
3901  const Quadrature<dim - 1> & face_quadrature,
3902  const UpdateFlags update_flags);
3903 
3910  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3911  void
3912  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3913  level_dof_access>> &cell,
3914  const unsigned int face_no,
3915  const unsigned int subface_no);
3916 
3921  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3922  void
3923  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3924  level_dof_access>> & cell,
3925  const typename Triangulation<dim, spacedim>::face_iterator &face,
3926  const typename Triangulation<dim, spacedim>::face_iterator &subface);
3927 
3941  void
3942  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3943  const unsigned int face_no,
3944  const unsigned int subface_no);
3945 
3965  void
3966  reinit(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3967  const typename Triangulation<dim, spacedim>::face_iterator &face,
3968  const typename Triangulation<dim, spacedim>::face_iterator &subface);
3969 
3985  get_present_fe_values() const;
3986 
3992  DeclException0(ExcReinitCalledWithBoundaryFace);
3993 
3999  DeclException0(ExcFaceHasNoSubfaces);
4000 
4001 private:
4005  void
4006  initialize(const UpdateFlags update_flags);
4007 
4014  void
4015  do_reinit(const unsigned int face_no, const unsigned int subface_no);
4016 };
4017 
4018 
4019 #ifndef DOXYGEN
4020 
4021 
4022 /*------------------------ Inline functions: namespace FEValuesViews --------*/
4023 
4024 namespace FEValuesViews
4025 {
4026  template <int dim, int spacedim>
4027  inline typename Scalar<dim, spacedim>::value_type
4028  Scalar<dim, spacedim>::value(const unsigned int shape_function,
4029  const unsigned int q_point) const
4030  {
4031  Assert(shape_function < fe_values->fe->dofs_per_cell,
4032  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4033  Assert(
4034  fe_values->update_flags & update_values,
4036  "update_values"))));
4037 
4038  // an adaptation of the FEValuesBase::shape_value_component function
4039  // except that here we know the component as fixed and we have
4040  // pre-computed and cached a bunch of information. See the comments there.
4041  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4042  return fe_values->finite_element_output.shape_values(
4043  shape_function_data[shape_function].row_index, q_point);
4044  else
4045  return 0;
4046  }
4047 
4048 
4049 
4050  template <int dim, int spacedim>
4051  inline typename Scalar<dim, spacedim>::gradient_type
4052  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
4053  const unsigned int q_point) const
4054  {
4055  Assert(shape_function < fe_values->fe->dofs_per_cell,
4056  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4057  Assert(fe_values->update_flags & update_gradients,
4059  "update_gradients")));
4060 
4061  // an adaptation of the FEValuesBase::shape_grad_component
4062  // function except that here we know the component as fixed and we have
4063  // pre-computed and cached a bunch of information. See the comments there.
4064  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4065  return fe_values->finite_element_output
4066  .shape_gradients[shape_function_data[shape_function].row_index]
4067  [q_point];
4068  else
4069  return gradient_type();
4070  }
4071 
4072 
4073 
4074  template <int dim, int spacedim>
4075  inline typename Scalar<dim, spacedim>::hessian_type
4076  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
4077  const unsigned int q_point) const
4078  {
4079  Assert(shape_function < fe_values->fe->dofs_per_cell,
4080  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4081  Assert(fe_values->update_flags & update_hessians,
4083  "update_hessians")));
4084 
4085  // an adaptation of the FEValuesBase::shape_hessian_component
4086  // function except that here we know the component as fixed and we have
4087  // pre-computed and cached a bunch of information. See the comments there.
4088  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4089  return fe_values->finite_element_output
4090  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
4091  else
4092  return hessian_type();
4093  }
4094 
4095 
4096 
4097  template <int dim, int spacedim>
4098  inline typename Scalar<dim, spacedim>::third_derivative_type
4099  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
4100  const unsigned int q_point) const
4101  {
4102  Assert(shape_function < fe_values->fe->dofs_per_cell,
4103  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4104  Assert(fe_values->update_flags & update_3rd_derivatives,
4106  "update_3rd_derivatives")));
4107 
4108  // an adaptation of the FEValuesBase::shape_3rdderivative_component
4109  // function except that here we know the component as fixed and we have
4110  // pre-computed and cached a bunch of information. See the comments there.
4111  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4112  return fe_values->finite_element_output
4113  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
4114  [q_point];
4115  else
4116  return third_derivative_type();
4117  }
4118 
4119 
4120 
4121  template <int dim, int spacedim>
4122  inline typename Vector<dim, spacedim>::value_type
4123  Vector<dim, spacedim>::value(const unsigned int shape_function,
4124  const unsigned int q_point) const
4125  {
4126  Assert(shape_function < fe_values->fe->dofs_per_cell,
4127  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4128  Assert(fe_values->update_flags & update_values,
4130  "update_values")));
4131 
4132  // same as for the scalar case except that we have one more index
4133  const int snc =
4134  shape_function_data[shape_function].single_nonzero_component;
4135  if (snc == -2)
4136  return value_type();
4137  else if (snc != -1)
4138  {
4139  value_type return_value;
4140  return_value[shape_function_data[shape_function]
4141  .single_nonzero_component_index] =
4142  fe_values->finite_element_output.shape_values(snc, q_point);
4143  return return_value;
4144  }
4145  else
4146  {
4147  value_type return_value;
4148  for (unsigned int d = 0; d < dim; ++d)
4149  if (shape_function_data[shape_function]
4150  .is_nonzero_shape_function_component[d])
4151  return_value[d] = fe_values->finite_element_output.shape_values(
4152  shape_function_data[shape_function].row_index[d], q_point);
4153 
4154  return return_value;
4155  }
4156  }
4157 
4158 
4159 
4160  template <int dim, int spacedim>
4161  inline typename Vector<dim, spacedim>::gradient_type
4162  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4163  const unsigned int q_point) const
4164  {
4165  Assert(shape_function < fe_values->fe->dofs_per_cell,
4166  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4167  Assert(fe_values->update_flags & update_gradients,
4169  "update_gradients")));
4170 
4171  // same as for the scalar case except that we have one more index
4172  const int snc =
4173  shape_function_data[shape_function].single_nonzero_component;
4174  if (snc == -2)
4175  return gradient_type();
4176  else if (snc != -1)
4177  {
4178  gradient_type return_value;
4179  return_value[shape_function_data[shape_function]
4180  .single_nonzero_component_index] =
4181  fe_values->finite_element_output.shape_gradients[snc][q_point];
4182  return return_value;
4183  }
4184  else
4185  {
4186  gradient_type return_value;
4187  for (unsigned int d = 0; d < dim; ++d)
4188  if (shape_function_data[shape_function]
4189  .is_nonzero_shape_function_component[d])
4190  return_value[d] =
4191  fe_values->finite_element_output.shape_gradients
4192  [shape_function_data[shape_function].row_index[d]][q_point];
4193 
4194  return return_value;
4195  }
4196  }
4197 
4198 
4199 
4200  template <int dim, int spacedim>
4201  inline typename Vector<dim, spacedim>::divergence_type
4202  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4203  const unsigned int q_point) const
4204  {
4205  // this function works like in the case above
4206  Assert(shape_function < fe_values->fe->dofs_per_cell,
4207  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4208  Assert(fe_values->update_flags & update_gradients,
4210  "update_gradients")));
4211 
4212  // same as for the scalar case except that we have one more index
4213  const int snc =
4214  shape_function_data[shape_function].single_nonzero_component;
4215  if (snc == -2)
4216  return divergence_type();
4217  else if (snc != -1)
4218  return fe_values->finite_element_output
4219  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4220  .single_nonzero_component_index];
4221  else
4222  {
4223  divergence_type return_value = 0;
4224  for (unsigned int d = 0; d < dim; ++d)
4225  if (shape_function_data[shape_function]
4226  .is_nonzero_shape_function_component[d])
4227  return_value +=
4228  fe_values->finite_element_output.shape_gradients
4229  [shape_function_data[shape_function].row_index[d]][q_point][d];
4230 
4231  return return_value;
4232  }
4233  }
4234 
4235 
4236 
4237  template <int dim, int spacedim>
4238  inline typename Vector<dim, spacedim>::curl_type
4239  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4240  const unsigned int q_point) const
4241  {
4242  // this function works like in the case above
4243 
4244  Assert(shape_function < fe_values->fe->dofs_per_cell,
4245  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4246  Assert(fe_values->update_flags & update_gradients,
4248  "update_gradients")));
4249  // same as for the scalar case except that we have one more index
4250  const int snc =
4251  shape_function_data[shape_function].single_nonzero_component;
4252 
4253  if (snc == -2)
4254  return curl_type();
4255 
4256  else
4257  switch (dim)
4258  {
4259  case 1:
4260  {
4261  Assert(false,
4262  ExcMessage(
4263  "Computing the curl in 1d is not a useful operation"));
4264  return curl_type();
4265  }
4266 
4267  case 2:
4268  {
4269  if (snc != -1)
4270  {
4271  curl_type return_value;
4272 
4273  // the single nonzero component can only be zero or one in 2d
4274  if (shape_function_data[shape_function]
4275  .single_nonzero_component_index == 0)
4276  return_value[0] =
4277  -1.0 * fe_values->finite_element_output
4278  .shape_gradients[snc][q_point][1];
4279  else
4280  return_value[0] = fe_values->finite_element_output
4281  .shape_gradients[snc][q_point][0];
4282 
4283  return return_value;
4284  }
4285 
4286  else
4287  {
4288  curl_type return_value;
4289 
4290  return_value[0] = 0.0;
4291 
4292  if (shape_function_data[shape_function]
4293  .is_nonzero_shape_function_component[0])
4294  return_value[0] -=
4295  fe_values->finite_element_output
4296  .shape_gradients[shape_function_data[shape_function]
4297  .row_index[0]][q_point][1];
4298 
4299  if (shape_function_data[shape_function]
4300  .is_nonzero_shape_function_component[1])
4301  return_value[0] +=
4302  fe_values->finite_element_output
4303  .shape_gradients[shape_function_data[shape_function]
4304  .row_index[1]][q_point][0];
4305 
4306  return return_value;
4307  }
4308  }
4309 
4310  case 3:
4311  {
4312  if (snc != -1)
4313  {
4314  curl_type return_value;
4315 
4316  switch (shape_function_data[shape_function]
4317  .single_nonzero_component_index)
4318  {
4319  case 0:
4320  {
4321  return_value[0] = 0;
4322  return_value[1] = fe_values->finite_element_output
4323  .shape_gradients[snc][q_point][2];
4324  return_value[2] =
4325  -1.0 * fe_values->finite_element_output
4326  .shape_gradients[snc][q_point][1];
4327  return return_value;
4328  }
4329 
4330  case 1:
4331  {
4332  return_value[0] =
4333  -1.0 * fe_values->finite_element_output
4334  .shape_gradients[snc][q_point][2];
4335  return_value[1] = 0;
4336  return_value[2] = fe_values->finite_element_output
4337  .shape_gradients[snc][q_point][0];
4338  return return_value;
4339  }
4340 
4341  default:
4342  {
4343  return_value[0] = fe_values->finite_element_output
4344  .shape_gradients[snc][q_point][1];
4345  return_value[1] =
4346  -1.0 * fe_values->finite_element_output
4347  .shape_gradients[snc][q_point][0];
4348  return_value[2] = 0;
4349  return return_value;
4350  }
4351  }
4352  }
4353 
4354  else
4355  {
4356  curl_type return_value;
4357 
4358  for (unsigned int i = 0; i < dim; ++i)
4359  return_value[i] = 0.0;
4360 
4361  if (shape_function_data[shape_function]
4362  .is_nonzero_shape_function_component[0])
4363  {
4364  return_value[1] +=
4365  fe_values->finite_element_output
4366  .shape_gradients[shape_function_data[shape_function]
4367  .row_index[0]][q_point][2];
4368  return_value[2] -=
4369  fe_values->finite_element_output
4370  .shape_gradients[shape_function_data[shape_function]
4371  .row_index[0]][q_point][1];
4372  }
4373 
4374  if (shape_function_data[shape_function]
4375  .is_nonzero_shape_function_component[1])
4376  {
4377  return_value[0] -=
4378  fe_values->finite_element_output
4379  .shape_gradients[shape_function_data[shape_function]
4380  .row_index[1]][q_point][2];
4381  return_value[2] +=
4382  fe_values->finite_element_output
4383  .shape_gradients[shape_function_data[shape_function]
4384  .row_index[1]][q_point][0];
4385  }
4386 
4387  if (shape_function_data[shape_function]
4388  .is_nonzero_shape_function_component[2])
4389  {
4390  return_value[0] +=
4391  fe_values->finite_element_output
4392  .shape_gradients[shape_function_data[shape_function]
4393  .row_index[2]][q_point][1];
4394  return_value[1] -=
4395  fe_values->finite_element_output
4396  .shape_gradients[shape_function_data[shape_function]
4397  .row_index[2]][q_point][0];
4398  }
4399 
4400  return return_value;
4401  }
4402  }
4403  }
4404  // should not end up here
4405  Assert(false, ExcInternalError());
4406  return curl_type();
4407  }
4408 
4409 
4410 
4411  template <int dim, int spacedim>
4412  inline typename Vector<dim, spacedim>::hessian_type
4413  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
4414  const unsigned int q_point) const
4415  {
4416  // this function works like in the case above
4417  Assert(shape_function < fe_values->fe->dofs_per_cell,
4418  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4419  Assert(fe_values->update_flags & update_hessians,
4421  "update_hessians")));
4422 
4423  // same as for the scalar case except that we have one more index
4424  const int snc =
4425  shape_function_data[shape_function].single_nonzero_component;
4426  if (snc == -2)
4427  return hessian_type();
4428  else if (snc != -1)
4429  {
4430  hessian_type return_value;
4431  return_value[shape_function_data[shape_function]
4432  .single_nonzero_component_index] =
4433  fe_values->finite_element_output.shape_hessians[snc][q_point];
4434  return return_value;
4435  }
4436  else
4437  {
4438  hessian_type return_value;
4439  for (unsigned int d = 0; d < dim; ++d)
4440  if (shape_function_data[shape_function]
4441  .is_nonzero_shape_function_component[d])
4442  return_value[d] =
4443  fe_values->finite_element_output.shape_hessians
4444  [shape_function_data[shape_function].row_index[d]][q_point];
4445 
4446  return return_value;
4447  }
4448  }
4449 
4450 
4451 
4452  template <int dim, int spacedim>
4453  inline typename Vector<dim, spacedim>::third_derivative_type
4454  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
4455  const unsigned int q_point) const
4456  {
4457  // this function works like in the case above
4458  Assert(shape_function < fe_values->fe->dofs_per_cell,
4459  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4460  Assert(fe_values->update_flags & update_3rd_derivatives,
4462  "update_3rd_derivatives")));
4463 
4464  // same as for the scalar case except that we have one more index
4465  const int snc =
4466  shape_function_data[shape_function].single_nonzero_component;
4467  if (snc == -2)
4468  return third_derivative_type();
4469  else if (snc != -1)
4470  {
4471  third_derivative_type return_value;
4472  return_value[shape_function_data[shape_function]
4473  .single_nonzero_component_index] =
4474  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
4475  return return_value;
4476  }
4477  else
4478  {
4479  third_derivative_type return_value;
4480  for (unsigned int d = 0; d < dim; ++d)
4481  if (shape_function_data[shape_function]
4482  .is_nonzero_shape_function_component[d])
4483  return_value[d] =
4484  fe_values->finite_element_output.shape_3rd_derivatives
4485  [shape_function_data[shape_function].row_index[d]][q_point];
4486 
4487  return return_value;
4488  }
4489  }
4490 
4491 
4492 
4493  namespace internal
4494  {
4499  inline ::SymmetricTensor<2, 1>
4500  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
4501  {
4502  Assert(n < 1, ExcIndexRange(n, 0, 1));
4503  (void)n;
4504 
4505  return {{t[0]}};
4506  }
4507 
4508 
4509 
4510  inline ::SymmetricTensor<2, 2>
4511  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
4512  {
4513  switch (n)
4514  {
4515  case 0:
4516  {
4517  return {{t[0], 0, t[1] / 2}};
4518  }
4519  case 1:
4520  {
4521  return {{0, t[1], t[0] / 2}};
4522  }
4523  default:
4524  {
4525  Assert(false, ExcIndexRange(n, 0, 2));
4526  return {};
4527  }
4528  }
4529  }
4530 
4531 
4532 
4533  inline ::SymmetricTensor<2, 3>
4534  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
4535  {
4536  switch (n)
4537  {
4538  case 0:
4539  {
4540  return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
4541  }
4542  case 1:
4543  {
4544  return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
4545  }
4546  case 2:
4547  {
4548  return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
4549  }
4550  default:
4551  {
4552  Assert(false, ExcIndexRange(n, 0, 3));
4553  return {};
4554  }
4555  }
4556  }
4557  } // namespace internal
4558 
4559 
4560 
4561  template <int dim, int spacedim>
4562  inline typename Vector<dim, spacedim>::symmetric_gradient_type
4563  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
4564  const unsigned int q_point) const
4565  {
4566  Assert(shape_function < fe_values->fe->dofs_per_cell,
4567  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4568  Assert(fe_values->update_flags & update_gradients,
4570  "update_gradients")));
4571 
4572  // same as for the scalar case except that we have one more index
4573  const int snc =
4574  shape_function_data[shape_function].single_nonzero_component;
4575  if (snc == -2)
4576  return symmetric_gradient_type();
4577  else if (snc != -1)
4578  return internal::symmetrize_single_row(
4579  shape_function_data[shape_function].single_nonzero_component_index,
4580  fe_values->finite_element_output.shape_gradients[snc][q_point]);
4581  else
4582  {
4583  gradient_type return_value;
4584  for (unsigned int d = 0; d < dim; ++d)
4585  if (shape_function_data[shape_function]
4586  .is_nonzero_shape_function_component[d])
4587  return_value[d] =
4588  fe_values->finite_element_output.shape_gradients
4589  [shape_function_data[shape_function].row_index[d]][q_point];
4590 
4591  return symmetrize(return_value);
4592  }
4593  }
4594 
4595 
4596 
4597  template <int dim, int spacedim>
4599  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
4600  const unsigned int q_point) const
4601  {
4602  Assert(shape_function < fe_values->fe->dofs_per_cell,
4603  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4604  Assert(fe_values->update_flags & update_values,
4606  "update_values")));
4607 
4608  // similar to the vector case where we have more then one index and we need
4609  // to convert between unrolled and component indexing for tensors
4610  const int snc =
4611  shape_function_data[shape_function].single_nonzero_component;
4612 
4613  if (snc == -2)
4614  {
4615  // shape function is zero for the selected components
4616  return value_type();
4617  }
4618  else if (snc != -1)
4619  {
4620  value_type return_value;
4621  const unsigned int comp =
4622  shape_function_data[shape_function].single_nonzero_component_index;
4623  return_value[value_type::unrolled_to_component_indices(comp)] =
4624  fe_values->finite_element_output.shape_values(snc, q_point);
4625  return return_value;
4626  }
4627  else
4628  {
4629  value_type return_value;
4630  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
4631  if (shape_function_data[shape_function]
4632  .is_nonzero_shape_function_component[d])
4633  return_value[value_type::unrolled_to_component_indices(d)] =
4634  fe_values->finite_element_output.shape_values(
4635  shape_function_data[shape_function].row_index[d], q_point);
4636  return return_value;
4637  }
4638  }
4639 
4640 
4641 
4642  template <int dim, int spacedim>
4645  const unsigned int shape_function,
4646  const unsigned int q_point) const
4647  {
4648  Assert(shape_function < fe_values->fe->dofs_per_cell,
4649  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4650  Assert(fe_values->update_flags & update_gradients,
4652  "update_gradients")));
4653 
4654  const int snc =
4655  shape_function_data[shape_function].single_nonzero_component;
4656 
4657  if (snc == -2)
4658  {
4659  // shape function is zero for the selected components
4660  return divergence_type();
4661  }
4662  else if (snc != -1)
4663  {
4664  // we have a single non-zero component when the symmetric tensor is
4665  // represented in unrolled form. this implies we potentially have
4666  // two non-zero components when represented in component form! we
4667  // will only have one non-zero entry if the non-zero component lies on
4668  // the diagonal of the tensor.
4669  //
4670  // the divergence of a second-order tensor is a first order tensor.
4671  //
4672  // assume the second-order tensor is A with components A_{ij}. then
4673  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
4674  // entries in the tensorial representation. define the
4675  // divergence as:
4676  // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
4677  // (which is incidentally also
4678  // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
4679  // In both cases, a sum is implied.
4680  //
4681  // Now, we know the nonzero component in unrolled form: it is indicated
4682  // by 'snc'. we can figure out which tensor components belong to this:
4683  const unsigned int comp =
4684  shape_function_data[shape_function].single_nonzero_component_index;
4685  const unsigned int ii =
4686  value_type::unrolled_to_component_indices(comp)[0];
4687  const unsigned int jj =
4688  value_type::unrolled_to_component_indices(comp)[1];
4689 
4690  // given the form of the divergence above, if ii=jj there is only a
4691  // single nonzero component of the full tensor and the gradient
4692  // equals
4693  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
4694  // all other entries of 'b' are zero
4695  //
4696  // on the other hand, if ii!=jj, then there are two nonzero entries in
4697  // the full tensor and
4698  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
4699  // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
4700  // again, all other entries of 'b' are zero
4701  const ::Tensor<1, spacedim> &phi_grad =
4702  fe_values->finite_element_output.shape_gradients[snc][q_point];
4703 
4704  divergence_type return_value;
4705  return_value[ii] = phi_grad[jj];
4706 
4707  if (ii != jj)
4708  return_value[jj] = phi_grad[ii];
4709 
4710  return return_value;
4711  }
4712  else
4713  {
4714  Assert(false, ExcNotImplemented());
4715  divergence_type return_value;
4716  return return_value;
4717  }
4718  }
4719 
4720 
4721 
4722  template <int dim, int spacedim>
4723  inline typename Tensor<2, dim, spacedim>::value_type
4724  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
4725  const unsigned int q_point) const
4726  {
4727  Assert(shape_function < fe_values->fe->dofs_per_cell,
4728  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4729  Assert(fe_values->update_flags & update_values,
4731  "update_values")));
4732 
4733  // similar to the vector case where we have more then one index and we need
4734  // to convert between unrolled and component indexing for tensors
4735  const int snc =
4736  shape_function_data[shape_function].single_nonzero_component;
4737 
4738  if (snc == -2)
4739  {
4740  // shape function is zero for the selected components
4741  return value_type();
4742  }
4743  else if (snc != -1)
4744  {
4745  value_type return_value;
4746  const unsigned int comp =
4747  shape_function_data[shape_function].single_nonzero_component_index;
4748  const TableIndices<2> indices =
4750  return_value[indices] =
4751  fe_values->finite_element_output.shape_values(snc, q_point);
4752  return return_value;
4753  }
4754  else
4755  {
4756  value_type return_value;
4757  for (unsigned int d = 0; d < dim * dim; ++d)
4758  if (shape_function_data[shape_function]
4759  .is_nonzero_shape_function_component[d])
4760  {
4761  const TableIndices<2> indices =
4763  return_value[indices] =
4764  fe_values->finite_element_output.shape_values(
4765  shape_function_data[shape_function].row_index[d], q_point);
4766  }
4767  return return_value;
4768  }
4769  }
4770 
4771 
4772 
4773  template <int dim, int spacedim>
4775  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
4776  const unsigned int q_point) const
4777  {
4778  Assert(shape_function < fe_values->fe->dofs_per_cell,
4779  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4780  Assert(fe_values->update_flags & update_gradients,
4782  "update_gradients")));
4783 
4784  const int snc =
4785  shape_function_data[shape_function].single_nonzero_component;
4786 
4787  if (snc == -2)
4788  {
4789  // shape function is zero for the selected components
4790  return divergence_type();
4791  }
4792  else if (snc != -1)
4793  {
4794  // we have a single non-zero component when the tensor is
4795  // represented in unrolled form.
4796  //
4797  // the divergence of a second-order tensor is a first order tensor.
4798  //
4799  // assume the second-order tensor is A with components A_{ij},
4800  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
4801  //
4802  // Now, we know the nonzero component in unrolled form: it is indicated
4803  // by 'snc'. we can figure out which tensor components belong to this:
4804  const unsigned int comp =
4805  shape_function_data[shape_function].single_nonzero_component_index;
4806  const TableIndices<2> indices =
4808  const unsigned int ii = indices[0];
4809  const unsigned int jj = indices[1];
4810 
4811  const ::Tensor<1, spacedim> &phi_grad =
4812  fe_values->finite_element_output.shape_gradients[snc][q_point];
4813 
4814  divergence_type return_value;
4815  // note that we contract \nabla from the right
4816  return_value[ii] = phi_grad[jj];
4817 
4818  return return_value;
4819  }
4820  else
4821  {
4822  Assert(false, ExcNotImplemented());
4823  divergence_type return_value;
4824  return return_value;
4825  }
4826  }
4827 
4828 
4829 
4830  template <int dim, int spacedim>
4832  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
4833  const unsigned int q_point) const
4834  {
4835  Assert(shape_function < fe_values->fe->dofs_per_cell,
4836  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4837  Assert(fe_values->update_flags & update_gradients,
4839  "update_gradients")));
4840 
4841  const int snc =
4842  shape_function_data[shape_function].single_nonzero_component;
4843 
4844  if (snc == -2)
4845  {
4846  // shape function is zero for the selected components
4847  return gradient_type();
4848  }
4849  else if (snc != -1)
4850  {
4851  // we have a single non-zero component when the tensor is
4852  // represented in unrolled form.
4853  //
4854  // the gradient of a second-order tensor is a third order tensor.
4855  //
4856  // assume the second-order tensor is A with components A_{ij},
4857  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
4858  //
4859  // Now, we know the nonzero component in unrolled form: it is indicated
4860  // by 'snc'. we can figure out which tensor components belong to this:
4861  const unsigned int comp =
4862  shape_function_data[shape_function].single_nonzero_component_index;
4863  const TableIndices<2> indices =
4865  const unsigned int ii = indices[0];
4866  const unsigned int jj = indices[1];
4867 
4868  const ::Tensor<1, spacedim> &phi_grad =
4869  fe_values->finite_element_output.shape_gradients[snc][q_point];
4870 
4871  gradient_type return_value;
4872  return_value[ii][jj] = phi_grad;
4873 
4874  return return_value;
4875  }
4876  else
4877  {
4878  Assert(false, ExcNotImplemented());
4879  gradient_type return_value;
4880  return return_value;
4881  }
4882  }
4883 
4884 } // namespace FEValuesViews
4885 
4886 
4887 
4888 /*---------------------- Inline functions: FEValuesBase ---------------------*/
4889 
4890 
4891 
4892 template <int dim, int spacedim>
4894  operator[](const FEValuesExtractors::Scalar &scalar) const
4895 {
4896  Assert(scalar.component < fe_values_views_cache.scalars.size(),
4897  ExcIndexRange(scalar.component,
4898  0,
4899  fe_values_views_cache.scalars.size()));
4900 
4901  return fe_values_views_cache.scalars[scalar.component];
4902 }
4903 
4904 
4905 
4906 template <int dim, int spacedim>
4908  operator[](const FEValuesExtractors::Vector &vector) const
4909 {
4910  Assert(vector.first_vector_component < fe_values_views_cache.vectors.size(),
4912  0,
4913  fe_values_views_cache.vectors.size()));
4914 
4915  return fe_values_views_cache.vectors[vector.first_vector_component];
4916 }
4917 
4918 
4919 
4920 template <int dim, int spacedim>
4924 {
4925  Assert(
4926  tensor.first_tensor_component <
4927  fe_values_views_cache.symmetric_second_order_tensors.size(),
4929  0,
4930  fe_values_views_cache.symmetric_second_order_tensors.size()));
4931 
4932  return fe_values_views_cache
4933  .symmetric_second_order_tensors[tensor.first_tensor_component];
4934 }
4935 
4936 
4937 
4938 template <int dim, int spacedim>
4941  operator[](const FEValuesExtractors::Tensor<2> &tensor) const
4942 {
4944  fe_values_views_cache.second_order_tensors.size(),
4946  0,
4947  fe_values_views_cache.second_order_tensors.size()));
4948 
4949  return fe_values_views_cache
4950  .second_order_tensors[tensor.first_tensor_component];
4951 }
4952 
4953 
4954 
4955 template <int dim, int spacedim>
4956 inline const double &
4957 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
4958  const unsigned int j) const
4959 {
4960  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4961  Assert(this->update_flags & update_values,
4962  ExcAccessToUninitializedField("update_values"));
4963  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4964  Assert(present_cell.get() != nullptr,
4965  ExcMessage("FEValues object is not reinit'ed to any cell"));
4966  // if the entire FE is primitive,
4967  // then we can take a short-cut:
4968  if (fe->is_primitive())
4969  return this->finite_element_output.shape_values(i, j);
4970  else
4971  {
4972  // otherwise, use the mapping
4973  // between shape function
4974  // numbers and rows. note that
4975  // by the assertions above, we
4976  // know that this particular
4977  // shape function is primitive,
4978  // so we can call
4979  // system_to_component_index
4980  const unsigned int row =
4981  this->finite_element_output
4982  .shape_function_to_row_table[i * fe->n_components() +
4983  fe->system_to_component_index(i).first];
4984  return this->finite_element_output.shape_values(row, j);
4985  }
4986 }
4987 
4988 
4989 
4990 template <int dim, int spacedim>
4991 inline double
4993  const unsigned int i,
4994  const unsigned int j,
4995  const unsigned int component) const
4996 {
4997  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4998  Assert(this->update_flags & update_values,
4999  ExcAccessToUninitializedField("update_values"));
5000  Assert(component < fe->n_components(),
5001  ExcIndexRange(component, 0, fe->n_components()));
5002  Assert(present_cell.get() != nullptr,
5003  ExcMessage("FEValues object is not reinit'ed to any cell"));
5004 
5005  // check whether the shape function
5006  // is non-zero at all within
5007  // this component:
5008  if (fe->get_nonzero_components(i)[component] == false)
5009  return 0;
5010 
5011  // look up the right row in the
5012  // table and take the data from
5013  // there
5014  const unsigned int row =
5015  this->finite_element_output
5016  .shape_function_to_row_table[i * fe->n_components() + component];
5017  return this->finite_element_output.shape_values(row, j);
5018 }
5019 
5020 
5021 
5022 template <int dim, int spacedim>
5023 inline const Tensor<1, spacedim> &
5024 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
5025  const unsigned int j) const
5026 {
5027  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5028  Assert(this->update_flags & update_gradients,
5029  ExcAccessToUninitializedField("update_gradients"));
5030  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5031  Assert(present_cell.get() != nullptr,
5032  ExcMessage("FEValues object is not reinit'ed to any cell"));
5033  // if the entire FE is primitive,
5034  // then we can take a short-cut:
5035  if (fe->is_primitive())
5036  return this->finite_element_output.shape_gradients[i][j];
5037  else
5038  {
5039  // otherwise, use the mapping
5040  // between shape function
5041  // numbers and rows. note that
5042  // by the assertions above, we
5043  // know that this particular
5044  // shape function is primitive,
5045  // so we can call
5046  // system_to_component_index
5047  const unsigned int row =
5048  this->finite_element_output
5049  .shape_function_to_row_table[i * fe->n_components() +
5050  fe->system_to_component_index(i).first];
5051  return this->finite_element_output.shape_gradients[row][j];
5052  }
5053 }
5054 
5055 
5056 
5057 template <int dim, int spacedim>
5058 inline Tensor<1, spacedim>
5060  const unsigned int i,
5061  const unsigned int j,
5062  const unsigned int component) const
5063 {
5064  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5065  Assert(this->update_flags & update_gradients,
5066  ExcAccessToUninitializedField("update_gradients"));
5067  Assert(component < fe->n_components(),
5068  ExcIndexRange(component, 0, fe->n_components()));
5069  Assert(present_cell.get() != nullptr,
5070  ExcMessage("FEValues object is not reinit'ed to any cell"));
5071  // check whether the shape function
5072  // is non-zero at all within
5073  // this component:
5074  if (fe->get_nonzero_components(i)[component] == false)
5075  return Tensor<1, spacedim>();
5076 
5077  // look up the right row in the
5078  // table and take the data from
5079  // there
5080  const unsigned int row =
5081  this->finite_element_output
5082  .shape_function_to_row_table[i * fe->n_components() + component];
5083  return this->finite_element_output.shape_gradients[row][j];
5084 }
5085 
5086 
5087 
5088 template <int dim, int spacedim>
5089 inline const Tensor<2, spacedim> &
5090 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
5091  const unsigned int j) const
5092 {
5093  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5094  Assert(this->update_flags & update_hessians,
5095  ExcAccessToUninitializedField("update_hessians"));
5096  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5097  Assert(present_cell.get() != nullptr,
5098  ExcMessage("FEValues object is not reinit'ed to any cell"));
5099  // if the entire FE is primitive,
5100  // then we can take a short-cut:
5101  if (fe->is_primitive())
5102  return this->finite_element_output.shape_hessians[i][j];
5103  else
5104  {
5105  // otherwise, use the mapping
5106  // between shape function
5107  // numbers and rows. note that
5108  // by the assertions above, we
5109  // know that this particular
5110  // shape function is primitive,
5111  // so we can call
5112  // system_to_component_index
5113  const unsigned int row =
5114  this->finite_element_output
5115  .shape_function_to_row_table[i * fe->n_components() +
5116  fe->system_to_component_index(i).first];
5117  return this->finite_element_output.shape_hessians[row][j];
5118  }
5119 }
5120 
5121 
5122 
5123 template <int dim, int spacedim>
5124 inline Tensor<2, spacedim>
5126  const unsigned int i,
5127  const unsigned int j,
5128  const unsigned int component) const
5129 {
5130  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5131  Assert(this->update_flags & update_hessians,
5132  ExcAccessToUninitializedField("update_hessians"));
5133  Assert(component < fe->n_components(),
5134  ExcIndexRange(component, 0, fe->n_components()));
5135  Assert(present_cell.get() != nullptr,
5136  ExcMessage("FEValues object is not reinit'ed to any cell"));
5137  // check whether the shape function
5138  // is non-zero at all within
5139  // this component:
5140  if (fe->get_nonzero_components(i)[component] == false)
5141  return Tensor<2, spacedim>();
5142 
5143  // look up the right row in the
5144  // table and take the data from
5145  // there
5146  const unsigned int row =
5147  this->finite_element_output
5148  .shape_function_to_row_table[i * fe->n_components() + component];
5149  return this->finite_element_output.shape_hessians[row][j];
5150 }
5151 
5152 
5153 
5154 template <int dim, int spacedim>
5155 inline const Tensor<3, spacedim> &
5157  const unsigned int j) const
5158 {
5159  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5160  Assert(this->update_flags & update_hessians,
5161  ExcAccessToUninitializedField("update_3rd_derivatives"));
5162  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5163  Assert(present_cell.get() != nullptr,
5164  ExcMessage("FEValues object is not reinit'ed to any cell"));
5165  // if the entire FE is primitive,
5166  // then we can take a short-cut:
5167  if (fe->is_primitive())
5168  return this->finite_element_output.shape_3rd_derivatives[i][j];
5169  else
5170  {
5171  // otherwise, use the mapping
5172  // between shape function
5173  // numbers and rows. note that
5174  // by the assertions above, we
5175  // know that this particular
5176  // shape function is primitive,
5177  // so we can call
5178  // system_to_component_index
5179  const unsigned int row =
5180  this->finite_element_output
5181  .shape_function_to_row_table[i * fe->n_components() +
5182  fe->system_to_component_index(i).first];
5183  return this->finite_element_output.shape_3rd_derivatives[row][j];
5184  }
5185 }
5186 
5187 
5188 
5189 template <int dim, int spacedim>
5190 inline Tensor<3, spacedim>
5192  const unsigned int i,
5193  const unsigned int j,
5194  const unsigned int component) const
5195 {
5196  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5197  Assert(this->update_flags & update_hessians,
5198  ExcAccessToUninitializedField("update_3rd_derivatives"));
5199  Assert(component < fe->n_components(),
5200  ExcIndexRange(component, 0, fe->n_components()));
5201  Assert(present_cell.get() != nullptr,
5202  ExcMessage("FEValues object is not reinit'ed to any cell"));
5203  // check whether the shape function
5204  // is non-zero at all within
5205  // this component:
5206  if (fe->get_nonzero_components(i)[component] == false)
5207  return Tensor<3, spacedim>();
5208 
5209  // look up the right row in the
5210  // table and take the data from
5211  // there
5212  const unsigned int row =
5213  this->finite_element_output
5214  .shape_function_to_row_table[i * fe->n_components() + component];
5215  return this->finite_element_output.shape_3rd_derivatives[row][j];
5216 }
5217 
5218 
5219 
5220 template <int dim, int spacedim>
5221 inline const FiniteElement<dim, spacedim> &
5223 {
5224  return *fe;
5225 }
5226 
5227 
5228 
5229 template <int dim, int spacedim>
5230 inline const Mapping<dim, spacedim> &
5232 {
5233  return *mapping;
5234 }
5235 
5236 
5237 
5238 template <int dim, int spacedim>
5239 inline UpdateFlags
5241 {
5242  return this->update_flags;
5243 }
5244 
5245 
5246 
5247 template <int dim, int spacedim>
5248 inline const std::vector<Point<spacedim>> &
5250 {
5251  Assert(this->update_flags & update_quadrature_points,
5252  ExcAccessToUninitializedField("update_quadrature_points"));
5253  Assert(present_cell.get() != nullptr,
5254  ExcMessage("FEValues object is not reinit'ed to any cell"));
5255  return this->mapping_output.quadrature_points;
5256 }
5257 
5258 
5259 
5260 template <int dim, int spacedim>
5261 inline const std::vector<double> &
5263 {
5264  Assert(this->update_flags & update_JxW_values,
5265  ExcAccessToUninitializedField("update_JxW_values"));
5266  Assert(present_cell.get() != nullptr,
5267  ExcMessage("FEValues object is not reinit'ed to any cell"));
5268  return this->mapping_output.JxW_values;
5269 }
5270 
5271 
5272 
5273 template <int dim, int spacedim>
5274 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5276 {
5277  Assert(this->update_flags & update_jacobians,
5278  ExcAccessToUninitializedField("update_jacobians"));
5279  Assert(present_cell.get() != nullptr,
5280  ExcMessage("FEValues object is not reinit'ed to any cell"));
5281  return this->mapping_output.jacobians;
5282 }
5283 
5284 
5285 
5286 template <int dim, int spacedim>
5287 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5289 {
5290  Assert(this->update_flags & update_jacobian_grads,
5291  ExcAccessToUninitializedField("update_jacobians_grads"));
5292  Assert(present_cell.get() != nullptr,
5293  ExcMessage("FEValues object is not reinit'ed to any cell"));
5294  return this->mapping_output.jacobian_grads;
5295 }
5296 
5297 
5298 
5299 template <int dim, int spacedim>
5300 inline const Tensor<3, spacedim> &
5302  const unsigned int i) const
5303 {
5304  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5305  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5306  Assert(present_cell.get() != nullptr,
5307  ExcMessage("FEValues object is not reinit'ed to any cell"));
5308  return this->mapping_output.jacobian_pushed_forward_grads[i];
5309 }
5310 
5311 
5312 
5313 template <int dim, int spacedim>
5314 inline const std::vector<Tensor<3, spacedim>> &
5316 {
5317  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5318  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5319  Assert(present_cell.get() != nullptr,
5320  ExcMessage("FEValues object is not reinit'ed to any cell"));
5321  return this->mapping_output.jacobian_pushed_forward_grads;
5322 }
5323 
5324 
5325 
5326 template <int dim, int spacedim>
5327 inline const DerivativeForm<3, dim, spacedim> &
5328 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5329 {
5330  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5331  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5332  Assert(present_cell.get() != nullptr,
5333  ExcMessage("FEValues object is not reinit'ed to any cell"));
5334  return this->mapping_output.jacobian_2nd_derivatives[i];
5335 }
5336 
5337 
5338 
5339 template <int dim, int spacedim>
5340 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5342 {
5343  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5344  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5345  Assert(present_cell.get() != nullptr,
5346  ExcMessage("FEValues object is not reinit'ed to any cell"));
5347  return this->mapping_output.jacobian_2nd_derivatives;
5348 }
5349 
5350 
5351 
5352 template <int dim, int spacedim>
5353 inline const Tensor<4, spacedim> &
5355  const unsigned int i) const
5356 {
5359  "update_jacobian_pushed_forward_2nd_derivatives"));
5360  Assert(present_cell.get() != nullptr,
5361  ExcMessage("FEValues object is not reinit'ed to any cell"));
5362  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5363 }
5364 
5365 
5366 
5367 template <int dim, int spacedim>
5368 inline const std::vector<Tensor<4, spacedim>> &
5370 {
5371  Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
5373  "update_jacobian_pushed_forward_2nd_derivatives"));
5374  Assert(present_cell.get() != nullptr,
5375  ExcMessage("FEValues object is not reinit'ed to any cell"));
5376  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5377 }
5378 
5379 
5380 
5381 template <int dim, int spacedim>
5382 inline const DerivativeForm<4, dim, spacedim> &
5383 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5384 {
5385  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5386  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5387  Assert(present_cell.get() != nullptr,
5388  ExcMessage("FEValues object is not reinit'ed to any cell"));
5389  return this->mapping_output.jacobian_3rd_derivatives[i];
5390 }
5391 
5392 
5393 
5394 template <int dim, int spacedim>
5395 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5397 {
5398  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5399  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5400  Assert(present_cell.get() != nullptr,
5401  ExcMessage("FEValues object is not reinit'ed to any cell"));
5402  return this->mapping_output.jacobian_3rd_derivatives;
5403 }
5404 
5405 
5406 
5407 template <int dim, int spacedim>
5408 inline const Tensor<5, spacedim> &
5410  const unsigned int i) const
5411 {
5414  "update_jacobian_pushed_forward_3rd_derivatives"));
5415  Assert(present_cell.get() != nullptr,
5416  ExcMessage("FEValues object is not reinit'ed to any cell"));
5417  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
5418 }
5419 
5420 
5421 
5422 template <int dim, int spacedim>
5423 inline const std::vector<Tensor<5, spacedim>> &
5425 {
5426  Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
5428  "update_jacobian_pushed_forward_3rd_derivatives"));
5429  Assert(present_cell.get() != nullptr,
5430  ExcMessage("FEValues object is not reinit'ed to any cell"));
5431  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
5432 }
5433 
5434 
5435 
5436 template <int dim, int spacedim>
5437 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
5439 {
5440  Assert(this->update_flags & update_inverse_jacobians,
5441  ExcAccessToUninitializedField("update_inverse_jacobians"));
5442  Assert(present_cell.get() != nullptr,
5443  ExcMessage("FEValues object is not reinit'ed to any cell"));
5444  return this->mapping_output.inverse_jacobians;
5445 }
5446 
5447 
5448 
5449 template <int dim, int spacedim>
5450 inline const Point<spacedim> &
5451 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
5452 {
5453  Assert(this->update_flags & update_quadrature_points,
5454  ExcAccessToUninitializedField("update_quadrature_points"));
5455  Assert(i < this->mapping_output.quadrature_points.size(),
5456  ExcIndexRange(i, 0, this->mapping_output.quadrature_points.size()));
5457  Assert(present_cell.get() != nullptr,
5458  ExcMessage("FEValues object is not reinit'ed to any cell"));
5459 
5460  return this->mapping_output.quadrature_points[i];
5461 }
5462 
5463 
5464 
5465 template <int dim, int spacedim>
5466 inline double
5467 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
5468 {
5469  Assert(this->update_flags & update_JxW_values,
5470  ExcAccessToUninitializedField("update_JxW_values"));
5471  Assert(i < this->mapping_output.JxW_values.size(),
5472  ExcIndexRange(i, 0, this->mapping_output.JxW_values.size()));
5473  Assert(present_cell.get() != nullptr,
5474  ExcMessage("FEValues object is not reinit'ed to any cell"));
5475 
5476  return this->mapping_output.JxW_values[i];
5477 }
5478 
5479 
5480 
5481 template <int dim, int spacedim>
5482 inline const DerivativeForm<1, dim, spacedim> &
5483 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
5484 {
5485  Assert(this->update_flags & update_jacobians,
5486  ExcAccessToUninitializedField("update_jacobians"));
5487  Assert(i < this->mapping_output.jacobians.size(),
5488  ExcIndexRange(i, 0, this->mapping_output.jacobians.size()));
5489  Assert(present_cell.get() != nullptr,
5490  ExcMessage("FEValues object is not reinit'ed to any cell"));
5491 
5492  return this->mapping_output.jacobians[i];
5493 }
5494 
5495 
5496 
5497 template <int dim, int spacedim>
5498 inline const DerivativeForm<2, dim, spacedim> &
5499 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
5500 {
5501  Assert(this->update_flags & update_jacobian_grads,
5502  ExcAccessToUninitializedField("update_jacobians_grads"));
5503  Assert(i < this->mapping_output.jacobian_grads.size(),
5504  ExcIndexRange(i, 0, this->mapping_output.jacobian_grads.size()));
5505  Assert(present_cell.get() != nullptr,
5506  ExcMessage("FEValues object is not reinit'ed to any cell"));
5507 
5508  return this->mapping_output.jacobian_grads[i];
5509 }
5510 
5511 
5512 
5513 template <int dim, int spacedim>
5514 inline const DerivativeForm<1, spacedim, dim> &
5515 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
5516 {
5517  Assert(this->update_flags & update_inverse_jacobians,
5518  ExcAccessToUninitializedField("update_inverse_jacobians"));
5519  Assert(i < this->mapping_output.inverse_jacobians.size(),
5520  ExcIndexRange(i, 0, this->mapping_output.inverse_jacobians.size()));
5521  Assert(present_cell.get() != nullptr,
5522  ExcMessage("FEValues object is not reinit'ed to any cell"));
5523 
5524  return this->mapping_output.inverse_jacobians[i];
5525 }
5526 
5527 
5528 
5529 template <int dim, int spacedim>
5530 inline const Tensor<1, spacedim> &
5531 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
5532 {
5533  Assert(this->update_flags & update_normal_vectors,
5535  "update_normal_vectors")));
5536  Assert(i < this->mapping_output.normal_vectors.size(),
5537  ExcIndexRange(i, 0, this->mapping_output.normal_vectors.size()));
5538  Assert(present_cell.get() != nullptr,
5539  ExcMessage("FEValues object is not reinit'ed to any cell"));
5540 
5541  return this->mapping_output.normal_vectors[i];
5542 }
5543 
5544 
5545 
5546 /*--------------------- Inline functions: FEValues --------------------------*/
5547 
5548 
5549 template <int dim, int spacedim>
5550 inline const Quadrature<dim> &
5552 {
5553  return quadrature;
5554 }
5555 
5556 
5557 
5558 template <int dim, int spacedim>
5559 inline const FEValues<dim, spacedim> &
5561 {
5562  return *this;
5563 }
5564 
5565 
5566 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
5567 
5568 
5569 template <int dim, int spacedim>
5570 inline unsigned int
5572 {
5573  return present_face_index;
5574 }
5575 
5576 
5577 /*----------------------- Inline functions: FE*FaceValues -------------------*/
5578 
5579 template <int dim, int spacedim>
5580 inline const Quadrature<dim - 1> &
5582 {
5583  return quadrature;
5584 }
5585 
5586 
5587 
5588 template <int dim, int spacedim>
5589 inline const FEFaceValues<dim, spacedim> &
5591 {
5592  return *this;
5593 }
5594 
5595 
5596 
5597 template <int dim, int spacedim>
5598 inline const FESubfaceValues<dim, spacedim> &
5600 {
5601  return *this;
5602 }
5603 
5604 
5605 
5606 template <int dim, int spacedim>
5607 inline const Tensor<1, spacedim> &
5608 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
5609 {
5610  Assert(i < this->mapping_output.boundary_forms.size(),
5611  ExcIndexRange(i, 0, this->mapping_output.boundary_forms.size()));
5612  Assert(this->update_flags & update_boundary_forms,
5614  "update_boundary_forms")));
5615 
5616  return this->mapping_output.boundary_forms[i];
5617 }
5618 
5619 #endif // DOXYGEN
5620 
5621 DEAL_II_NAMESPACE_CLOSE
5622 
5623 #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:3431
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:3463
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:3696
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:3701
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:3347
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:3399
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:3478
LinearAlgebra::distributed::Vector< Number > Vector
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
static ::ExceptionBase & ExcFENotPrimitive()
Tensor< 2, spacedim > shape_hessian_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
typename ProductType< Number, typename Vector< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:664
const Quadrature< dim - 1 > & get_quadrature() const
const std::vector< Point< spacedim > > & get_quadrature_points() const
const Point< spacedim > & quadrature_point(const unsigned int q) const
const std::vector< Tensor< 4, spacedim > > & get_jacobian_pushed_forward_2nd_derivatives() const
typename ProductType< Number, typename Scalar< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:197
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1581
static ::ExceptionBase & ExcMessage(std::string arg1)
constexpr SymmetricTensor()=default
Number value_type
Definition: vector.h:125
typename ProductType< Number, typename Vector< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:712
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:518
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
Tensor< 3, spacedim > shape_3rd_derivative_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
Third derivatives of shape functions.
const FEValues< dim, spacedim > & get_present_fe_values() const
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1867
unsigned int get_face_index() const
#define Assert(cond, exc)
Definition: exceptions.h:1411
UpdateFlags
const Tensor< 2, spacedim > & shape_hessian(const unsigned int function_no, const unsigned int point_no) const
Abstract base class for mapping classes.
Definition: mapping.h:302
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1227
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
const Quadrature< dim > quadrature
Definition: fe_values.h:3594
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:3407
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:3372
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:1525
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:3363
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3439
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:3414
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:3445
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:3423
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