Reference documentation for deal.II version Git 6a05540 2018-06-21 10:59:02 -0400
fe_values.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2018 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 #include <deal.II/base/table.h>
29 #include <deal.II/base/vector_slice.h>
30 
31 #include <deal.II/dofs/dof_accessor.h>
32 #include <deal.II/dofs/dof_handler.h>
33 
34 #include <deal.II/fe/fe.h>
35 #include <deal.II/fe/fe_update_flags.h>
36 #include <deal.II/fe/fe_values_extractors.h>
37 #include <deal.II/fe/mapping.h>
38 
39 #include <deal.II/grid/tria.h>
40 #include <deal.II/grid/tria_iterator.h>
41 
42 #include <deal.II/hp/dof_handler.h>
43 
44 #include <algorithm>
45 #include <memory>
46 #include <type_traits>
47 
48 
49 // dummy include in order to have the
50 // definition of PetscScalar available
51 // without including other PETSc stuff
52 #ifdef DEAL_II_WITH_PETSC
53 # include <petsc.h>
54 #endif
55 
56 DEAL_II_NAMESPACE_OPEN
57 
58 template <int dim, int spacedim = dim>
59 class FEValuesBase;
60 
61 namespace internal
62 {
67  template <int dim, class NumberType = double>
68  struct CurlType;
69 
76  template <class NumberType>
77  struct CurlType<1, NumberType>
78  {
80  };
81 
88  template <class NumberType>
89  struct CurlType<2, NumberType>
90  {
92  };
93 
100  template <class NumberType>
101  struct CurlType<3, NumberType>
102  {
104  };
105 } // namespace internal
106 
107 
108 
130 namespace FEValuesViews
131 {
143  template <int dim, int spacedim = dim>
144  class Scalar
145  {
146  public:
152  typedef double value_type;
153 
159  typedef ::Tensor<1, spacedim> gradient_type;
160 
166  typedef ::Tensor<2, spacedim> hessian_type;
167 
173  typedef ::Tensor<3, spacedim> third_derivative_type;
174 
179  template <typename Number>
180  struct OutputType
181  {
186  typedef
187  typename ProductType<Number,
188  typename Scalar<dim, spacedim>::value_type>::type
190 
195  typedef typename ProductType<
196  Number,
198 
203  typedef
204  typename ProductType<Number,
205  typename Scalar<dim, spacedim>::value_type>::type
207 
212  typedef
213  typename ProductType<Number,
216 
221  typedef typename ProductType<
222  Number,
225  };
226 
232  {
242 
251  unsigned int row_index;
252  };
253 
257  Scalar();
258 
264  Scalar(const FEValuesBase<dim, spacedim> &fe_values_base,
265  const unsigned int component);
266 
271  Scalar &
273 
287  value_type
288  value(const unsigned int shape_function, const unsigned int q_point) const;
289 
300  gradient_type
301  gradient(const unsigned int shape_function,
302  const unsigned int q_point) const;
303 
314  hessian_type
315  hessian(const unsigned int shape_function,
316  const unsigned int q_point) const;
317 
328  third_derivative_type
329  third_derivative(const unsigned int shape_function,
330  const unsigned int q_point) const;
331 
349  template <class InputVector>
350  void
352  const InputVector &fe_function,
353  std::vector<typename ProductType<value_type,
354  typename InputVector::value_type>::type>
355  &values) const;
356 
378  template <class InputVector>
379  void
381  const InputVector &dof_values,
382  std::vector<
384  &values) const;
385 
403  template <class InputVector>
404  void
406  const InputVector &fe_function,
407  std::vector<typename ProductType<gradient_type,
408  typename InputVector::value_type>::type>
409  &gradients) const;
410 
414  template <class InputVector>
415  void
417  const InputVector &dof_values,
418  std::vector<
420  &gradients) const;
421 
439  template <class InputVector>
440  void
442  const InputVector &fe_function,
443  std::vector<typename ProductType<hessian_type,
444  typename InputVector::value_type>::type>
445  &hessians) const;
446 
450  template <class InputVector>
451  void
453  const InputVector &dof_values,
454  std::vector<
456  &hessians) const;
457 
458 
477  template <class InputVector>
478  void
480  const InputVector &fe_function,
481  std::vector<typename ProductType<value_type,
482  typename InputVector::value_type>::type>
483  &laplacians) const;
484 
488  template <class InputVector>
489  void
491  const InputVector &dof_values,
492  std::vector<
494  &laplacians) const;
495 
496 
515  template <class InputVector>
516  void
518  const InputVector &fe_function,
519  std::vector<typename ProductType<third_derivative_type,
520  typename InputVector::value_type>::type>
521  &third_derivatives) const;
522 
526  template <class InputVector>
527  void
529  const InputVector & dof_values,
530  std::vector<typename OutputType<typename InputVector::value_type>::
531  third_derivative_type> &third_derivatives) const;
532 
533 
534  private:
539 
544  const unsigned int component;
545 
549  std::vector<ShapeFunctionData> shape_function_data;
550  };
551 
552 
553 
583  template <int dim, int spacedim = dim>
584  class Vector
585  {
586  public:
592  typedef ::Tensor<1, spacedim> value_type;
593 
602  typedef ::Tensor<2, spacedim> gradient_type;
603 
614  typedef ::SymmetricTensor<2, spacedim> symmetric_gradient_type;
615 
621  typedef double divergence_type;
622 
629  typedef typename ::internal::CurlType<spacedim>::type curl_type;
630 
636  typedef ::Tensor<3, spacedim> hessian_type;
637 
643  typedef ::Tensor<4, spacedim> third_derivative_type;
644 
649  template <typename Number>
650  struct OutputType
651  {
656  typedef
657  typename ProductType<Number,
658  typename Vector<dim, spacedim>::value_type>::type
660 
665  typedef typename ProductType<
666  Number,
668 
673  typedef typename ProductType<
674  Number,
677 
682  typedef typename ProductType<
683  Number,
685 
690  typedef
691  typename ProductType<Number,
692  typename Vector<dim, spacedim>::value_type>::type
694 
699  typedef
700  typename ProductType<Number,
701  typename Vector<dim, spacedim>::curl_type>::type
703 
708  typedef
709  typename ProductType<Number,
712 
717  typedef typename ProductType<
718  Number,
721  };
722 
728  {
738 
748  unsigned int row_index[spacedim];
749 
759  unsigned int single_nonzero_component_index;
760  };
761 
765  Vector();
766 
775  Vector(const FEValuesBase<dim, spacedim> &fe_values_base,
776  const unsigned int first_vector_component);
777 
782  Vector &
784 
801  value_type
802  value(const unsigned int shape_function, const unsigned int q_point) const;
803 
817  gradient_type
818  gradient(const unsigned int shape_function,
819  const unsigned int q_point) const;
820 
836  symmetric_gradient_type
837  symmetric_gradient(const unsigned int shape_function,
838  const unsigned int q_point) const;
839 
850  divergence_type
851  divergence(const unsigned int shape_function,
852  const unsigned int q_point) const;
853 
874  curl_type
875  curl(const unsigned int shape_function, const unsigned int q_point) const;
876 
887  hessian_type
888  hessian(const unsigned int shape_function,
889  const unsigned int q_point) const;
890 
901  third_derivative_type
902  third_derivative(const unsigned int shape_function,
903  const unsigned int q_point) const;
904 
922  template <class InputVector>
923  void
925  const InputVector &fe_function,
926  std::vector<typename ProductType<value_type,
927  typename InputVector::value_type>::type>
928  &values) const;
929 
951  template <class InputVector>
952  void
954  const InputVector &dof_values,
955  std::vector<
957  &values) const;
958 
976  template <class InputVector>
977  void
979  const InputVector &fe_function,
980  std::vector<typename ProductType<gradient_type,
981  typename InputVector::value_type>::type>
982  &gradients) const;
983 
987  template <class InputVector>
988  void
990  const InputVector &dof_values,
991  std::vector<
993  &gradients) const;
994 
1018  template <class InputVector>
1019  void
1021  const InputVector &fe_function,
1022  std::vector<typename ProductType<symmetric_gradient_type,
1023  typename InputVector::value_type>::type>
1024  &symmetric_gradients) const;
1025 
1029  template <class InputVector>
1030  void
1032  const InputVector & dof_values,
1033  std::vector<typename OutputType<typename InputVector::value_type>::
1034  symmetric_gradient_type> &symmetric_gradients) const;
1035 
1054  template <class InputVector>
1055  void
1057  const InputVector &fe_function,
1058  std::vector<typename ProductType<divergence_type,
1059  typename InputVector::value_type>::type>
1060  &divergences) const;
1061 
1065  template <class InputVector>
1066  void
1068  const InputVector &dof_values,
1069  std::vector<
1071  &divergences) const;
1072 
1091  template <class InputVector>
1092  void
1094  const InputVector &fe_function,
1095  std::vector<
1096  typename ProductType<curl_type, typename InputVector::value_type>::type>
1097  &curls) const;
1098 
1102  template <class InputVector>
1103  void
1105  const InputVector &dof_values,
1106  std::vector<
1108  &curls) const;
1109 
1127  template <class InputVector>
1128  void
1130  const InputVector &fe_function,
1131  std::vector<typename ProductType<hessian_type,
1132  typename InputVector::value_type>::type>
1133  &hessians) const;
1134 
1138  template <class InputVector>
1139  void
1141  const InputVector &dof_values,
1142  std::vector<
1144  &hessians) const;
1145 
1164  template <class InputVector>
1165  void
1167  const InputVector &fe_function,
1168  std::vector<typename ProductType<value_type,
1169  typename InputVector::value_type>::type>
1170  &laplacians) const;
1171 
1175  template <class InputVector>
1176  void
1178  const InputVector &dof_values,
1179  std::vector<
1181  &laplacians) const;
1182 
1201  template <class InputVector>
1202  void
1204  const InputVector &fe_function,
1205  std::vector<typename ProductType<third_derivative_type,
1206  typename InputVector::value_type>::type>
1207  &third_derivatives) const;
1208 
1212  template <class InputVector>
1213  void
1215  const InputVector & dof_values,
1216  std::vector<typename OutputType<typename InputVector::value_type>::
1217  third_derivative_type> &third_derivatives) const;
1218 
1219  private:
1224 
1229  const unsigned int first_vector_component;
1230 
1234  std::vector<ShapeFunctionData> shape_function_data;
1235  };
1236 
1237 
1238  template <int rank, int dim, int spacedim = dim>
1239  class SymmetricTensor;
1240 
1265  template <int dim, int spacedim>
1266  class SymmetricTensor<2, dim, spacedim>
1267  {
1268  public:
1275  typedef ::SymmetricTensor<2, spacedim> value_type;
1276 
1286  typedef ::Tensor<1, spacedim> divergence_type;
1287 
1292  template <typename Number>
1293  struct OutputType
1294  {
1299  typedef typename ProductType<
1300  Number,
1303 
1308  typedef typename ProductType<
1309  Number,
1312  };
1313 
1318  struct ShapeFunctionData
1319  {
1328  bool is_nonzero_shape_function_component
1329  [value_type::n_independent_components];
1330 
1340  unsigned int row_index[value_type::n_independent_components];
1341 
1351 
1356  };
1357 
1361  SymmetricTensor();
1362 
1372  SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1373  const unsigned int first_tensor_component);
1374 
1379  SymmetricTensor &
1380  operator=(const SymmetricTensor<2, dim, spacedim> &);
1381 
1399  value_type
1400  value(const unsigned int shape_function, const unsigned int q_point) const;
1401 
1415  divergence_type
1416  divergence(const unsigned int shape_function,
1417  const unsigned int q_point) const;
1418 
1436  template <class InputVector>
1437  void
1438  get_function_values(
1439  const InputVector &fe_function,
1440  std::vector<typename ProductType<value_type,
1441  typename InputVector::value_type>::type>
1442  &values) const;
1443 
1465  template <class InputVector>
1466  void
1467  get_function_values_from_local_dof_values(
1468  const InputVector &dof_values,
1469  std::vector<
1471  &values) const;
1472 
1494  template <class InputVector>
1495  void
1496  get_function_divergences(
1497  const InputVector &fe_function,
1498  std::vector<typename ProductType<divergence_type,
1499  typename InputVector::value_type>::type>
1500  &divergences) const;
1501 
1505  template <class InputVector>
1506  void
1507  get_function_divergences_from_local_dof_values(
1508  const InputVector &dof_values,
1509  std::vector<
1511  &divergences) const;
1512 
1513  private:
1518 
1523  const unsigned int first_tensor_component;
1524 
1528  std::vector<ShapeFunctionData> shape_function_data;
1529  };
1530 
1531 
1532  template <int rank, int dim, int spacedim = dim>
1533  class Tensor;
1534 
1555  template <int dim, int spacedim>
1556  class Tensor<2, dim, spacedim>
1557  {
1558  public:
1563  typedef ::Tensor<2, spacedim> value_type;
1564 
1568  typedef ::Tensor<1, spacedim> divergence_type;
1569 
1574  typedef ::Tensor<3, spacedim> gradient_type;
1575 
1580  template <typename Number>
1581  struct OutputType
1582  {
1587  typedef typename ProductType<
1588  Number,
1590 
1595  typedef typename ProductType<
1596  Number,
1599 
1604  typedef typename ProductType<
1605  Number,
1607  };
1608 
1613  struct ShapeFunctionData
1614  {
1623  bool is_nonzero_shape_function_component
1624  [value_type::n_independent_components];
1625 
1635  unsigned int row_index[value_type::n_independent_components];
1636 
1646 
1651  };
1652 
1656  Tensor();
1657 
1667  Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1668  const unsigned int first_tensor_component);
1669 
1670 
1675  Tensor &
1676  operator=(const Tensor<2, dim, spacedim> &);
1677 
1694  value_type
1695  value(const unsigned int shape_function, const unsigned int q_point) const;
1696 
1710  divergence_type
1711  divergence(const unsigned int shape_function,
1712  const unsigned int q_point) const;
1713 
1727  gradient_type
1728  gradient(const unsigned int shape_function,
1729  const unsigned int q_point) const;
1730 
1748  template <class InputVector>
1749  void
1750  get_function_values(
1751  const InputVector &fe_function,
1752  std::vector<typename ProductType<value_type,
1753  typename InputVector::value_type>::type>
1754  &values) const;
1755 
1777  template <class InputVector>
1778  void
1779  get_function_values_from_local_dof_values(
1780  const InputVector &dof_values,
1781  std::vector<
1783  &values) const;
1784 
1806  template <class InputVector>
1807  void
1808  get_function_divergences(
1809  const InputVector &fe_function,
1810  std::vector<typename ProductType<divergence_type,
1811  typename InputVector::value_type>::type>
1812  &divergences) const;
1813 
1817  template <class InputVector>
1818  void
1819  get_function_divergences_from_local_dof_values(
1820  const InputVector &dof_values,
1821  std::vector<
1823  &divergences) const;
1824 
1841  template <class InputVector>
1842  void
1843  get_function_gradients(
1844  const InputVector &fe_function,
1845  std::vector<typename ProductType<gradient_type,
1846  typename InputVector::value_type>::type>
1847  &gradients) const;
1848 
1852  template <class InputVector>
1853  void
1854  get_function_gradients_from_local_dof_values(
1855  const InputVector &dof_values,
1856  std::vector<
1858  &gradients) const;
1859 
1860  private:
1865 
1870  const unsigned int first_tensor_component;
1871 
1875  std::vector<ShapeFunctionData> shape_function_data;
1876  };
1877 
1878 } // namespace FEValuesViews
1879 
1880 
1881 namespace internal
1882 {
1883  namespace FEValuesViews
1884  {
1892  template <int dim, int spacedim>
1893  struct Cache
1894  {
1899  std::vector<::FEValuesViews::Scalar<dim, spacedim>> scalars;
1900  std::vector<::FEValuesViews::Vector<dim, spacedim>> vectors;
1901  std::vector<::FEValuesViews::SymmetricTensor<2, dim, spacedim>>
1902  symmetric_second_order_tensors;
1903  std::vector<::FEValuesViews::Tensor<2, dim, spacedim>>
1904  second_order_tensors;
1905 
1909  Cache(const FEValuesBase<dim, spacedim> &fe_values);
1910  };
1911  } // namespace FEValuesViews
1912 } // namespace internal
1913 
1914 
1915 
2018 template <int dim, int spacedim>
2019 class FEValuesBase : public Subscriptor
2020 {
2021 public:
2025  static const unsigned int dimension = dim;
2026 
2030  static const unsigned int space_dimension = spacedim;
2031 
2035  const unsigned int n_quadrature_points;
2036 
2042  const unsigned int dofs_per_cell;
2043 
2044 
2052  FEValuesBase(const unsigned int n_q_points,
2053  const unsigned int dofs_per_cell,
2054  const UpdateFlags update_flags,
2057 
2058 
2062  ~FEValuesBase() override;
2063 
2064 
2067 
2068 
2089  const double &
2090  shape_value(const unsigned int function_no,
2091  const unsigned int point_no) const;
2092 
2113  double
2114  shape_value_component(const unsigned int function_no,
2115  const unsigned int point_no,
2116  const unsigned int component) const;
2117 
2143  const Tensor<1, spacedim> &
2144  shape_grad(const unsigned int function_no,
2145  const unsigned int quadrature_point) const;
2146 
2164  shape_grad_component(const unsigned int function_no,
2165  const unsigned int point_no,
2166  const unsigned int component) const;
2167 
2187  const Tensor<2, spacedim> &
2188  shape_hessian(const unsigned int function_no,
2189  const unsigned int point_no) const;
2190 
2208  shape_hessian_component(const unsigned int function_no,
2209  const unsigned int point_no,
2210  const unsigned int component) const;
2211 
2231  const Tensor<3, spacedim> &
2232  shape_3rd_derivative(const unsigned int function_no,
2233  const unsigned int point_no) const;
2234 
2252  shape_3rd_derivative_component(const unsigned int function_no,
2253  const unsigned int point_no,
2254  const unsigned int component) const;
2255 
2257 
2259 
2298  template <class InputVector>
2299  void
2301  const InputVector & fe_function,
2302  std::vector<typename InputVector::value_type> &values) const;
2303 
2317  template <class InputVector>
2318  void
2320  const InputVector & fe_function,
2321  std::vector<Vector<typename InputVector::value_type>> &values) const;
2322 
2341  template <class InputVector>
2342  void
2344  const InputVector & fe_function,
2345  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2346  std::vector<typename InputVector::value_type> &values) const;
2347 
2369  template <class InputVector>
2370  void
2372  const InputVector & fe_function,
2373  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2374  std::vector<Vector<typename InputVector::value_type>> &values) const;
2375 
2376 
2407  template <class InputVector>
2408  void
2410  const InputVector & fe_function,
2411  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2412  VectorSlice<std::vector<std::vector<typename InputVector::value_type>>>
2413  values,
2414  const bool quadrature_points_fastest) const;
2415 
2417 
2419 
2458  template <class InputVector>
2459  void
2461  const InputVector &fe_function,
2463  &gradients) const;
2464 
2481  template <class InputVector>
2482  void
2484  const InputVector &fe_function,
2485  std::vector<
2487  &gradients) const;
2488 
2495  template <class InputVector>
2496  void
2498  const InputVector & fe_function,
2499  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2501  &gradients) const;
2502 
2509  template <class InputVector>
2510  void
2512  const InputVector & fe_function,
2513  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2514  VectorSlice<std::vector<
2516  gradients,
2517  bool quadrature_points_fastest = false) const;
2518 
2520 
2523 
2563  template <class InputVector>
2564  void
2566  const InputVector &fe_function,
2568  &hessians) const;
2569 
2587  template <class InputVector>
2588  void
2590  const InputVector &fe_function,
2591  std::vector<
2593  & hessians,
2594  bool quadrature_points_fastest = false) const;
2595 
2600  template <class InputVector>
2601  void
2603  const InputVector & fe_function,
2604  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2606  &hessians) const;
2607 
2614  template <class InputVector>
2615  void
2617  const InputVector & fe_function,
2618  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2619  VectorSlice<std::vector<
2621  hessians,
2622  bool quadrature_points_fastest = false) const;
2623 
2666  template <class InputVector>
2667  void
2669  const InputVector & fe_function,
2670  std::vector<typename InputVector::value_type> &laplacians) const;
2671 
2691  template <class InputVector>
2692  void
2694  const InputVector & fe_function,
2695  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2696 
2703  template <class InputVector>
2704  void
2706  const InputVector & fe_function,
2707  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2708  std::vector<typename InputVector::value_type> &laplacians) const;
2709 
2716  template <class InputVector>
2717  void
2719  const InputVector & fe_function,
2720  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2721  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2722 
2729  template <class InputVector>
2730  void
2732  const InputVector & fe_function,
2733  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2734  std::vector<std::vector<typename InputVector::value_type>> & laplacians,
2735  bool quadrature_points_fastest = false) const;
2736 
2738 
2740 
2781  template <class InputVector>
2782  void
2784  const InputVector &fe_function,
2786  &third_derivatives) const;
2787 
2806  template <class InputVector>
2807  void
2809  const InputVector &fe_function,
2810  std::vector<
2812  & third_derivatives,
2813  bool quadrature_points_fastest = false) const;
2814 
2819  template <class InputVector>
2820  void
2822  const InputVector & fe_function,
2823  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2825  &third_derivatives) const;
2826 
2833  template <class InputVector>
2834  void
2836  const InputVector & fe_function,
2837  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2838  VectorSlice<std::vector<
2840  third_derivatives,
2841  bool quadrature_points_fastest = false) const;
2843 
2845 
2846 
2852  const Point<spacedim> &
2853  quadrature_point(const unsigned int q) const;
2854 
2860  const std::vector<Point<spacedim>> &
2861  get_quadrature_points() const;
2862 
2878  double
2879  JxW(const unsigned int quadrature_point) const;
2880 
2884  const std::vector<double> &
2885  get_JxW_values() const;
2886 
2894  jacobian(const unsigned int quadrature_point) const;
2895 
2902  const std::vector<DerivativeForm<1, dim, spacedim>> &
2903  get_jacobians() const;
2904 
2913  jacobian_grad(const unsigned int quadrature_point) const;
2914 
2921  const std::vector<DerivativeForm<2, dim, spacedim>> &
2922  get_jacobian_grads() const;
2923 
2932  const Tensor<3, spacedim> &
2933  jacobian_pushed_forward_grad(const unsigned int quadrature_point) const;
2934 
2941  const std::vector<Tensor<3, spacedim>> &
2943 
2952  jacobian_2nd_derivative(const unsigned int quadrature_point) const;
2953 
2960  const std::vector<DerivativeForm<3, dim, spacedim>> &
2962 
2972  const Tensor<4, spacedim> &
2974  const unsigned int quadrature_point) const;
2975 
2982  const std::vector<Tensor<4, spacedim>> &
2984 
2994  jacobian_3rd_derivative(const unsigned int quadrature_point) const;
2995 
3002  const std::vector<DerivativeForm<4, dim, spacedim>> &
3004 
3014  const Tensor<5, spacedim> &
3016  const unsigned int quadrature_point) const;
3017 
3024  const std::vector<Tensor<5, spacedim>> &
3026 
3034  inverse_jacobian(const unsigned int quadrature_point) const;
3035 
3042  const std::vector<DerivativeForm<1, spacedim, dim>> &
3043  get_inverse_jacobians() const;
3044 
3058  const Tensor<1, spacedim> &
3059  normal_vector(const unsigned int i) const;
3060 
3071  DEAL_II_DEPRECATED
3072  const std::vector<Tensor<1, spacedim>> &
3073  get_all_normal_vectors() const;
3074 
3082  const std::vector<Tensor<1, spacedim>> &
3083  get_normal_vectors() const;
3084 
3086 
3088 
3089 
3099  operator[](const FEValuesExtractors::Scalar &scalar) const;
3100 
3110  operator[](const FEValuesExtractors::Vector &vector) const;
3111 
3123 
3124 
3134  operator[](const FEValuesExtractors::Tensor<2> &tensor) const;
3135 
3137 
3139 
3140 
3144  const Mapping<dim, spacedim> &
3145  get_mapping() const;
3146 
3151  get_fe() const;
3152 
3156  UpdateFlags
3157  get_update_flags() const;
3158 
3163  get_cell() const;
3164 
3171  get_cell_similarity() const;
3172 
3177  std::size_t
3178  memory_consumption() const;
3180 
3181 
3190  std::string,
3191  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3192  << "object for which this kind of information has not been computed. What "
3193  << "information these objects compute is determined by the update_* flags you "
3194  << "pass to the constructor. Here, the operation you are attempting requires "
3195  << "the <" << arg1
3196  << "> flag to be set, but it was apparently not specified "
3197  << "upon construction.");
3198 
3207  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3208  "to the DoFHandler that provided the cell iterator do not match.");
3215  int,
3216  << "The shape function with index " << arg1
3217  << " is not primitive, i.e. it is vector-valued and "
3218  << "has more than one non-zero vector component. This "
3219  << "function cannot be called for these shape functions. "
3220  << "Maybe you want to use the same function with the "
3221  << "_component suffix?");
3222 
3231  "The given FiniteElement is not a primitive element but the requested operation "
3232  "only works for those. See FiniteElement::is_primitive() for more information.");
3233 
3234 protected:
3263  class CellIteratorBase;
3264 
3269  template <typename CI>
3272 
3278  std::unique_ptr<const CellIteratorBase> present_cell;
3279 
3287  boost::signals2::connection tria_listener_refinement;
3288 
3296  boost::signals2::connection tria_listener_mesh_transform;
3297 
3303  void
3305 
3315  void
3317  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3318 
3324 
3330  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3332 
3339 
3340 
3348 
3354  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3356 
3362  spacedim>
3364 
3365 
3370 
3379  UpdateFlags
3380  compute_update_flags(const UpdateFlags update_flags) const;
3381 
3388 
3394  void
3396  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3397 
3398 private:
3403  FEValuesBase(const FEValuesBase &);
3404 
3409  FEValuesBase &
3410  operator=(const FEValuesBase &);
3411 
3416 
3421  template <int, int>
3423  template <int, int>
3424  friend class FEValuesViews::Vector;
3425  template <int, int, int>
3426  friend class FEValuesViews::SymmetricTensor;
3427  template <int, int, int>
3428  friend class FEValuesViews::Tensor;
3429 };
3430 
3431 
3432 
3443 template <int dim, int spacedim = dim>
3444 class FEValues : public FEValuesBase<dim, spacedim>
3445 {
3446 public:
3451  static const unsigned int integral_dimension = dim;
3452 
3459  const Quadrature<dim> & quadrature,
3460  const UpdateFlags update_flags);
3461 
3468  const Quadrature<dim> & quadrature,
3469  const UpdateFlags update_flags);
3470 
3477  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3478  void
3479  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3480  level_dof_access>> &cell);
3481 
3495  void
3497 
3502  const Quadrature<dim> &
3503  get_quadrature() const;
3504 
3509  std::size_t
3510  memory_consumption() const;
3511 
3526  const FEValues<dim, spacedim> &
3527  get_present_fe_values() const;
3528 
3529 private:
3534 
3538  void
3540 
3547  void
3548  do_reinit();
3549 };
3550 
3551 
3562 template <int dim, int spacedim = dim>
3563 class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
3564 {
3565 public:
3570  static const unsigned int integral_dimension = dim - 1;
3571 
3583  FEFaceValuesBase(const unsigned int n_q_points,
3584  const unsigned int dofs_per_cell,
3585  const UpdateFlags update_flags,
3589 
3597  const Tensor<1, spacedim> &
3598  boundary_form(const unsigned int i) const;
3599 
3606  const std::vector<Tensor<1, spacedim>> &
3607  get_boundary_forms() const;
3608 
3613  unsigned int
3614  get_face_index() const;
3615 
3620  const Quadrature<dim - 1> &
3621  get_quadrature() const;
3622 
3627  std::size_t
3628  memory_consumption() const;
3629 
3630 protected:
3635  unsigned int present_face_index;
3636 
3640  const Quadrature<dim - 1> quadrature;
3641 };
3642 
3643 
3644 
3659 template <int dim, int spacedim = dim>
3660 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
3661 {
3662 public:
3667  static const unsigned int dimension = dim;
3668 
3669  static const unsigned int space_dimension = spacedim;
3670 
3675  static const unsigned int integral_dimension = dim - 1;
3676 
3684  const UpdateFlags update_flags);
3685 
3693  const UpdateFlags update_flags);
3694 
3699  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3700  void
3701  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3702  level_dof_access>> &cell,
3703  const unsigned int face_no);
3704 
3718  void
3720  const unsigned int face_no);
3721 
3737  get_present_fe_values() const;
3738 
3739 private:
3743  void
3745 
3752  void
3753  do_reinit(const unsigned int face_no);
3754 };
3755 
3756 
3774 template <int dim, int spacedim = dim>
3775 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
3776 {
3777 public:
3781  static const unsigned int dimension = dim;
3782 
3786  static const unsigned int space_dimension = spacedim;
3787 
3792  static const unsigned int integral_dimension = dim - 1;
3793 
3800  const Quadrature<dim - 1> & face_quadrature,
3801  const UpdateFlags update_flags);
3802 
3809  const Quadrature<dim - 1> & face_quadrature,
3810  const UpdateFlags update_flags);
3811 
3818  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3819  void
3820  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3821  level_dof_access>> &cell,
3822  const unsigned int face_no,
3823  const unsigned int subface_no);
3824 
3838  void
3840  const unsigned int face_no,
3841  const unsigned int subface_no);
3842 
3858  get_present_fe_values() const;
3859 
3866 
3873 
3874 private:
3878  void
3880 
3887  void
3888  do_reinit(const unsigned int face_no, const unsigned int subface_no);
3889 };
3890 
3891 
3892 #ifndef DOXYGEN
3893 
3894 
3895 /*------------------------ Inline functions: namespace FEValuesViews --------*/
3896 
3897 namespace FEValuesViews
3898 {
3899  template <int dim, int spacedim>
3900  inline typename Scalar<dim, spacedim>::value_type
3901  Scalar<dim, spacedim>::value(const unsigned int shape_function,
3902  const unsigned int q_point) const
3903  {
3904  Assert(shape_function < fe_values->fe->dofs_per_cell,
3905  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3906  Assert(
3907  fe_values->update_flags & update_values,
3909  "update_values"))));
3910 
3911  // an adaptation of the FEValuesBase::shape_value_component function
3912  // except that here we know the component as fixed and we have
3913  // pre-computed and cached a bunch of information. See the comments there.
3914  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3915  return fe_values->finite_element_output.shape_values(
3916  shape_function_data[shape_function].row_index, q_point);
3917  else
3918  return 0;
3919  }
3920 
3921 
3922 
3923  template <int dim, int spacedim>
3924  inline typename Scalar<dim, spacedim>::gradient_type
3925  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
3926  const unsigned int q_point) const
3927  {
3928  Assert(shape_function < fe_values->fe->dofs_per_cell,
3929  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3930  Assert(fe_values->update_flags & update_gradients,
3932  "update_gradients")));
3933 
3934  // an adaptation of the FEValuesBase::shape_grad_component
3935  // function except that here we know the component as fixed and we have
3936  // pre-computed and cached a bunch of information. See the comments there.
3937  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3938  return fe_values->finite_element_output
3939  .shape_gradients[shape_function_data[shape_function].row_index]
3940  [q_point];
3941  else
3942  return gradient_type();
3943  }
3944 
3945 
3946 
3947  template <int dim, int spacedim>
3948  inline typename Scalar<dim, spacedim>::hessian_type
3949  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
3950  const unsigned int q_point) const
3951  {
3952  Assert(shape_function < fe_values->fe->dofs_per_cell,
3953  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3954  Assert(fe_values->update_flags & update_hessians,
3956  "update_hessians")));
3957 
3958  // an adaptation of the FEValuesBase::shape_hessian_component
3959  // function except that here we know the component as fixed and we have
3960  // pre-computed and cached a bunch of information. See the comments there.
3961  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3962  return fe_values->finite_element_output
3963  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
3964  else
3965  return hessian_type();
3966  }
3967 
3968 
3969 
3970  template <int dim, int spacedim>
3971  inline typename Scalar<dim, spacedim>::third_derivative_type
3972  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
3973  const unsigned int q_point) const
3974  {
3975  Assert(shape_function < fe_values->fe->dofs_per_cell,
3976  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3977  Assert(fe_values->update_flags & update_3rd_derivatives,
3979  "update_3rd_derivatives")));
3980 
3981  // an adaptation of the FEValuesBase::shape_3rdderivative_component
3982  // function except that here we know the component as fixed and we have
3983  // pre-computed and cached a bunch of information. See the comments there.
3984  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3985  return fe_values->finite_element_output
3986  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
3987  [q_point];
3988  else
3989  return third_derivative_type();
3990  }
3991 
3992 
3993 
3994  template <int dim, int spacedim>
3995  inline typename Vector<dim, spacedim>::value_type
3996  Vector<dim, spacedim>::value(const unsigned int shape_function,
3997  const unsigned int q_point) const
3998  {
3999  Assert(shape_function < fe_values->fe->dofs_per_cell,
4000  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4001  Assert(fe_values->update_flags & update_values,
4003  "update_values")));
4004 
4005  // same as for the scalar case except that we have one more index
4006  const int snc =
4007  shape_function_data[shape_function].single_nonzero_component;
4008  if (snc == -2)
4009  return value_type();
4010  else if (snc != -1)
4011  {
4012  value_type return_value;
4013  return_value[shape_function_data[shape_function]
4014  .single_nonzero_component_index] =
4015  fe_values->finite_element_output.shape_values(snc, q_point);
4016  return return_value;
4017  }
4018  else
4019  {
4020  value_type return_value;
4021  for (unsigned int d = 0; d < dim; ++d)
4022  if (shape_function_data[shape_function]
4023  .is_nonzero_shape_function_component[d])
4024  return_value[d] = fe_values->finite_element_output.shape_values(
4025  shape_function_data[shape_function].row_index[d], q_point);
4026 
4027  return return_value;
4028  }
4029  }
4030 
4031 
4032 
4033  template <int dim, int spacedim>
4034  inline typename Vector<dim, spacedim>::gradient_type
4035  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4036  const unsigned int q_point) const
4037  {
4038  Assert(shape_function < fe_values->fe->dofs_per_cell,
4039  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4040  Assert(fe_values->update_flags & update_gradients,
4042  "update_gradients")));
4043 
4044  // same as for the scalar case except that we have one more index
4045  const int snc =
4046  shape_function_data[shape_function].single_nonzero_component;
4047  if (snc == -2)
4048  return gradient_type();
4049  else if (snc != -1)
4050  {
4051  gradient_type return_value;
4052  return_value[shape_function_data[shape_function]
4053  .single_nonzero_component_index] =
4054  fe_values->finite_element_output.shape_gradients[snc][q_point];
4055  return return_value;
4056  }
4057  else
4058  {
4059  gradient_type return_value;
4060  for (unsigned int d = 0; d < dim; ++d)
4061  if (shape_function_data[shape_function]
4062  .is_nonzero_shape_function_component[d])
4063  return_value[d] =
4064  fe_values->finite_element_output.shape_gradients
4065  [shape_function_data[shape_function].row_index[d]][q_point];
4066 
4067  return return_value;
4068  }
4069  }
4070 
4071 
4072 
4073  template <int dim, int spacedim>
4074  inline typename Vector<dim, spacedim>::divergence_type
4075  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4076  const unsigned int q_point) const
4077  {
4078  // this function works like in the case above
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_gradients,
4083  "update_gradients")));
4084 
4085  // same as for the scalar case except that we have one more index
4086  const int snc =
4087  shape_function_data[shape_function].single_nonzero_component;
4088  if (snc == -2)
4089  return divergence_type();
4090  else if (snc != -1)
4091  return fe_values->finite_element_output
4092  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4093  .single_nonzero_component_index];
4094  else
4095  {
4096  divergence_type return_value = 0;
4097  for (unsigned int d = 0; d < dim; ++d)
4098  if (shape_function_data[shape_function]
4099  .is_nonzero_shape_function_component[d])
4100  return_value +=
4101  fe_values->finite_element_output.shape_gradients
4102  [shape_function_data[shape_function].row_index[d]][q_point][d];
4103 
4104  return return_value;
4105  }
4106  }
4107 
4108 
4109 
4110  template <int dim, int spacedim>
4111  inline typename Vector<dim, spacedim>::curl_type
4112  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4113  const unsigned int q_point) const
4114  {
4115  // this function works like in the case above
4116 
4117  Assert(shape_function < fe_values->fe->dofs_per_cell,
4118  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4119  Assert(fe_values->update_flags & update_gradients,
4121  "update_gradients")));
4122  // same as for the scalar case except that we have one more index
4123  const int snc =
4124  shape_function_data[shape_function].single_nonzero_component;
4125 
4126  if (snc == -2)
4127  return curl_type();
4128 
4129  else
4130  switch (dim)
4131  {
4132  case 1:
4133  {
4134  Assert(false,
4135  ExcMessage(
4136  "Computing the curl in 1d is not a useful operation"));
4137  return curl_type();
4138  }
4139 
4140  case 2:
4141  {
4142  if (snc != -1)
4143  {
4144  curl_type return_value;
4145 
4146  // the single nonzero component can only be zero or one in 2d
4147  if (shape_function_data[shape_function]
4148  .single_nonzero_component_index == 0)
4149  return_value[0] =
4150  -1.0 * fe_values->finite_element_output
4151  .shape_gradients[snc][q_point][1];
4152  else
4153  return_value[0] = fe_values->finite_element_output
4154  .shape_gradients[snc][q_point][0];
4155 
4156  return return_value;
4157  }
4158 
4159  else
4160  {
4161  curl_type return_value;
4162 
4163  return_value[0] = 0.0;
4164 
4165  if (shape_function_data[shape_function]
4166  .is_nonzero_shape_function_component[0])
4167  return_value[0] -=
4168  fe_values->finite_element_output
4169  .shape_gradients[shape_function_data[shape_function]
4170  .row_index[0]][q_point][1];
4171 
4172  if (shape_function_data[shape_function]
4173  .is_nonzero_shape_function_component[1])
4174  return_value[0] +=
4175  fe_values->finite_element_output
4176  .shape_gradients[shape_function_data[shape_function]
4177  .row_index[1]][q_point][0];
4178 
4179  return return_value;
4180  }
4181  }
4182 
4183  case 3:
4184  {
4185  if (snc != -1)
4186  {
4187  curl_type return_value;
4188 
4189  switch (shape_function_data[shape_function]
4190  .single_nonzero_component_index)
4191  {
4192  case 0:
4193  {
4194  return_value[0] = 0;
4195  return_value[1] = fe_values->finite_element_output
4196  .shape_gradients[snc][q_point][2];
4197  return_value[2] =
4198  -1.0 * fe_values->finite_element_output
4199  .shape_gradients[snc][q_point][1];
4200  return return_value;
4201  }
4202 
4203  case 1:
4204  {
4205  return_value[0] =
4206  -1.0 * fe_values->finite_element_output
4207  .shape_gradients[snc][q_point][2];
4208  return_value[1] = 0;
4209  return_value[2] = fe_values->finite_element_output
4210  .shape_gradients[snc][q_point][0];
4211  return return_value;
4212  }
4213 
4214  default:
4215  {
4216  return_value[0] = fe_values->finite_element_output
4217  .shape_gradients[snc][q_point][1];
4218  return_value[1] =
4219  -1.0 * fe_values->finite_element_output
4220  .shape_gradients[snc][q_point][0];
4221  return_value[2] = 0;
4222  return return_value;
4223  }
4224  }
4225  }
4226 
4227  else
4228  {
4229  curl_type return_value;
4230 
4231  for (unsigned int i = 0; i < dim; ++i)
4232  return_value[i] = 0.0;
4233 
4234  if (shape_function_data[shape_function]
4235  .is_nonzero_shape_function_component[0])
4236  {
4237  return_value[1] +=
4238  fe_values->finite_element_output
4239  .shape_gradients[shape_function_data[shape_function]
4240  .row_index[0]][q_point][2];
4241  return_value[2] -=
4242  fe_values->finite_element_output
4243  .shape_gradients[shape_function_data[shape_function]
4244  .row_index[0]][q_point][1];
4245  }
4246 
4247  if (shape_function_data[shape_function]
4248  .is_nonzero_shape_function_component[1])
4249  {
4250  return_value[0] -=
4251  fe_values->finite_element_output
4252  .shape_gradients[shape_function_data[shape_function]
4253  .row_index[1]][q_point][2];
4254  return_value[2] +=
4255  fe_values->finite_element_output
4256  .shape_gradients[shape_function_data[shape_function]
4257  .row_index[1]][q_point][0];
4258  }
4259 
4260  if (shape_function_data[shape_function]
4261  .is_nonzero_shape_function_component[2])
4262  {
4263  return_value[0] +=
4264  fe_values->finite_element_output
4265  .shape_gradients[shape_function_data[shape_function]
4266  .row_index[2]][q_point][1];
4267  return_value[1] -=
4268  fe_values->finite_element_output
4269  .shape_gradients[shape_function_data[shape_function]
4270  .row_index[2]][q_point][0];
4271  }
4272 
4273  return return_value;
4274  }
4275  }
4276  }
4277  // should not end up here
4278  Assert(false, ExcInternalError());
4279  return curl_type();
4280  }
4281 
4282 
4283 
4284  template <int dim, int spacedim>
4285  inline typename Vector<dim, spacedim>::hessian_type
4286  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
4287  const unsigned int q_point) const
4288  {
4289  // this function works like in the case above
4290  Assert(shape_function < fe_values->fe->dofs_per_cell,
4291  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4292  Assert(fe_values->update_flags & update_hessians,
4294  "update_hessians")));
4295 
4296  // same as for the scalar case except that we have one more index
4297  const int snc =
4298  shape_function_data[shape_function].single_nonzero_component;
4299  if (snc == -2)
4300  return hessian_type();
4301  else if (snc != -1)
4302  {
4303  hessian_type return_value;
4304  return_value[shape_function_data[shape_function]
4305  .single_nonzero_component_index] =
4306  fe_values->finite_element_output.shape_hessians[snc][q_point];
4307  return return_value;
4308  }
4309  else
4310  {
4311  hessian_type return_value;
4312  for (unsigned int d = 0; d < dim; ++d)
4313  if (shape_function_data[shape_function]
4314  .is_nonzero_shape_function_component[d])
4315  return_value[d] =
4316  fe_values->finite_element_output.shape_hessians
4317  [shape_function_data[shape_function].row_index[d]][q_point];
4318 
4319  return return_value;
4320  }
4321  }
4322 
4323 
4324 
4325  template <int dim, int spacedim>
4326  inline typename Vector<dim, spacedim>::third_derivative_type
4327  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
4328  const unsigned int q_point) const
4329  {
4330  // this function works like in the case above
4331  Assert(shape_function < fe_values->fe->dofs_per_cell,
4332  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4333  Assert(fe_values->update_flags & update_3rd_derivatives,
4335  "update_3rd_derivatives")));
4336 
4337  // same as for the scalar case except that we have one more index
4338  const int snc =
4339  shape_function_data[shape_function].single_nonzero_component;
4340  if (snc == -2)
4341  return third_derivative_type();
4342  else if (snc != -1)
4343  {
4344  third_derivative_type return_value;
4345  return_value[shape_function_data[shape_function]
4346  .single_nonzero_component_index] =
4347  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
4348  return return_value;
4349  }
4350  else
4351  {
4352  third_derivative_type return_value;
4353  for (unsigned int d = 0; d < dim; ++d)
4354  if (shape_function_data[shape_function]
4355  .is_nonzero_shape_function_component[d])
4356  return_value[d] =
4357  fe_values->finite_element_output.shape_3rd_derivatives
4358  [shape_function_data[shape_function].row_index[d]][q_point];
4359 
4360  return return_value;
4361  }
4362  }
4363 
4364 
4365 
4366  namespace internal
4367  {
4372  inline ::SymmetricTensor<2, 1>
4373  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
4374  {
4375  Assert(n < 1, ExcIndexRange(n, 0, 1));
4376  (void)n;
4377 
4378  const double array[1] = {t[0]};
4379  return ::SymmetricTensor<2, 1>(array);
4380  }
4381 
4382 
4383 
4384  inline ::SymmetricTensor<2, 2>
4385  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
4386  {
4387  switch (n)
4388  {
4389  case 0:
4390  {
4391  const double array[3] = {t[0], 0, t[1] / 2};
4392  return ::SymmetricTensor<2, 2>(array);
4393  }
4394  case 1:
4395  {
4396  const double array[3] = {0, t[1], t[0] / 2};
4397  return ::SymmetricTensor<2, 2>(array);
4398  }
4399  default:
4400  {
4401  Assert(false, ExcIndexRange(n, 0, 2));
4402  return ::SymmetricTensor<2, 2>();
4403  }
4404  }
4405  }
4406 
4407 
4408 
4409  inline ::SymmetricTensor<2, 3>
4410  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
4411  {
4412  switch (n)
4413  {
4414  case 0:
4415  {
4416  const double array[6] = {t[0], 0, 0, t[1] / 2, t[2] / 2, 0};
4417  return ::SymmetricTensor<2, 3>(array);
4418  }
4419  case 1:
4420  {
4421  const double array[6] = {0, t[1], 0, t[0] / 2, 0, t[2] / 2};
4422  return ::SymmetricTensor<2, 3>(array);
4423  }
4424  case 2:
4425  {
4426  const double array[6] = {0, 0, t[2], 0, t[0] / 2, t[1] / 2};
4427  return ::SymmetricTensor<2, 3>(array);
4428  }
4429  default:
4430  {
4431  Assert(false, ExcIndexRange(n, 0, 3));
4432  return ::SymmetricTensor<2, 3>();
4433  }
4434  }
4435  }
4436  } // namespace internal
4437 
4438 
4439 
4440  template <int dim, int spacedim>
4441  inline typename Vector<dim, spacedim>::symmetric_gradient_type
4442  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
4443  const unsigned int q_point) const
4444  {
4445  Assert(shape_function < fe_values->fe->dofs_per_cell,
4446  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4447  Assert(fe_values->update_flags & update_gradients,
4449  "update_gradients")));
4450 
4451  // same as for the scalar case except that we have one more index
4452  const int snc =
4453  shape_function_data[shape_function].single_nonzero_component;
4454  if (snc == -2)
4455  return symmetric_gradient_type();
4456  else if (snc != -1)
4457  return internal::symmetrize_single_row(
4458  shape_function_data[shape_function].single_nonzero_component_index,
4459  fe_values->finite_element_output.shape_gradients[snc][q_point]);
4460  else
4461  {
4462  gradient_type return_value;
4463  for (unsigned int d = 0; d < dim; ++d)
4464  if (shape_function_data[shape_function]
4465  .is_nonzero_shape_function_component[d])
4466  return_value[d] =
4467  fe_values->finite_element_output.shape_gradients
4468  [shape_function_data[shape_function].row_index[d]][q_point];
4469 
4470  return symmetrize(return_value);
4471  }
4472  }
4473 
4474 
4475 
4476  template <int dim, int spacedim>
4478  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
4479  const unsigned int q_point) const
4480  {
4481  Assert(shape_function < fe_values->fe->dofs_per_cell,
4482  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4483  Assert(fe_values->update_flags & update_values,
4485  "update_values")));
4486 
4487  // similar to the vector case where we have more then one index and we need
4488  // to convert between unrolled and component indexing for tensors
4489  const int snc =
4490  shape_function_data[shape_function].single_nonzero_component;
4491 
4492  if (snc == -2)
4493  {
4494  // shape function is zero for the selected components
4495  return value_type();
4496  }
4497  else if (snc != -1)
4498  {
4499  value_type return_value;
4500  const unsigned int comp =
4501  shape_function_data[shape_function].single_nonzero_component_index;
4502  return_value[value_type::unrolled_to_component_indices(comp)] =
4503  fe_values->finite_element_output.shape_values(snc, q_point);
4504  return return_value;
4505  }
4506  else
4507  {
4508  value_type return_value;
4509  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
4510  if (shape_function_data[shape_function]
4511  .is_nonzero_shape_function_component[d])
4512  return_value[value_type::unrolled_to_component_indices(d)] =
4513  fe_values->finite_element_output.shape_values(
4514  shape_function_data[shape_function].row_index[d], q_point);
4515  return return_value;
4516  }
4517  }
4518 
4519 
4520 
4521  template <int dim, int spacedim>
4524  const unsigned int shape_function,
4525  const unsigned int q_point) const
4526  {
4527  Assert(shape_function < fe_values->fe->dofs_per_cell,
4528  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4529  Assert(fe_values->update_flags & update_gradients,
4531  "update_gradients")));
4532 
4533  const int snc =
4534  shape_function_data[shape_function].single_nonzero_component;
4535 
4536  if (snc == -2)
4537  {
4538  // shape function is zero for the selected components
4539  return divergence_type();
4540  }
4541  else if (snc != -1)
4542  {
4543  // we have a single non-zero component when the symmetric tensor is
4544  // represented in unrolled form. this implies we potentially have
4545  // two non-zero components when represented in component form! we
4546  // will only have one non-zero entry if the non-zero component lies on
4547  // the diagonal of the tensor.
4548  //
4549  // the divergence of a second-order tensor is a first order tensor.
4550  //
4551  // assume the second-order tensor is A with components A_{ij}. then
4552  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
4553  // entries in the tensorial representation. define the
4554  // divergence as:
4555  // b_i := \dfrac{\partial phi_{ij}}{\partial x_j}.
4556  // (which is incidentally also
4557  // b_j := \dfrac{\partial phi_{ij}}{\partial x_i}).
4558  // In both cases, a sum is implied.
4559  //
4560  // Now, we know the nonzero component in unrolled form: it is indicated
4561  // by 'snc'. we can figure out which tensor components belong to this:
4562  const unsigned int comp =
4563  shape_function_data[shape_function].single_nonzero_component_index;
4564  const unsigned int ii =
4565  value_type::unrolled_to_component_indices(comp)[0];
4566  const unsigned int jj =
4567  value_type::unrolled_to_component_indices(comp)[1];
4568 
4569  // given the form of the divergence above, if ii=jj there is only a
4570  // single nonzero component of the full tensor and the gradient
4571  // equals
4572  // b_ii := \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
4573  // all other entries of 'b' are zero
4574  //
4575  // on the other hand, if ii!=jj, then there are two nonzero entries in
4576  // the full tensor and
4577  // b_ii := \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
4578  // b_jj := \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
4579  // again, all other entries of 'b' are zero
4580  const ::Tensor<1, spacedim> &phi_grad =
4581  fe_values->finite_element_output.shape_gradients[snc][q_point];
4582 
4583  divergence_type return_value;
4584  return_value[ii] = phi_grad[jj];
4585 
4586  if (ii != jj)
4587  return_value[jj] = phi_grad[ii];
4588 
4589  return return_value;
4590  }
4591  else
4592  {
4593  Assert(false, ExcNotImplemented());
4594  divergence_type return_value;
4595  return return_value;
4596  }
4597  }
4598 
4599 
4600 
4601  template <int dim, int spacedim>
4602  inline typename Tensor<2, dim, spacedim>::value_type
4603  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
4604  const unsigned int q_point) const
4605  {
4606  Assert(shape_function < fe_values->fe->dofs_per_cell,
4607  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4608  Assert(fe_values->update_flags & update_values,
4610  "update_values")));
4611 
4612  // similar to the vector case where we have more then one index and we need
4613  // to convert between unrolled and component indexing for tensors
4614  const int snc =
4615  shape_function_data[shape_function].single_nonzero_component;
4616 
4617  if (snc == -2)
4618  {
4619  // shape function is zero for the selected components
4620  return value_type();
4621  }
4622  else if (snc != -1)
4623  {
4624  value_type return_value;
4625  const unsigned int comp =
4626  shape_function_data[shape_function].single_nonzero_component_index;
4627  const TableIndices<2> indices =
4629  return_value[indices] =
4630  fe_values->finite_element_output.shape_values(snc, q_point);
4631  return return_value;
4632  }
4633  else
4634  {
4635  value_type return_value;
4636  for (unsigned int d = 0; d < dim * dim; ++d)
4637  if (shape_function_data[shape_function]
4638  .is_nonzero_shape_function_component[d])
4639  {
4640  const TableIndices<2> indices =
4642  return_value[indices] =
4643  fe_values->finite_element_output.shape_values(
4644  shape_function_data[shape_function].row_index[d], q_point);
4645  }
4646  return return_value;
4647  }
4648  }
4649 
4650 
4651 
4652  template <int dim, int spacedim>
4654  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
4655  const unsigned int q_point) const
4656  {
4657  Assert(shape_function < fe_values->fe->dofs_per_cell,
4658  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4659  Assert(fe_values->update_flags & update_gradients,
4661  "update_gradients")));
4662 
4663  const int snc =
4664  shape_function_data[shape_function].single_nonzero_component;
4665 
4666  if (snc == -2)
4667  {
4668  // shape function is zero for the selected components
4669  return divergence_type();
4670  }
4671  else if (snc != -1)
4672  {
4673  // we have a single non-zero component when the tensor is
4674  // represented in unrolled form.
4675  //
4676  // the divergence of a second-order tensor is a first order tensor.
4677  //
4678  // assume the second-order tensor is A with components A_{ij},
4679  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
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 TableIndices<2> indices =
4687  const unsigned int ii = indices[0];
4688  const unsigned int jj = indices[1];
4689 
4690  const ::Tensor<1, spacedim> &phi_grad =
4691  fe_values->finite_element_output.shape_gradients[snc][q_point];
4692 
4693  divergence_type return_value;
4694  // note that we contract \nabla from the right
4695  return_value[ii] = phi_grad[jj];
4696 
4697  return return_value;
4698  }
4699  else
4700  {
4701  Assert(false, ExcNotImplemented());
4702  divergence_type return_value;
4703  return return_value;
4704  }
4705  }
4706 
4707 
4708 
4709  template <int dim, int spacedim>
4711  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
4712  const unsigned int q_point) const
4713  {
4714  Assert(shape_function < fe_values->fe->dofs_per_cell,
4715  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4716  Assert(fe_values->update_flags & update_gradients,
4718  "update_gradients")));
4719 
4720  const int snc =
4721  shape_function_data[shape_function].single_nonzero_component;
4722 
4723  if (snc == -2)
4724  {
4725  // shape function is zero for the selected components
4726  return gradient_type();
4727  }
4728  else if (snc != -1)
4729  {
4730  // we have a single non-zero component when the tensor is
4731  // represented in unrolled form.
4732  //
4733  // the gradient of a second-order tensor is a third order tensor.
4734  //
4735  // assume the second-order tensor is A with components A_{ij},
4736  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
4737  //
4738  // Now, we know the nonzero component in unrolled form: it is indicated
4739  // by 'snc'. we can figure out which tensor components belong to this:
4740  const unsigned int comp =
4741  shape_function_data[shape_function].single_nonzero_component_index;
4742  const TableIndices<2> indices =
4744  const unsigned int ii = indices[0];
4745  const unsigned int jj = indices[1];
4746 
4747  const ::Tensor<1, spacedim> &phi_grad =
4748  fe_values->finite_element_output.shape_gradients[snc][q_point];
4749 
4750  gradient_type return_value;
4751  return_value[ii][jj] = phi_grad;
4752 
4753  return return_value;
4754  }
4755  else
4756  {
4757  Assert(false, ExcNotImplemented());
4758  gradient_type return_value;
4759  return return_value;
4760  }
4761  }
4762 
4763 } // namespace FEValuesViews
4764 
4765 
4766 
4767 /*---------------------- Inline functions: FEValuesBase ---------------------*/
4768 
4769 
4770 
4771 template <int dim, int spacedim>
4773  operator[](const FEValuesExtractors::Scalar &scalar) const
4774 {
4775  Assert(scalar.component < fe_values_views_cache.scalars.size(),
4776  ExcIndexRange(scalar.component,
4777  0,
4778  fe_values_views_cache.scalars.size()));
4779 
4780  return fe_values_views_cache.scalars[scalar.component];
4781 }
4782 
4783 
4784 
4785 template <int dim, int spacedim>
4787  operator[](const FEValuesExtractors::Vector &vector) const
4788 {
4789  Assert(vector.first_vector_component < fe_values_views_cache.vectors.size(),
4791  0,
4792  fe_values_views_cache.vectors.size()));
4793 
4794  return fe_values_views_cache.vectors[vector.first_vector_component];
4795 }
4796 
4797 
4798 
4799 template <int dim, int spacedim>
4803 {
4804  Assert(
4805  tensor.first_tensor_component <
4806  fe_values_views_cache.symmetric_second_order_tensors.size(),
4808  0,
4809  fe_values_views_cache.symmetric_second_order_tensors.size()));
4810 
4811  return fe_values_views_cache
4812  .symmetric_second_order_tensors[tensor.first_tensor_component];
4813 }
4814 
4815 
4816 
4817 template <int dim, int spacedim>
4820  operator[](const FEValuesExtractors::Tensor<2> &tensor) const
4821 {
4823  fe_values_views_cache.second_order_tensors.size(),
4825  0,
4826  fe_values_views_cache.second_order_tensors.size()));
4827 
4828  return fe_values_views_cache
4829  .second_order_tensors[tensor.first_tensor_component];
4830 }
4831 
4832 
4833 
4834 template <int dim, int spacedim>
4835 inline const double &
4836 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
4837  const unsigned int j) const
4838 {
4839  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4840  Assert(this->update_flags & update_values,
4841  ExcAccessToUninitializedField("update_values"));
4842  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4843  Assert(present_cell.get() != nullptr,
4844  ExcMessage("FEValues object is not reinit'ed to any cell"));
4845  // if the entire FE is primitive,
4846  // then we can take a short-cut:
4847  if (fe->is_primitive())
4848  return this->finite_element_output.shape_values(i, j);
4849  else
4850  {
4851  // otherwise, use the mapping
4852  // between shape function
4853  // numbers and rows. note that
4854  // by the assertions above, we
4855  // know that this particular
4856  // shape function is primitive,
4857  // so we can call
4858  // system_to_component_index
4859  const unsigned int row =
4860  this->finite_element_output
4861  .shape_function_to_row_table[i * fe->n_components() +
4862  fe->system_to_component_index(i).first];
4863  return this->finite_element_output.shape_values(row, j);
4864  }
4865 }
4866 
4867 
4868 
4869 template <int dim, int spacedim>
4870 inline double
4872  const unsigned int i,
4873  const unsigned int j,
4874  const unsigned int component) const
4875 {
4876  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4877  Assert(this->update_flags & update_values,
4878  ExcAccessToUninitializedField("update_values"));
4879  Assert(component < fe->n_components(),
4880  ExcIndexRange(component, 0, fe->n_components()));
4881  Assert(present_cell.get() != nullptr,
4882  ExcMessage("FEValues object is not reinit'ed to any cell"));
4883 
4884  // check whether the shape function
4885  // is non-zero at all within
4886  // this component:
4887  if (fe->get_nonzero_components(i)[component] == false)
4888  return 0;
4889 
4890  // look up the right row in the
4891  // table and take the data from
4892  // there
4893  const unsigned int row =
4894  this->finite_element_output
4895  .shape_function_to_row_table[i * fe->n_components() + component];
4896  return this->finite_element_output.shape_values(row, j);
4897 }
4898 
4899 
4900 
4901 template <int dim, int spacedim>
4902 inline const Tensor<1, spacedim> &
4903 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
4904  const unsigned int j) const
4905 {
4906  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4907  Assert(this->update_flags & update_gradients,
4908  ExcAccessToUninitializedField("update_gradients"));
4909  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4910  Assert(present_cell.get() != nullptr,
4911  ExcMessage("FEValues object is not reinit'ed to any cell"));
4912  // if the entire FE is primitive,
4913  // then we can take a short-cut:
4914  if (fe->is_primitive())
4915  return this->finite_element_output.shape_gradients[i][j];
4916  else
4917  {
4918  // otherwise, use the mapping
4919  // between shape function
4920  // numbers and rows. note that
4921  // by the assertions above, we
4922  // know that this particular
4923  // shape function is primitive,
4924  // so we can call
4925  // system_to_component_index
4926  const unsigned int row =
4927  this->finite_element_output
4928  .shape_function_to_row_table[i * fe->n_components() +
4929  fe->system_to_component_index(i).first];
4930  return this->finite_element_output.shape_gradients[row][j];
4931  }
4932 }
4933 
4934 
4935 
4936 template <int dim, int spacedim>
4937 inline Tensor<1, spacedim>
4939  const unsigned int i,
4940  const unsigned int j,
4941  const unsigned int component) const
4942 {
4943  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4944  Assert(this->update_flags & update_gradients,
4945  ExcAccessToUninitializedField("update_gradients"));
4946  Assert(component < fe->n_components(),
4947  ExcIndexRange(component, 0, fe->n_components()));
4948  Assert(present_cell.get() != nullptr,
4949  ExcMessage("FEValues object is not reinit'ed to any cell"));
4950  // check whether the shape function
4951  // is non-zero at all within
4952  // this component:
4953  if (fe->get_nonzero_components(i)[component] == false)
4954  return Tensor<1, spacedim>();
4955 
4956  // look up the right row in the
4957  // table and take the data from
4958  // there
4959  const unsigned int row =
4960  this->finite_element_output
4961  .shape_function_to_row_table[i * fe->n_components() + component];
4962  return this->finite_element_output.shape_gradients[row][j];
4963 }
4964 
4965 
4966 
4967 template <int dim, int spacedim>
4968 inline const Tensor<2, spacedim> &
4969 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
4970  const unsigned int j) const
4971 {
4972  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4973  Assert(this->update_flags & update_hessians,
4974  ExcAccessToUninitializedField("update_hessians"));
4975  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4976  Assert(present_cell.get() != nullptr,
4977  ExcMessage("FEValues object is not reinit'ed to any cell"));
4978  // if the entire FE is primitive,
4979  // then we can take a short-cut:
4980  if (fe->is_primitive())
4981  return this->finite_element_output.shape_hessians[i][j];
4982  else
4983  {
4984  // otherwise, use the mapping
4985  // between shape function
4986  // numbers and rows. note that
4987  // by the assertions above, we
4988  // know that this particular
4989  // shape function is primitive,
4990  // so we can call
4991  // system_to_component_index
4992  const unsigned int row =
4993  this->finite_element_output
4994  .shape_function_to_row_table[i * fe->n_components() +
4995  fe->system_to_component_index(i).first];
4996  return this->finite_element_output.shape_hessians[row][j];
4997  }
4998 }
4999 
5000 
5001 
5002 template <int dim, int spacedim>
5003 inline Tensor<2, spacedim>
5005  const unsigned int i,
5006  const unsigned int j,
5007  const unsigned int component) const
5008 {
5009  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5010  Assert(this->update_flags & update_hessians,
5011  ExcAccessToUninitializedField("update_hessians"));
5012  Assert(component < fe->n_components(),
5013  ExcIndexRange(component, 0, fe->n_components()));
5014  Assert(present_cell.get() != nullptr,
5015  ExcMessage("FEValues object is not reinit'ed to any cell"));
5016  // check whether the shape function
5017  // is non-zero at all within
5018  // this component:
5019  if (fe->get_nonzero_components(i)[component] == false)
5020  return Tensor<2, spacedim>();
5021 
5022  // look up the right row in the
5023  // table and take the data from
5024  // there
5025  const unsigned int row =
5026  this->finite_element_output
5027  .shape_function_to_row_table[i * fe->n_components() + component];
5028  return this->finite_element_output.shape_hessians[row][j];
5029 }
5030 
5031 
5032 
5033 template <int dim, int spacedim>
5034 inline const Tensor<3, spacedim> &
5036  const unsigned int j) const
5037 {
5038  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5039  Assert(this->update_flags & update_hessians,
5040  ExcAccessToUninitializedField("update_3rd_derivatives"));
5041  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5042  Assert(present_cell.get() != nullptr,
5043  ExcMessage("FEValues object is not reinit'ed to any cell"));
5044  // if the entire FE is primitive,
5045  // then we can take a short-cut:
5046  if (fe->is_primitive())
5047  return this->finite_element_output.shape_3rd_derivatives[i][j];
5048  else
5049  {
5050  // otherwise, use the mapping
5051  // between shape function
5052  // numbers and rows. note that
5053  // by the assertions above, we
5054  // know that this particular
5055  // shape function is primitive,
5056  // so we can call
5057  // system_to_component_index
5058  const unsigned int row =
5059  this->finite_element_output
5060  .shape_function_to_row_table[i * fe->n_components() +
5061  fe->system_to_component_index(i).first];
5062  return this->finite_element_output.shape_3rd_derivatives[row][j];
5063  }
5064 }
5065 
5066 
5067 
5068 template <int dim, int spacedim>
5069 inline Tensor<3, spacedim>
5071  const unsigned int i,
5072  const unsigned int j,
5073  const unsigned int component) const
5074 {
5075  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5076  Assert(this->update_flags & update_hessians,
5077  ExcAccessToUninitializedField("update_3rd_derivatives"));
5078  Assert(component < fe->n_components(),
5079  ExcIndexRange(component, 0, fe->n_components()));
5080  Assert(present_cell.get() != nullptr,
5081  ExcMessage("FEValues object is not reinit'ed to any cell"));
5082  // check whether the shape function
5083  // is non-zero at all within
5084  // this component:
5085  if (fe->get_nonzero_components(i)[component] == false)
5086  return Tensor<3, spacedim>();
5087 
5088  // look up the right row in the
5089  // table and take the data from
5090  // there
5091  const unsigned int row =
5092  this->finite_element_output
5093  .shape_function_to_row_table[i * fe->n_components() + component];
5094  return this->finite_element_output.shape_3rd_derivatives[row][j];
5095 }
5096 
5097 
5098 
5099 template <int dim, int spacedim>
5100 inline const FiniteElement<dim, spacedim> &
5102 {
5103  return *fe;
5104 }
5105 
5106 
5107 
5108 template <int dim, int spacedim>
5109 inline const Mapping<dim, spacedim> &
5111 {
5112  return *mapping;
5113 }
5114 
5115 
5116 
5117 template <int dim, int spacedim>
5118 inline UpdateFlags
5120 {
5121  return this->update_flags;
5122 }
5123 
5124 
5125 
5126 template <int dim, int spacedim>
5127 inline const std::vector<Point<spacedim>> &
5129 {
5130  Assert(this->update_flags & update_quadrature_points,
5131  ExcAccessToUninitializedField("update_quadrature_points"));
5132  Assert(present_cell.get() != nullptr,
5133  ExcMessage("FEValues object is not reinit'ed to any cell"));
5134  return this->mapping_output.quadrature_points;
5135 }
5136 
5137 
5138 
5139 template <int dim, int spacedim>
5140 inline const std::vector<double> &
5142 {
5143  Assert(this->update_flags & update_JxW_values,
5144  ExcAccessToUninitializedField("update_JxW_values"));
5145  Assert(present_cell.get() != nullptr,
5146  ExcMessage("FEValues object is not reinit'ed to any cell"));
5147  return this->mapping_output.JxW_values;
5148 }
5149 
5150 
5151 
5152 template <int dim, int spacedim>
5153 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5155 {
5156  Assert(this->update_flags & update_jacobians,
5157  ExcAccessToUninitializedField("update_jacobians"));
5158  Assert(present_cell.get() != nullptr,
5159  ExcMessage("FEValues object is not reinit'ed to any cell"));
5160  return this->mapping_output.jacobians;
5161 }
5162 
5163 
5164 
5165 template <int dim, int spacedim>
5166 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5168 {
5169  Assert(this->update_flags & update_jacobian_grads,
5170  ExcAccessToUninitializedField("update_jacobians_grads"));
5171  Assert(present_cell.get() != nullptr,
5172  ExcMessage("FEValues object is not reinit'ed to any cell"));
5173  return this->mapping_output.jacobian_grads;
5174 }
5175 
5176 
5177 
5178 template <int dim, int spacedim>
5179 inline const Tensor<3, spacedim> &
5181  const unsigned int i) const
5182 {
5183  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5184  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5185  Assert(present_cell.get() != nullptr,
5186  ExcMessage("FEValues object is not reinit'ed to any cell"));
5187  return this->mapping_output.jacobian_pushed_forward_grads[i];
5188 }
5189 
5190 
5191 
5192 template <int dim, int spacedim>
5193 inline const std::vector<Tensor<3, spacedim>> &
5195 {
5196  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5197  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5198  Assert(present_cell.get() != nullptr,
5199  ExcMessage("FEValues object is not reinit'ed to any cell"));
5200  return this->mapping_output.jacobian_pushed_forward_grads;
5201 }
5202 
5203 
5204 
5205 template <int dim, int spacedim>
5206 inline const DerivativeForm<3, dim, spacedim> &
5207 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5208 {
5209  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5210  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5211  Assert(present_cell.get() != nullptr,
5212  ExcMessage("FEValues object is not reinit'ed to any cell"));
5213  return this->mapping_output.jacobian_2nd_derivatives[i];
5214 }
5215 
5216 
5217 
5218 template <int dim, int spacedim>
5219 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5221 {
5222  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5223  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5224  Assert(present_cell.get() != nullptr,
5225  ExcMessage("FEValues object is not reinit'ed to any cell"));
5226  return this->mapping_output.jacobian_2nd_derivatives;
5227 }
5228 
5229 
5230 
5231 template <int dim, int spacedim>
5232 inline const Tensor<4, spacedim> &
5234  const unsigned int i) const
5235 {
5238  "update_jacobian_pushed_forward_2nd_derivatives"));
5239  Assert(present_cell.get() != nullptr,
5240  ExcMessage("FEValues object is not reinit'ed to any cell"));
5241  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5242 }
5243 
5244 
5245 
5246 template <int dim, int spacedim>
5247 inline const std::vector<Tensor<4, spacedim>> &
5249 {
5250  Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
5252  "update_jacobian_pushed_forward_2nd_derivatives"));
5253  Assert(present_cell.get() != nullptr,
5254  ExcMessage("FEValues object is not reinit'ed to any cell"));
5255  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5256 }
5257 
5258 
5259 
5260 template <int dim, int spacedim>
5261 inline const DerivativeForm<4, dim, spacedim> &
5262 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5263 {
5264  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5265  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5266  Assert(present_cell.get() != nullptr,
5267  ExcMessage("FEValues object is not reinit'ed to any cell"));
5268  return this->mapping_output.jacobian_3rd_derivatives[i];
5269 }
5270 
5271 
5272 
5273 template <int dim, int spacedim>
5274 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5276 {
5277  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5278  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5279  Assert(present_cell.get() != nullptr,
5280  ExcMessage("FEValues object is not reinit'ed to any cell"));
5281  return this->mapping_output.jacobian_3rd_derivatives;
5282 }
5283 
5284 
5285 
5286 template <int dim, int spacedim>
5287 inline const Tensor<5, spacedim> &
5289  const unsigned int i) const
5290 {
5293  "update_jacobian_pushed_forward_3rd_derivatives"));
5294  Assert(present_cell.get() != nullptr,
5295  ExcMessage("FEValues object is not reinit'ed to any cell"));
5296  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
5297 }
5298 
5299 
5300 
5301 template <int dim, int spacedim>
5302 inline const std::vector<Tensor<5, spacedim>> &
5304 {
5305  Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
5307  "update_jacobian_pushed_forward_3rd_derivatives"));
5308  Assert(present_cell.get() != nullptr,
5309  ExcMessage("FEValues object is not reinit'ed to any cell"));
5310  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
5311 }
5312 
5313 
5314 
5315 template <int dim, int spacedim>
5316 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
5318 {
5319  Assert(this->update_flags & update_inverse_jacobians,
5320  ExcAccessToUninitializedField("update_inverse_jacobians"));
5321  Assert(present_cell.get() != nullptr,
5322  ExcMessage("FEValues object is not reinit'ed to any cell"));
5323  return this->mapping_output.inverse_jacobians;
5324 }
5325 
5326 
5327 
5328 template <int dim, int spacedim>
5329 inline const Point<spacedim> &
5330 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
5331 {
5332  Assert(this->update_flags & update_quadrature_points,
5333  ExcAccessToUninitializedField("update_quadrature_points"));
5334  Assert(i < this->mapping_output.quadrature_points.size(),
5335  ExcIndexRange(i, 0, this->mapping_output.quadrature_points.size()));
5336  Assert(present_cell.get() != nullptr,
5337  ExcMessage("FEValues object is not reinit'ed to any cell"));
5338 
5339  return this->mapping_output.quadrature_points[i];
5340 }
5341 
5342 
5343 
5344 template <int dim, int spacedim>
5345 inline double
5346 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
5347 {
5348  Assert(this->update_flags & update_JxW_values,
5349  ExcAccessToUninitializedField("update_JxW_values"));
5350  Assert(i < this->mapping_output.JxW_values.size(),
5351  ExcIndexRange(i, 0, this->mapping_output.JxW_values.size()));
5352  Assert(present_cell.get() != nullptr,
5353  ExcMessage("FEValues object is not reinit'ed to any cell"));
5354 
5355  return this->mapping_output.JxW_values[i];
5356 }
5357 
5358 
5359 
5360 template <int dim, int spacedim>
5361 inline const DerivativeForm<1, dim, spacedim> &
5362 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
5363 {
5364  Assert(this->update_flags & update_jacobians,
5365  ExcAccessToUninitializedField("update_jacobians"));
5366  Assert(i < this->mapping_output.jacobians.size(),
5367  ExcIndexRange(i, 0, this->mapping_output.jacobians.size()));
5368  Assert(present_cell.get() != nullptr,
5369  ExcMessage("FEValues object is not reinit'ed to any cell"));
5370 
5371  return this->mapping_output.jacobians[i];
5372 }
5373 
5374 
5375 
5376 template <int dim, int spacedim>
5377 inline const DerivativeForm<2, dim, spacedim> &
5378 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
5379 {
5380  Assert(this->update_flags & update_jacobian_grads,
5381  ExcAccessToUninitializedField("update_jacobians_grads"));
5382  Assert(i < this->mapping_output.jacobian_grads.size(),
5383  ExcIndexRange(i, 0, this->mapping_output.jacobian_grads.size()));
5384  Assert(present_cell.get() != nullptr,
5385  ExcMessage("FEValues object is not reinit'ed to any cell"));
5386 
5387  return this->mapping_output.jacobian_grads[i];
5388 }
5389 
5390 
5391 
5392 template <int dim, int spacedim>
5393 inline const DerivativeForm<1, spacedim, dim> &
5394 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
5395 {
5396  Assert(this->update_flags & update_inverse_jacobians,
5397  ExcAccessToUninitializedField("update_inverse_jacobians"));
5398  Assert(i < this->mapping_output.inverse_jacobians.size(),
5399  ExcIndexRange(i, 0, this->mapping_output.inverse_jacobians.size()));
5400  Assert(present_cell.get() != nullptr,
5401  ExcMessage("FEValues object is not reinit'ed to any cell"));
5402 
5403  return this->mapping_output.inverse_jacobians[i];
5404 }
5405 
5406 
5407 
5408 template <int dim, int spacedim>
5409 inline const Tensor<1, spacedim> &
5410 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
5411 {
5412  Assert(this->update_flags & update_normal_vectors,
5414  "update_normal_vectors")));
5415  Assert(i < this->mapping_output.normal_vectors.size(),
5416  ExcIndexRange(i, 0, this->mapping_output.normal_vectors.size()));
5417  Assert(present_cell.get() != nullptr,
5418  ExcMessage("FEValues object is not reinit'ed to any cell"));
5419 
5420  return this->mapping_output.normal_vectors[i];
5421 }
5422 
5423 
5424 
5425 /*--------------------- Inline functions: FEValues --------------------------*/
5426 
5427 
5428 template <int dim, int spacedim>
5429 inline const Quadrature<dim> &
5431 {
5432  return quadrature;
5433 }
5434 
5435 
5436 
5437 template <int dim, int spacedim>
5438 inline const FEValues<dim, spacedim> &
5440 {
5441  return *this;
5442 }
5443 
5444 
5445 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
5446 
5447 
5448 template <int dim, int spacedim>
5449 inline unsigned int
5451 {
5452  return present_face_index;
5453 }
5454 
5455 
5456 /*----------------------- Inline functions: FE*FaceValues -------------------*/
5457 
5458 template <int dim, int spacedim>
5459 inline const Quadrature<dim - 1> &
5461 {
5462  return quadrature;
5463 }
5464 
5465 
5466 
5467 template <int dim, int spacedim>
5468 inline const FEFaceValues<dim, spacedim> &
5470 {
5471  return *this;
5472 }
5473 
5474 
5475 
5476 template <int dim, int spacedim>
5477 inline const FESubfaceValues<dim, spacedim> &
5479 {
5480  return *this;
5481 }
5482 
5483 
5484 
5485 template <int dim, int spacedim>
5486 inline const Tensor<1, spacedim> &
5487 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
5488 {
5489  Assert(i < this->mapping_output.boundary_forms.size(),
5490  ExcIndexRange(i, 0, this->mapping_output.boundary_forms.size()));
5491  Assert(this->update_flags & update_boundary_forms,
5493  "update_boundary_forms")));
5494 
5495  return this->mapping_output.boundary_forms[i];
5496 }
5497 
5498 #endif // DOXYGEN
5499 
5500 DEAL_II_NAMESPACE_CLOSE
5501 
5502 #endif
ProductType< Number, typename Scalar< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:215
Transformed quadrature weights.
const std::vector< Tensor< 4, spacedim > > & get_jacobian_pushed_forward_2nd_derivatives() const
void get_function_hessians_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::hessian_type > &hessians) const
Definition: fe_values.cc:2210
~FEValuesBase() override
Definition: fe_values.cc:3184
void get_function_values(const InputVector &fe_function, std::vector< typename InputVector::value_type > &values) const
Definition: fe_values.cc:3657
const Tensor< 3, spacedim > & shape_3rd_derivative(const unsigned int function_no, const unsigned int point_no) const
void get_function_curls(const InputVector &fe_function, std::vector< typename ProductType< curl_type, typename InputVector::value_type >::type > &curls) const
Definition: fe_values.cc:2123
Shape function values.
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int quadrature_point) const
const Tensor< 1, spacedim > & shape_grad(const unsigned int function_no, const unsigned int quadrature_point) const
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int quadrature_point) const
void get_function_laplacians(const InputVector &fe_function, std::vector< typename ProductType< value_type, typename InputVector::value_type >::type > &laplacians) const
Definition: fe_values.cc:1788
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access >> &cell, const unsigned int face_no, const unsigned int subface_no)
Definition: fe_values.cc:5023
void get_function_values(const InputVector &fe_function, std::vector< typename ProductType< value_type, typename InputVector::value_type >::type > &values) const
Definition: fe_values.cc:1619
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > fe_data
Definition: fe_values.h:3355
Number value_type
Definition: vector.h:115
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:3387
::Tensor< 1, spacedim > divergence_type
Definition: fe_values.h:1568
static const unsigned int dimension
Definition: fe_values.h:3781
const FEFaceValues< dim, spacedim > & get_present_fe_values() const
Cache(const FEValuesBase< dim, spacedim > &fe_values)
Definition: fe_values.cc:2641
std::size_t memory_consumption() const
Definition: fe_values.cc:4742
unsigned int present_face_index
Definition: fe_values.h:3635
ProductType< Number, typename Vector< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:693
void get_function_hessians(const InputVector &fe_function, std::vector< typename ProductType< hessian_type, typename InputVector::value_type >::type > &hessians) const
Definition: fe_values.cc:1732
void get_function_third_derivatives(const InputVector &fe_function, std::vector< Tensor< 3, spacedim, typename InputVector::value_type >> &third_derivatives) const
Definition: fe_values.cc:4164
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1528
void get_function_third_derivatives_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::third_derivative_type > &third_derivatives) const
Definition: fe_values.cc:1875
static::ExceptionBase & ExcAccessToUninitializedField()
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:538
curl_type curl(const unsigned int shape_function, const unsigned int q_point) const
void get_function_symmetric_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::symmetric_gradient_type > &symmetric_gradients) const
Definition: fe_values.cc:2042
static::ExceptionBase & ExcFaceHasNoSubfaces()
const DerivativeForm< 3, dim, spacedim > & jacobian_2nd_derivative(const unsigned int quadrature_point) const
const unsigned int dofs_per_cell
Definition: fe_values.h:2042
const unsigned int component
Definition: fe_values.h:544
::Tensor< 3, spacedim > hessian_type
Definition: fe_values.h:636
const Quadrature< dim-1 > quadrature
Definition: fe_values.h:3640
::SymmetricTensor< 2, spacedim > symmetric_gradient_type
Definition: fe_values.h:614
FEValuesBase(const unsigned int n_q_points, const unsigned int dofs_per_cell, const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe)
Definition: fe_values.cc:3161
ProductType< Number, typename Vector< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:659
const Tensor< 1, spacedim > & normal_vector(const unsigned int i) const
void get_function_gradients(const InputVector &fe_function, std::vector< Tensor< 1, spacedim, typename InputVector::value_type >> &gradients) const
Definition: fe_values.cc:3799
Volume element.
const Tensor< 3, spacedim > & jacobian_pushed_forward_grad(const unsigned int quadrature_point) const
static const unsigned int dimension
Definition: fe_values.h:3667
unsigned int get_face_index() const
static::ExceptionBase & ExcAccessToUninitializedField(std::string arg1)
ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1302
::internal::CurlType< spacedim >::type curl_type
Definition: fe_values.h:629
static::ExceptionBase & ExcReinitCalledWithBoundaryFace()
Outer normal vector, not normalized.
std::unique_ptr< const CellIteratorBase > present_cell
Definition: fe_values.h:3271
third_derivative_type third_derivative(const unsigned int shape_function, const unsigned int q_point) const
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
static::ExceptionBase & ExcFEDontMatch()
Scalar & operator=(const Scalar< dim, spacedim > &)
Definition: fe_values.cc:188
ProductType< Number, typename Tensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1598
Transformed quadrature points.
void get_function_curls_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::curl_type > &curls) const
Definition: fe_values.cc:2154
void do_reinit(const unsigned int face_no)
Definition: fe_values.cc:4897
const Triangulation< dim, spacedim >::cell_iterator get_cell() const
Definition: fe_values.cc:4279
ProductType< Number, typename Vector< dim, spacedim >::symmetric_gradient_type >::type symmetric_gradient_type
Definition: fe_values.h:676
const std::vector< DerivativeForm< 4, dim, spacedim > > & get_jacobian_3rd_derivatives() const
value_type value(const unsigned int shape_function, const unsigned int q_point) const
ProductType< Number, typename Vector< dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:684
::SymmetricTensor< 2, spacedim > value_type
Definition: fe_values.h:1275
unsigned int row_index[spacedim]
Definition: fe_values.h:748
void get_function_values_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::value_type > &values) const
Definition: fe_values.cc:1930
void check_cell_similarity(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
Definition: fe_values.cc:4417
static const unsigned int space_dimension
Definition: fe_values.h:2030
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:3323
::Tensor< 1, spacedim > value_type
Definition: fe_values.h:592
::Tensor< 3, spacedim > third_derivative_type
Definition: fe_values.h:173
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
void get_function_laplacians_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::laplacian_type > &laplacians) const
Definition: fe_values.cc:1819
void get_function_third_derivatives(const InputVector &fe_function, std::vector< typename ProductType< third_derivative_type, typename InputVector::value_type >::type > &third_derivatives) const
Definition: fe_values.cc:1844
static::ExceptionBase & ExcShapeFunctionNotPrimitive(int arg1)
void get_function_symmetric_gradients(const InputVector &fe_function, std::vector< typename ProductType< symmetric_gradient_type, typename InputVector::value_type >::type > &symmetric_gradients) const
Definition: fe_values.cc:2011
void get_function_gradients(const InputVector &fe_function, std::vector< typename ProductType< gradient_type, typename InputVector::value_type >::type > &gradients) const
Definition: fe_values.cc:1955
const std::vector< double > & get_JxW_values() const
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition: fe_values.h:3415
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
void do_reinit(const unsigned int face_no, const unsigned int subface_no)
Definition: fe_values.cc:5096
UpdateFlags get_update_flags() const
static const unsigned int integral_dimension
Definition: fe_values.h:3792
void get_function_laplacians_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::laplacian_type > &laplacians) const
Definition: fe_values.cc:2271
const Quadrature< dim > & get_quadrature() const
UpdateFlags compute_update_flags(const UpdateFlags update_flags) const
Definition: fe_values.cc:4333
::Tensor< 2, spacedim > hessian_type
Definition: fe_values.h:166
static const unsigned int integral_dimension
Definition: fe_values.h:3675
void get_function_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::gradient_type > &gradients) const
Definition: fe_values.cc:1986
const std::vector< DerivativeForm< 3, dim, spacedim > > & get_jacobian_2nd_derivatives() const
void get_function_laplacians(const InputVector &fe_function, std::vector< typename InputVector::value_type > &laplacians) const
Definition: fe_values.cc:4027
static::ExceptionBase & ExcMessage(std::string arg1)
const std::vector< Point< spacedim > > & get_quadrature_points() const
void get_function_hessians(const InputVector &fe_function, std::vector< typename ProductType< hessian_type, typename InputVector::value_type >::type > &hessians) const
Definition: fe_values.cc:2179
symmetric_gradient_type symmetric_gradient(const unsigned int shape_function, const unsigned int q_point) const
const std::vector< Tensor< 5, spacedim > > & get_jacobian_pushed_forward_3rd_derivatives() const
const Quadrature< dim-1 > & get_quadrature() const
divergence_type divergence(const unsigned int shape_function, const unsigned int q_point) const
const std::vector< Tensor< 1, spacedim > > & get_boundary_forms() const
Definition: fe_values.cc:4730
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:408
third_derivative_type third_derivative(const unsigned int shape_function, const unsigned int q_point) const
const Tensor< 2, spacedim > & shape_hessian(const unsigned int function_no, const unsigned int point_no) const
Third derivatives of shape functions.
void get_function_gradients(const InputVector &fe_function, std::vector< typename ProductType< gradient_type, typename InputVector::value_type >::type > &gradients) const
Definition: fe_values.cc:1676
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1875
void get_function_hessians_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::hessian_type > &hessians) const
Definition: fe_values.cc:1763
#define Assert(cond, exc)
Definition: exceptions.h:1227
::Tensor< 2, spacedim > value_type
Definition: fe_values.h:1563
UpdateFlags
void get_function_third_derivatives(const InputVector &fe_function, std::vector< typename ProductType< third_derivative_type, typename InputVector::value_type >::type > &third_derivatives) const
Definition: fe_values.cc:2299
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
ProductType< Number, typename Scalar< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:224
Tensor< 3, spacedim > shape_3rd_derivative_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
Abstract base class for mapping classes.
Definition: dof_tools.h:57
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1234
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:397
const Quadrature< dim > quadrature
Definition: fe_values.h:3533
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
const unsigned int first_vector_component
Definition: fe_values.h:1229
std::size_t memory_consumption() const
Definition: fe_values.cc:4313
#define DeclException0(Exception0)
Definition: exceptions.h:385
const FiniteElement< dim, spacedim > & get_fe() const
void get_function_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::gradient_type > &gradients) const
Definition: fe_values.cc:1707
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
void invalidate_present_cell()
Definition: fe_values.cc:4352
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition: fe_values.h:3331
static const unsigned int integral_dimension
Definition: fe_values.h:3451
Vector & operator=(const Vector< dim, spacedim > &)
Definition: fe_values.cc:280
SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
Tensor()
Definition: tensor.h:1017
ProductType< Number, typename Scalar< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:197
const std::vector< Tensor< 1, spacedim > > & get_normal_vectors() const
Definition: fe_values.cc:4300
void get_function_divergences_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::divergence_type > &divergences) const
Definition: fe_values.cc:2098
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
bool is_nonzero_shape_function_component[spacedim]
Definition: fe_values.h:737
const FEValues< dim, spacedim > & get_present_fe_values() const
static const unsigned int integral_dimension
Definition: fe_values.h:3570
Second derivatives of shape functions.
static::ExceptionBase & ExcFENotPrimitive()
Gradient of volume element.
::Tensor< 2, spacedim > gradient_type
Definition: fe_values.h:602
FEValuesBase & operator=(const FEValuesBase &)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
Definition: tensor.h:1364
std::vector<::FEValuesViews::Scalar< dim, spacedim > > scalars
Definition: fe_values.h:1899
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access >> &cell)
Definition: fe_values.cc:4638
void get_function_values(const InputVector &fe_function, std::vector< typename ProductType< value_type, typename InputVector::value_type >::type > &values) const
Definition: fe_values.cc:1899
const unsigned int n_quadrature_points
Definition: fe_values.h:2035
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
::Tensor< 1, spacedim > gradient_type
Definition: fe_values.h:159
const FESubfaceValues< dim, spacedim > & get_present_fe_values() const
static const unsigned int dimension
Definition: fe_values.h:2025
boost::signals2::connection tria_listener_mesh_transform
Definition: fe_values.h:3296
ProductType< Number, typename Vector< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:667
ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1311
static const unsigned int space_dimension
Definition: fe_values.h:3786
const Tensor< 5, spacedim > & jacobian_pushed_forward_3rd_derivative(const unsigned int quadrature_point) const
const Tensor< 4, spacedim > & jacobian_pushed_forward_2nd_derivative(const unsigned int quadrature_point) const
Definition: mpi.h:55
ProductType< Number, typename Tensor< 2, dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:1606
FEFaceValuesBase(const unsigned int n_q_points, const unsigned int dofs_per_cell, const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim-1 > &quadrature)
Definition: fe_values.cc:4710
void initialize(const UpdateFlags update_flags)
Definition: fe_values.cc:4529
ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:189
Shape function gradients.
Normal vectors.
void maybe_invalidate_previous_present_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
Definition: fe_values.cc:4370
void get_function_divergences(const InputVector &fe_function, std::vector< typename ProductType< divergence_type, typename InputVector::value_type >::type > &divergences) const
Definition: fe_values.cc:2066
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1864
Definition: fe.h:36
void initialize(const UpdateFlags update_flags)
Definition: fe_values.cc:4974
const std::vector< DerivativeForm< 1, spacedim, dim > > & get_inverse_jacobians() const
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1517
const std::vector< Tensor< 1, spacedim > > & get_all_normal_vectors() const
Definition: fe_values.cc:4288
value_type value(const unsigned int shape_function, const unsigned int q_point) const
void get_function_laplacians(const InputVector &fe_function, std::vector< typename ProductType< value_type, typename InputVector::value_type >::type > &laplacians) const
Definition: fe_values.cc:2235
static::ExceptionBase & ExcNotImplemented()
ProductType< Number, typename Vector< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:720
ProductType< Number, typename Vector< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:711
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access >> &cell, const unsigned int face_no)
Definition: fe_values.cc:4845
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1223
boost::signals2::connection tria_listener_refinement
Definition: fe_values.h:3287
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3363
::Tensor< 3, spacedim > gradient_type
Definition: fe_values.h:1574
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int quadrature_point) const
const std::vector< DerivativeForm< 2, dim, spacedim > > & get_jacobian_grads() const
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)
const Point< spacedim > & quadrature_point(const unsigned int q) const
void initialize(const UpdateFlags update_flags)
Definition: fe_values.cc:4798
::Tensor< 4, spacedim > third_derivative_type
Definition: fe_values.h:643
::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:3338
FEValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4495
const Mapping< dim, spacedim > & get_mapping() const
void do_reinit()
Definition: fe_values.cc:4665
ProductType< Number, typename Vector< dim, spacedim >::curl_type >::type curl_type
Definition: fe_values.h:702
FESubfaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim-1 > &face_quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4938
Tensor< 1, spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
double JxW(const unsigned int quadrature_point) const
void get_function_values_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::value_type > &values) const
Definition: fe_values.cc:1651
ProductType< Number, typename Tensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1589
FEFaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim-1 > &quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4762
ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:206
CellSimilarity::Similarity get_cell_similarity() const
Definition: fe_values.cc:4472
void get_function_third_derivatives_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::third_derivative_type > &third_derivatives) const
Definition: fe_values.cc:2330
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:549
const DerivativeForm< 2, dim, spacedim > & jacobian_grad(const unsigned int quadrature_point) const
UpdateFlags update_flags
Definition: fe_values.h:3369
std::size_t memory_consumption() const
Definition: fe_values.cc:4699
static::ExceptionBase & ExcInternalError()
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
Definition: fe_values.h:3347
const Tensor< 1, spacedim > & boundary_form(const unsigned int i) const
Tensor< 2, spacedim > shape_hessian_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
void get_function_hessians(const InputVector &fe_function, std::vector< Tensor< 2, spacedim, typename InputVector::value_type >> &hessians) const
Definition: fe_values.cc:3913