Reference documentation for deal.II version Git e7c9a24 2019-02-08 08:26:19 +0100
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
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/vector_slice.h>
29 
30 #include <deal.II/dofs/dof_accessor.h>
31 #include <deal.II/dofs/dof_handler.h>
32 
33 #include <deal.II/fe/fe.h>
34 #include <deal.II/fe/fe_update_flags.h>
35 #include <deal.II/fe/fe_values_extractors.h>
36 #include <deal.II/fe/mapping.h>
37 
38 #include <deal.II/grid/tria.h>
39 #include <deal.II/grid/tria_iterator.h>
40 
41 #include <deal.II/hp/dof_handler.h>
42 
43 #include <algorithm>
44 #include <memory>
45 #include <type_traits>
46 
47 
48 // dummy include in order to have the
49 // definition of PetscScalar available
50 // without including other PETSc stuff
51 #ifdef DEAL_II_WITH_PETSC
52 # include <petsc.h>
53 #endif
54 
55 DEAL_II_NAMESPACE_OPEN
56 
57 template <int dim, int spacedim = dim>
58 class FEValuesBase;
59 
60 namespace internal
61 {
66  template <int dim, class NumberType = double>
67  struct CurlType;
68 
75  template <class NumberType>
76  struct CurlType<1, NumberType>
77  {
79  };
80 
87  template <class NumberType>
88  struct CurlType<2, NumberType>
89  {
91  };
92 
99  template <class NumberType>
100  struct CurlType<3, NumberType>
101  {
103  };
104 } // namespace internal
105 
106 
107 
129 namespace FEValuesViews
130 {
142  template <int dim, int spacedim = dim>
143  class Scalar
144  {
145  public:
151  using value_type = double;
152 
159 
166 
173 
178  template <typename Number>
179  struct OutputType
180  {
185  using value_type =
186  typename ProductType<Number,
188 
193  using gradient_type = typename ProductType<
194  Number,
196 
201  using laplacian_type =
202  typename ProductType<Number,
204 
209  using hessian_type = typename ProductType<
210  Number,
212 
217  using third_derivative_type = typename ProductType<
218  Number,
220  };
221 
227  {
237 
246  unsigned int row_index;
247  };
248 
252  Scalar();
253 
259  Scalar(const FEValuesBase<dim, spacedim> &fe_values_base,
260  const unsigned int component);
261 
266  Scalar &
267  operator=(const Scalar<dim, spacedim> &) = delete;
268 
282  value_type
283  value(const unsigned int shape_function, const unsigned int q_point) const;
284 
296  gradient(const unsigned int shape_function,
297  const unsigned int q_point) const;
298 
310  hessian(const unsigned int shape_function,
311  const unsigned int q_point) const;
312 
324  third_derivative(const unsigned int shape_function,
325  const unsigned int q_point) const;
326 
344  template <class InputVector>
345  void
347  const InputVector &fe_function,
348  std::vector<typename ProductType<value_type,
349  typename InputVector::value_type>::type>
350  &values) const;
351 
374  template <class InputVector>
375  void
377  const InputVector &dof_values,
378  std::vector<
380  &values) const;
381 
399  template <class InputVector>
400  void
402  const InputVector &fe_function,
403  std::vector<typename ProductType<gradient_type,
404  typename InputVector::value_type>::type>
405  &gradients) const;
406 
410  template <class InputVector>
411  void
413  const InputVector &dof_values,
414  std::vector<
416  &gradients) const;
417 
435  template <class InputVector>
436  void
438  const InputVector &fe_function,
439  std::vector<typename ProductType<hessian_type,
440  typename InputVector::value_type>::type>
441  &hessians) const;
442 
446  template <class InputVector>
447  void
449  const InputVector &dof_values,
450  std::vector<
452  &hessians) const;
453 
454 
473  template <class InputVector>
474  void
476  const InputVector &fe_function,
477  std::vector<typename ProductType<value_type,
478  typename InputVector::value_type>::type>
479  &laplacians) const;
480 
484  template <class InputVector>
485  void
487  const InputVector &dof_values,
488  std::vector<
490  &laplacians) const;
491 
492 
511  template <class InputVector>
512  void
514  const InputVector &fe_function,
515  std::vector<typename ProductType<third_derivative_type,
516  typename InputVector::value_type>::type>
517  &third_derivatives) const;
518 
522  template <class InputVector>
523  void
525  const InputVector & dof_values,
526  std::vector<typename OutputType<typename InputVector::value_type>::
527  third_derivative_type> &third_derivatives) const;
528 
529 
530  private:
535 
540  const unsigned int component;
541 
545  std::vector<ShapeFunctionData> shape_function_data;
546  };
547 
548 
549 
579  template <int dim, int spacedim = dim>
580  class Vector
581  {
582  public:
589 
599 
611 
617  using divergence_type = double;
618 
625  using curl_type = typename ::internal::CurlType<spacedim>::type;
626 
633 
640 
645  template <typename Number>
646  struct OutputType
647  {
652  using value_type =
653  typename ProductType<Number,
655 
660  using gradient_type = typename ProductType<
661  Number,
663 
668  using symmetric_gradient_type = typename ProductType<
669  Number,
671 
676  using divergence_type = typename ProductType<
677  Number,
679 
684  using laplacian_type =
685  typename ProductType<Number,
687 
692  using curl_type =
693  typename ProductType<Number,
695 
700  using hessian_type = typename ProductType<
701  Number,
703 
708  using third_derivative_type = typename ProductType<
709  Number,
711  };
712 
718  {
728 
738  unsigned int row_index[spacedim];
739 
749  unsigned int single_nonzero_component_index;
750  };
751 
755  Vector();
756 
765  Vector(const FEValuesBase<dim, spacedim> &fe_values_base,
766  const unsigned int first_vector_component);
767 
772  Vector &
773  operator=(const Vector<dim, spacedim> &) = delete;
774 
791  value_type
792  value(const unsigned int shape_function, const unsigned int q_point) const;
793 
808  gradient(const unsigned int shape_function,
809  const unsigned int q_point) const;
810 
827  symmetric_gradient(const unsigned int shape_function,
828  const unsigned int q_point) const;
829 
841  divergence(const unsigned int shape_function,
842  const unsigned int q_point) const;
843 
864  curl_type
865  curl(const unsigned int shape_function, const unsigned int q_point) const;
866 
878  hessian(const unsigned int shape_function,
879  const unsigned int q_point) const;
880 
892  third_derivative(const unsigned int shape_function,
893  const unsigned int q_point) const;
894 
912  template <class InputVector>
913  void
915  const InputVector &fe_function,
916  std::vector<typename ProductType<value_type,
917  typename InputVector::value_type>::type>
918  &values) const;
919 
942  template <class InputVector>
943  void
945  const InputVector &dof_values,
946  std::vector<
948  &values) const;
949 
967  template <class InputVector>
968  void
970  const InputVector &fe_function,
971  std::vector<typename ProductType<gradient_type,
972  typename InputVector::value_type>::type>
973  &gradients) const;
974 
978  template <class InputVector>
979  void
981  const InputVector &dof_values,
982  std::vector<
984  &gradients) const;
985 
1009  template <class InputVector>
1010  void
1012  const InputVector &fe_function,
1013  std::vector<typename ProductType<symmetric_gradient_type,
1014  typename InputVector::value_type>::type>
1015  &symmetric_gradients) const;
1016 
1020  template <class InputVector>
1021  void
1023  const InputVector & dof_values,
1024  std::vector<typename OutputType<typename InputVector::value_type>::
1025  symmetric_gradient_type> &symmetric_gradients) const;
1026 
1045  template <class InputVector>
1046  void
1048  const InputVector &fe_function,
1049  std::vector<typename ProductType<divergence_type,
1050  typename InputVector::value_type>::type>
1051  &divergences) const;
1052 
1056  template <class InputVector>
1057  void
1059  const InputVector &dof_values,
1060  std::vector<
1062  &divergences) const;
1063 
1082  template <class InputVector>
1083  void
1085  const InputVector &fe_function,
1086  std::vector<
1087  typename ProductType<curl_type, typename InputVector::value_type>::type>
1088  &curls) const;
1089 
1093  template <class InputVector>
1094  void
1096  const InputVector &dof_values,
1097  std::vector<
1099  &curls) const;
1100 
1118  template <class InputVector>
1119  void
1121  const InputVector &fe_function,
1122  std::vector<typename ProductType<hessian_type,
1123  typename InputVector::value_type>::type>
1124  &hessians) const;
1125 
1129  template <class InputVector>
1130  void
1132  const InputVector &dof_values,
1133  std::vector<
1135  &hessians) const;
1136 
1155  template <class InputVector>
1156  void
1158  const InputVector &fe_function,
1159  std::vector<typename ProductType<value_type,
1160  typename InputVector::value_type>::type>
1161  &laplacians) const;
1162 
1166  template <class InputVector>
1167  void
1169  const InputVector &dof_values,
1170  std::vector<
1172  &laplacians) const;
1173 
1192  template <class InputVector>
1193  void
1195  const InputVector &fe_function,
1196  std::vector<typename ProductType<third_derivative_type,
1197  typename InputVector::value_type>::type>
1198  &third_derivatives) const;
1199 
1203  template <class InputVector>
1204  void
1206  const InputVector & dof_values,
1207  std::vector<typename OutputType<typename InputVector::value_type>::
1208  third_derivative_type> &third_derivatives) const;
1209 
1210  private:
1215 
1220  const unsigned int first_vector_component;
1221 
1225  std::vector<ShapeFunctionData> shape_function_data;
1226  };
1227 
1228 
1229  template <int rank, int dim, int spacedim = dim>
1230  class SymmetricTensor;
1231 
1256  template <int dim, int spacedim>
1257  class SymmetricTensor<2, dim, spacedim>
1258  {
1259  public:
1267 
1278 
1283  template <typename Number>
1284  struct OutputType
1285  {
1290  using value_type = typename ProductType<
1291  Number,
1293 
1298  using divergence_type = typename ProductType<
1299  Number,
1301  };
1302 
1307  struct ShapeFunctionData
1308  {
1317  bool is_nonzero_shape_function_component
1318  [value_type::n_independent_components];
1319 
1329  unsigned int row_index[value_type::n_independent_components];
1330 
1340 
1345  };
1346 
1350  SymmetricTensor();
1351 
1361  SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1362  const unsigned int first_tensor_component);
1363 
1368  SymmetricTensor &
1369  operator=(const SymmetricTensor<2, dim, spacedim> &) = delete;
1370 
1388  value_type
1389  value(const unsigned int shape_function, const unsigned int q_point) const;
1390 
1405  divergence(const unsigned int shape_function,
1406  const unsigned int q_point) const;
1407 
1425  template <class InputVector>
1426  void
1427  get_function_values(
1428  const InputVector &fe_function,
1429  std::vector<typename ProductType<value_type,
1430  typename InputVector::value_type>::type>
1431  &values) const;
1432 
1455  template <class InputVector>
1456  void
1457  get_function_values_from_local_dof_values(
1458  const InputVector &dof_values,
1459  std::vector<
1461  &values) const;
1462 
1484  template <class InputVector>
1485  void
1486  get_function_divergences(
1487  const InputVector &fe_function,
1488  std::vector<typename ProductType<divergence_type,
1489  typename InputVector::value_type>::type>
1490  &divergences) const;
1491 
1495  template <class InputVector>
1496  void
1497  get_function_divergences_from_local_dof_values(
1498  const InputVector &dof_values,
1499  std::vector<
1501  &divergences) const;
1502 
1503  private:
1508 
1513  const unsigned int first_tensor_component;
1514 
1518  std::vector<ShapeFunctionData> shape_function_data;
1519  };
1520 
1521 
1522  template <int rank, int dim, int spacedim = dim>
1523  class Tensor;
1524 
1545  template <int dim, int spacedim>
1546  class Tensor<2, dim, spacedim>
1547  {
1548  public:
1554 
1559 
1565 
1570  template <typename Number>
1571  struct OutputType
1572  {
1577  using value_type = typename ProductType<
1578  Number,
1580 
1585  using divergence_type = typename ProductType<
1586  Number,
1588 
1593  using gradient_type = typename ProductType<
1594  Number,
1596  };
1597 
1602  struct ShapeFunctionData
1603  {
1612  bool is_nonzero_shape_function_component
1613  [value_type::n_independent_components];
1614 
1624  unsigned int row_index[value_type::n_independent_components];
1625 
1635 
1640  };
1641 
1645  Tensor();
1646 
1656  Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1657  const unsigned int first_tensor_component);
1658 
1659 
1664  Tensor &
1665  operator=(const Tensor<2, dim, spacedim> &) = delete;
1666 
1683  value_type
1684  value(const unsigned int shape_function, const unsigned int q_point) const;
1685 
1700  divergence(const unsigned int shape_function,
1701  const unsigned int q_point) const;
1702 
1717  gradient(const unsigned int shape_function,
1718  const unsigned int q_point) const;
1719 
1737  template <class InputVector>
1738  void
1739  get_function_values(
1740  const InputVector &fe_function,
1741  std::vector<typename ProductType<value_type,
1742  typename InputVector::value_type>::type>
1743  &values) const;
1744 
1767  template <class InputVector>
1768  void
1769  get_function_values_from_local_dof_values(
1770  const InputVector &dof_values,
1771  std::vector<
1773  &values) const;
1774 
1796  template <class InputVector>
1797  void
1798  get_function_divergences(
1799  const InputVector &fe_function,
1800  std::vector<typename ProductType<divergence_type,
1801  typename InputVector::value_type>::type>
1802  &divergences) const;
1803 
1807  template <class InputVector>
1808  void
1809  get_function_divergences_from_local_dof_values(
1810  const InputVector &dof_values,
1811  std::vector<
1813  &divergences) const;
1814 
1831  template <class InputVector>
1832  void
1833  get_function_gradients(
1834  const InputVector &fe_function,
1835  std::vector<typename ProductType<gradient_type,
1836  typename InputVector::value_type>::type>
1837  &gradients) const;
1838 
1842  template <class InputVector>
1843  void
1844  get_function_gradients_from_local_dof_values(
1845  const InputVector &dof_values,
1846  std::vector<
1848  &gradients) const;
1849 
1850  private:
1855 
1860  const unsigned int first_tensor_component;
1861 
1865  std::vector<ShapeFunctionData> shape_function_data;
1866  };
1867 
1868 } // namespace FEValuesViews
1869 
1870 
1871 namespace internal
1872 {
1873  namespace FEValuesViews
1874  {
1882  template <int dim, int spacedim>
1883  struct Cache
1884  {
1889  std::vector<::FEValuesViews::Scalar<dim, spacedim>> scalars;
1890  std::vector<::FEValuesViews::Vector<dim, spacedim>> vectors;
1891  std::vector<::FEValuesViews::SymmetricTensor<2, dim, spacedim>>
1892  symmetric_second_order_tensors;
1893  std::vector<::FEValuesViews::Tensor<2, dim, spacedim>>
1894  second_order_tensors;
1895 
1899  Cache(const FEValuesBase<dim, spacedim> &fe_values);
1900  };
1901  } // namespace FEValuesViews
1902 } // namespace internal
1903 
1904 
1905 
2008 template <int dim, int spacedim>
2009 class FEValuesBase : public Subscriptor
2010 {
2011 public:
2015  static const unsigned int dimension = dim;
2016 
2020  static const unsigned int space_dimension = spacedim;
2021 
2025  const unsigned int n_quadrature_points;
2026 
2032  const unsigned int dofs_per_cell;
2033 
2034 
2042  FEValuesBase(const unsigned int n_q_points,
2043  const unsigned int dofs_per_cell,
2044  const UpdateFlags update_flags,
2047 
2052  FEValuesBase &
2053  operator=(const FEValuesBase &) = delete;
2054 
2059  FEValuesBase(const FEValuesBase &) = delete;
2060 
2064  virtual ~FEValuesBase() override;
2065 
2066 
2069 
2070 
2091  const double &
2092  shape_value(const unsigned int function_no,
2093  const unsigned int point_no) const;
2094 
2115  double
2116  shape_value_component(const unsigned int function_no,
2117  const unsigned int point_no,
2118  const unsigned int component) const;
2119 
2145  const Tensor<1, spacedim> &
2146  shape_grad(const unsigned int function_no,
2147  const unsigned int quadrature_point) const;
2148 
2166  shape_grad_component(const unsigned int function_no,
2167  const unsigned int point_no,
2168  const unsigned int component) const;
2169 
2189  const Tensor<2, spacedim> &
2190  shape_hessian(const unsigned int function_no,
2191  const unsigned int point_no) const;
2192 
2210  shape_hessian_component(const unsigned int function_no,
2211  const unsigned int point_no,
2212  const unsigned int component) const;
2213 
2233  const Tensor<3, spacedim> &
2234  shape_3rd_derivative(const unsigned int function_no,
2235  const unsigned int point_no) const;
2236 
2254  shape_3rd_derivative_component(const unsigned int function_no,
2255  const unsigned int point_no,
2256  const unsigned int component) const;
2257 
2259 
2261 
2300  template <class InputVector>
2301  void
2303  const InputVector & fe_function,
2304  std::vector<typename InputVector::value_type> &values) const;
2305 
2319  template <class InputVector>
2320  void
2322  const InputVector & fe_function,
2323  std::vector<Vector<typename InputVector::value_type>> &values) const;
2324 
2343  template <class InputVector>
2344  void
2346  const InputVector & fe_function,
2348  std::vector<typename InputVector::value_type> & values) const;
2349 
2371  template <class InputVector>
2372  void
2374  const InputVector & fe_function,
2376  std::vector<Vector<typename InputVector::value_type>> &values) const;
2377 
2378 
2409  template <class InputVector>
2410  void
2412  const InputVector & fe_function,
2414  ArrayView<std::vector<typename InputVector::value_type>> values,
2415  const bool quadrature_points_fastest) const;
2416 
2418 
2420 
2459  template <class InputVector>
2460  void
2462  const InputVector &fe_function,
2464  &gradients) const;
2465 
2482  template <class InputVector>
2483  void
2485  const InputVector &fe_function,
2486  std::vector<
2488  &gradients) const;
2489 
2496  template <class InputVector>
2497  void
2499  const InputVector & fe_function,
2502  &gradients) const;
2503 
2510  template <class InputVector>
2511  void
2513  const InputVector & fe_function,
2515  ArrayView<
2517  gradients,
2518  bool quadrature_points_fastest = false) const;
2519 
2521 
2524 
2564  template <class InputVector>
2565  void
2567  const InputVector &fe_function,
2569  &hessians) const;
2570 
2588  template <class InputVector>
2589  void
2591  const InputVector &fe_function,
2592  std::vector<
2594  & hessians,
2595  bool quadrature_points_fastest = false) const;
2596 
2601  template <class InputVector>
2602  void
2604  const InputVector & fe_function,
2607  &hessians) const;
2608 
2615  template <class InputVector>
2616  void
2618  const InputVector & fe_function,
2620  ArrayView<
2622  hessians,
2623  bool quadrature_points_fastest = false) const;
2624 
2667  template <class InputVector>
2668  void
2670  const InputVector & fe_function,
2671  std::vector<typename InputVector::value_type> &laplacians) const;
2672 
2692  template <class InputVector>
2693  void
2695  const InputVector & fe_function,
2696  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2697 
2704  template <class InputVector>
2705  void
2707  const InputVector & fe_function,
2709  std::vector<typename InputVector::value_type> & laplacians) const;
2710 
2717  template <class InputVector>
2718  void
2720  const InputVector & fe_function,
2722  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2723 
2730  template <class InputVector>
2731  void
2733  const InputVector & fe_function,
2735  std::vector<std::vector<typename InputVector::value_type>> &laplacians,
2736  bool quadrature_points_fastest = false) const;
2737 
2739 
2741 
2782  template <class InputVector>
2783  void
2785  const InputVector &fe_function,
2787  &third_derivatives) const;
2788 
2807  template <class InputVector>
2808  void
2810  const InputVector &fe_function,
2811  std::vector<
2813  & third_derivatives,
2814  bool quadrature_points_fastest = false) const;
2815 
2820  template <class InputVector>
2821  void
2823  const InputVector & fe_function,
2826  &third_derivatives) const;
2827 
2834  template <class InputVector>
2835  void
2837  const InputVector & fe_function,
2839  ArrayView<
2841  third_derivatives,
2842  bool quadrature_points_fastest = false) const;
2844 
2846 
2847 
2853  const Point<spacedim> &
2854  quadrature_point(const unsigned int q) const;
2855 
2861  const std::vector<Point<spacedim>> &
2862  get_quadrature_points() const;
2863 
2879  double
2880  JxW(const unsigned int quadrature_point) const;
2881 
2885  const std::vector<double> &
2886  get_JxW_values() const;
2887 
2895  jacobian(const unsigned int quadrature_point) const;
2896 
2903  const std::vector<DerivativeForm<1, dim, spacedim>> &
2904  get_jacobians() const;
2905 
2914  jacobian_grad(const unsigned int quadrature_point) const;
2915 
2922  const std::vector<DerivativeForm<2, dim, spacedim>> &
2923  get_jacobian_grads() const;
2924 
2933  const Tensor<3, spacedim> &
2934  jacobian_pushed_forward_grad(const unsigned int quadrature_point) const;
2935 
2942  const std::vector<Tensor<3, spacedim>> &
2944 
2953  jacobian_2nd_derivative(const unsigned int quadrature_point) const;
2954 
2961  const std::vector<DerivativeForm<3, dim, spacedim>> &
2963 
2973  const Tensor<4, spacedim> &
2975  const unsigned int quadrature_point) const;
2976 
2983  const std::vector<Tensor<4, spacedim>> &
2985 
2995  jacobian_3rd_derivative(const unsigned int quadrature_point) const;
2996 
3003  const std::vector<DerivativeForm<4, dim, spacedim>> &
3005 
3015  const Tensor<5, spacedim> &
3017  const unsigned int quadrature_point) const;
3018 
3025  const std::vector<Tensor<5, spacedim>> &
3027 
3035  inverse_jacobian(const unsigned int quadrature_point) const;
3036 
3043  const std::vector<DerivativeForm<1, spacedim, dim>> &
3044  get_inverse_jacobians() const;
3045 
3059  const Tensor<1, spacedim> &
3060  normal_vector(const unsigned int i) const;
3061 
3072  DEAL_II_DEPRECATED
3073  const std::vector<Tensor<1, spacedim>> &
3074  get_all_normal_vectors() const;
3075 
3083  const std::vector<Tensor<1, spacedim>> &
3084  get_normal_vectors() const;
3085 
3087 
3089 
3090 
3100  operator[](const FEValuesExtractors::Scalar &scalar) const;
3101 
3111  operator[](const FEValuesExtractors::Vector &vector) const;
3112 
3124 
3125 
3135  operator[](const FEValuesExtractors::Tensor<2> &tensor) const;
3136 
3138 
3140 
3141 
3145  const Mapping<dim, spacedim> &
3146  get_mapping() const;
3147 
3152  get_fe() const;
3153 
3157  UpdateFlags
3158  get_update_flags() const;
3159 
3164  get_cell() const;
3165 
3172  get_cell_similarity() const;
3173 
3178  std::size_t
3179  memory_consumption() const;
3181 
3182 
3191  std::string,
3192  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3193  << "object for which this kind of information has not been computed. What "
3194  << "information these objects compute is determined by the update_* flags you "
3195  << "pass to the constructor. Here, the operation you are attempting requires "
3196  << "the <" << arg1
3197  << "> flag to be set, but it was apparently not specified "
3198  << "upon construction.");
3199 
3208  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3209  "to the DoFHandler that provided the cell iterator do not match.");
3216  int,
3217  << "The shape function with index " << arg1
3218  << " is not primitive, i.e. it is vector-valued and "
3219  << "has more than one non-zero vector component. This "
3220  << "function cannot be called for these shape functions. "
3221  << "Maybe you want to use the same function with the "
3222  << "_component suffix?");
3223 
3232  "The given FiniteElement is not a primitive element but the requested operation "
3233  "only works for those. See FiniteElement::is_primitive() for more information.");
3234 
3235 protected:
3264  class CellIteratorBase;
3265 
3270  template <typename CI>
3273 
3279  std::unique_ptr<const CellIteratorBase> present_cell;
3280 
3288  boost::signals2::connection tria_listener_refinement;
3289 
3297  boost::signals2::connection tria_listener_mesh_transform;
3298 
3304  void
3306 
3316  void
3318  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3319 
3325 
3331  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3333 
3340 
3341 
3349 
3355  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3357 
3363  spacedim>
3365 
3366 
3371 
3380  UpdateFlags
3381  compute_update_flags(const UpdateFlags update_flags) const;
3382 
3389 
3395  void
3397  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3398 
3399 private:
3404 
3409  template <int, int>
3411  template <int, int>
3412  friend class FEValuesViews::Vector;
3413  template <int, int, int>
3414  friend class FEValuesViews::SymmetricTensor;
3415  template <int, int, int>
3416  friend class FEValuesViews::Tensor;
3417 };
3418 
3419 
3420 
3431 template <int dim, int spacedim = dim>
3432 class FEValues : public FEValuesBase<dim, spacedim>
3433 {
3434 public:
3439  static const unsigned int integral_dimension = dim;
3440 
3447  const Quadrature<dim> & quadrature,
3448  const UpdateFlags update_flags);
3449 
3456  const Quadrature<dim> & quadrature,
3457  const UpdateFlags update_flags);
3458 
3465  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3466  void
3467  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3468  level_dof_access>> &cell);
3469 
3483  void
3485 
3490  const Quadrature<dim> &
3491  get_quadrature() const;
3492 
3497  std::size_t
3498  memory_consumption() const;
3499 
3514  const FEValues<dim, spacedim> &
3515  get_present_fe_values() const;
3516 
3517 private:
3522 
3526  void
3528 
3535  void
3536  do_reinit();
3537 };
3538 
3539 
3550 template <int dim, int spacedim = dim>
3551 class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
3552 {
3553 public:
3558  static const unsigned int integral_dimension = dim - 1;
3559 
3571  FEFaceValuesBase(const unsigned int n_q_points,
3572  const unsigned int dofs_per_cell,
3573  const UpdateFlags update_flags,
3577 
3585  const Tensor<1, spacedim> &
3586  boundary_form(const unsigned int i) const;
3587 
3594  const std::vector<Tensor<1, spacedim>> &
3595  get_boundary_forms() const;
3596 
3601  unsigned int
3602  get_face_index() const;
3603 
3608  const Quadrature<dim - 1> &
3609  get_quadrature() const;
3610 
3615  std::size_t
3616  memory_consumption() const;
3617 
3618 protected:
3623  unsigned int present_face_index;
3624 
3628  const Quadrature<dim - 1> quadrature;
3629 };
3630 
3631 
3632 
3647 template <int dim, int spacedim = dim>
3648 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
3649 {
3650 public:
3655  static const unsigned int dimension = dim;
3656 
3657  static const unsigned int space_dimension = spacedim;
3658 
3663  static const unsigned int integral_dimension = dim - 1;
3664 
3672  const UpdateFlags update_flags);
3673 
3681  const UpdateFlags update_flags);
3682 
3687  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3688  void
3689  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3690  level_dof_access>> &cell,
3691  const unsigned int face_no);
3692 
3706  void
3708  const unsigned int face_no);
3709 
3725  get_present_fe_values() const;
3726 
3727 private:
3731  void
3733 
3740  void
3741  do_reinit(const unsigned int face_no);
3742 };
3743 
3744 
3762 template <int dim, int spacedim = dim>
3763 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
3764 {
3765 public:
3769  static const unsigned int dimension = dim;
3770 
3774  static const unsigned int space_dimension = spacedim;
3775 
3780  static const unsigned int integral_dimension = dim - 1;
3781 
3788  const Quadrature<dim - 1> & face_quadrature,
3789  const UpdateFlags update_flags);
3790 
3797  const Quadrature<dim - 1> & face_quadrature,
3798  const UpdateFlags update_flags);
3799 
3806  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3807  void
3808  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3809  level_dof_access>> &cell,
3810  const unsigned int face_no,
3811  const unsigned int subface_no);
3812 
3826  void
3828  const unsigned int face_no,
3829  const unsigned int subface_no);
3830 
3846  get_present_fe_values() const;
3847 
3854 
3861 
3862 private:
3866  void
3868 
3875  void
3876  do_reinit(const unsigned int face_no, const unsigned int subface_no);
3877 };
3878 
3879 
3880 #ifndef DOXYGEN
3881 
3882 
3883 /*------------------------ Inline functions: namespace FEValuesViews --------*/
3884 
3885 namespace FEValuesViews
3886 {
3887  template <int dim, int spacedim>
3888  inline typename Scalar<dim, spacedim>::value_type
3889  Scalar<dim, spacedim>::value(const unsigned int shape_function,
3890  const unsigned int q_point) const
3891  {
3892  Assert(shape_function < fe_values->fe->dofs_per_cell,
3893  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3894  Assert(
3895  fe_values->update_flags & update_values,
3897  "update_values"))));
3898 
3899  // an adaptation of the FEValuesBase::shape_value_component function
3900  // except that here we know the component as fixed and we have
3901  // pre-computed and cached a bunch of information. See the comments there.
3902  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3903  return fe_values->finite_element_output.shape_values(
3904  shape_function_data[shape_function].row_index, q_point);
3905  else
3906  return 0;
3907  }
3908 
3909 
3910 
3911  template <int dim, int spacedim>
3912  inline typename Scalar<dim, spacedim>::gradient_type
3913  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
3914  const unsigned int q_point) const
3915  {
3916  Assert(shape_function < fe_values->fe->dofs_per_cell,
3917  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3918  Assert(fe_values->update_flags & update_gradients,
3920  "update_gradients")));
3921 
3922  // an adaptation of the FEValuesBase::shape_grad_component
3923  // function except that here we know the component as fixed and we have
3924  // pre-computed and cached a bunch of information. See the comments there.
3925  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3926  return fe_values->finite_element_output
3927  .shape_gradients[shape_function_data[shape_function].row_index]
3928  [q_point];
3929  else
3930  return gradient_type();
3931  }
3932 
3933 
3934 
3935  template <int dim, int spacedim>
3936  inline typename Scalar<dim, spacedim>::hessian_type
3937  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
3938  const unsigned int q_point) const
3939  {
3940  Assert(shape_function < fe_values->fe->dofs_per_cell,
3941  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3942  Assert(fe_values->update_flags & update_hessians,
3944  "update_hessians")));
3945 
3946  // an adaptation of the FEValuesBase::shape_hessian_component
3947  // function except that here we know the component as fixed and we have
3948  // pre-computed and cached a bunch of information. See the comments there.
3949  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3950  return fe_values->finite_element_output
3951  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
3952  else
3953  return hessian_type();
3954  }
3955 
3956 
3957 
3958  template <int dim, int spacedim>
3959  inline typename Scalar<dim, spacedim>::third_derivative_type
3960  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
3961  const unsigned int q_point) const
3962  {
3963  Assert(shape_function < fe_values->fe->dofs_per_cell,
3964  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3965  Assert(fe_values->update_flags & update_3rd_derivatives,
3967  "update_3rd_derivatives")));
3968 
3969  // an adaptation of the FEValuesBase::shape_3rdderivative_component
3970  // function except that here we know the component as fixed and we have
3971  // pre-computed and cached a bunch of information. See the comments there.
3972  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3973  return fe_values->finite_element_output
3974  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
3975  [q_point];
3976  else
3977  return third_derivative_type();
3978  }
3979 
3980 
3981 
3982  template <int dim, int spacedim>
3983  inline typename Vector<dim, spacedim>::value_type
3984  Vector<dim, spacedim>::value(const unsigned int shape_function,
3985  const unsigned int q_point) const
3986  {
3987  Assert(shape_function < fe_values->fe->dofs_per_cell,
3988  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3989  Assert(fe_values->update_flags & update_values,
3991  "update_values")));
3992 
3993  // same as for the scalar case except that we have one more index
3994  const int snc =
3995  shape_function_data[shape_function].single_nonzero_component;
3996  if (snc == -2)
3997  return value_type();
3998  else if (snc != -1)
3999  {
4000  value_type return_value;
4001  return_value[shape_function_data[shape_function]
4002  .single_nonzero_component_index] =
4003  fe_values->finite_element_output.shape_values(snc, q_point);
4004  return return_value;
4005  }
4006  else
4007  {
4008  value_type return_value;
4009  for (unsigned int d = 0; d < dim; ++d)
4010  if (shape_function_data[shape_function]
4011  .is_nonzero_shape_function_component[d])
4012  return_value[d] = fe_values->finite_element_output.shape_values(
4013  shape_function_data[shape_function].row_index[d], q_point);
4014 
4015  return return_value;
4016  }
4017  }
4018 
4019 
4020 
4021  template <int dim, int spacedim>
4022  inline typename Vector<dim, spacedim>::gradient_type
4023  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4024  const unsigned int q_point) const
4025  {
4026  Assert(shape_function < fe_values->fe->dofs_per_cell,
4027  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4028  Assert(fe_values->update_flags & update_gradients,
4030  "update_gradients")));
4031 
4032  // same as for the scalar case except that we have one more index
4033  const int snc =
4034  shape_function_data[shape_function].single_nonzero_component;
4035  if (snc == -2)
4036  return gradient_type();
4037  else if (snc != -1)
4038  {
4039  gradient_type return_value;
4040  return_value[shape_function_data[shape_function]
4041  .single_nonzero_component_index] =
4042  fe_values->finite_element_output.shape_gradients[snc][q_point];
4043  return return_value;
4044  }
4045  else
4046  {
4047  gradient_type return_value;
4048  for (unsigned int d = 0; d < dim; ++d)
4049  if (shape_function_data[shape_function]
4050  .is_nonzero_shape_function_component[d])
4051  return_value[d] =
4052  fe_values->finite_element_output.shape_gradients
4053  [shape_function_data[shape_function].row_index[d]][q_point];
4054 
4055  return return_value;
4056  }
4057  }
4058 
4059 
4060 
4061  template <int dim, int spacedim>
4062  inline typename Vector<dim, spacedim>::divergence_type
4063  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4064  const unsigned int q_point) const
4065  {
4066  // this function works like in the case above
4067  Assert(shape_function < fe_values->fe->dofs_per_cell,
4068  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4069  Assert(fe_values->update_flags & update_gradients,
4071  "update_gradients")));
4072 
4073  // same as for the scalar case except that we have one more index
4074  const int snc =
4075  shape_function_data[shape_function].single_nonzero_component;
4076  if (snc == -2)
4077  return divergence_type();
4078  else if (snc != -1)
4079  return fe_values->finite_element_output
4080  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4081  .single_nonzero_component_index];
4082  else
4083  {
4084  divergence_type return_value = 0;
4085  for (unsigned int d = 0; d < dim; ++d)
4086  if (shape_function_data[shape_function]
4087  .is_nonzero_shape_function_component[d])
4088  return_value +=
4089  fe_values->finite_element_output.shape_gradients
4090  [shape_function_data[shape_function].row_index[d]][q_point][d];
4091 
4092  return return_value;
4093  }
4094  }
4095 
4096 
4097 
4098  template <int dim, int spacedim>
4099  inline typename Vector<dim, spacedim>::curl_type
4100  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4101  const unsigned int q_point) const
4102  {
4103  // this function works like in the case above
4104 
4105  Assert(shape_function < fe_values->fe->dofs_per_cell,
4106  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4107  Assert(fe_values->update_flags & update_gradients,
4109  "update_gradients")));
4110  // same as for the scalar case except that we have one more index
4111  const int snc =
4112  shape_function_data[shape_function].single_nonzero_component;
4113 
4114  if (snc == -2)
4115  return curl_type();
4116 
4117  else
4118  switch (dim)
4119  {
4120  case 1:
4121  {
4122  Assert(false,
4123  ExcMessage(
4124  "Computing the curl in 1d is not a useful operation"));
4125  return curl_type();
4126  }
4127 
4128  case 2:
4129  {
4130  if (snc != -1)
4131  {
4132  curl_type return_value;
4133 
4134  // the single nonzero component can only be zero or one in 2d
4135  if (shape_function_data[shape_function]
4136  .single_nonzero_component_index == 0)
4137  return_value[0] =
4138  -1.0 * fe_values->finite_element_output
4139  .shape_gradients[snc][q_point][1];
4140  else
4141  return_value[0] = fe_values->finite_element_output
4142  .shape_gradients[snc][q_point][0];
4143 
4144  return return_value;
4145  }
4146 
4147  else
4148  {
4149  curl_type return_value;
4150 
4151  return_value[0] = 0.0;
4152 
4153  if (shape_function_data[shape_function]
4154  .is_nonzero_shape_function_component[0])
4155  return_value[0] -=
4156  fe_values->finite_element_output
4157  .shape_gradients[shape_function_data[shape_function]
4158  .row_index[0]][q_point][1];
4159 
4160  if (shape_function_data[shape_function]
4161  .is_nonzero_shape_function_component[1])
4162  return_value[0] +=
4163  fe_values->finite_element_output
4164  .shape_gradients[shape_function_data[shape_function]
4165  .row_index[1]][q_point][0];
4166 
4167  return return_value;
4168  }
4169  }
4170 
4171  case 3:
4172  {
4173  if (snc != -1)
4174  {
4175  curl_type return_value;
4176 
4177  switch (shape_function_data[shape_function]
4178  .single_nonzero_component_index)
4179  {
4180  case 0:
4181  {
4182  return_value[0] = 0;
4183  return_value[1] = fe_values->finite_element_output
4184  .shape_gradients[snc][q_point][2];
4185  return_value[2] =
4186  -1.0 * fe_values->finite_element_output
4187  .shape_gradients[snc][q_point][1];
4188  return return_value;
4189  }
4190 
4191  case 1:
4192  {
4193  return_value[0] =
4194  -1.0 * fe_values->finite_element_output
4195  .shape_gradients[snc][q_point][2];
4196  return_value[1] = 0;
4197  return_value[2] = fe_values->finite_element_output
4198  .shape_gradients[snc][q_point][0];
4199  return return_value;
4200  }
4201 
4202  default:
4203  {
4204  return_value[0] = fe_values->finite_element_output
4205  .shape_gradients[snc][q_point][1];
4206  return_value[1] =
4207  -1.0 * fe_values->finite_element_output
4208  .shape_gradients[snc][q_point][0];
4209  return_value[2] = 0;
4210  return return_value;
4211  }
4212  }
4213  }
4214 
4215  else
4216  {
4217  curl_type return_value;
4218 
4219  for (unsigned int i = 0; i < dim; ++i)
4220  return_value[i] = 0.0;
4221 
4222  if (shape_function_data[shape_function]
4223  .is_nonzero_shape_function_component[0])
4224  {
4225  return_value[1] +=
4226  fe_values->finite_element_output
4227  .shape_gradients[shape_function_data[shape_function]
4228  .row_index[0]][q_point][2];
4229  return_value[2] -=
4230  fe_values->finite_element_output
4231  .shape_gradients[shape_function_data[shape_function]
4232  .row_index[0]][q_point][1];
4233  }
4234 
4235  if (shape_function_data[shape_function]
4236  .is_nonzero_shape_function_component[1])
4237  {
4238  return_value[0] -=
4239  fe_values->finite_element_output
4240  .shape_gradients[shape_function_data[shape_function]
4241  .row_index[1]][q_point][2];
4242  return_value[2] +=
4243  fe_values->finite_element_output
4244  .shape_gradients[shape_function_data[shape_function]
4245  .row_index[1]][q_point][0];
4246  }
4247 
4248  if (shape_function_data[shape_function]
4249  .is_nonzero_shape_function_component[2])
4250  {
4251  return_value[0] +=
4252  fe_values->finite_element_output
4253  .shape_gradients[shape_function_data[shape_function]
4254  .row_index[2]][q_point][1];
4255  return_value[1] -=
4256  fe_values->finite_element_output
4257  .shape_gradients[shape_function_data[shape_function]
4258  .row_index[2]][q_point][0];
4259  }
4260 
4261  return return_value;
4262  }
4263  }
4264  }
4265  // should not end up here
4266  Assert(false, ExcInternalError());
4267  return curl_type();
4268  }
4269 
4270 
4271 
4272  template <int dim, int spacedim>
4273  inline typename Vector<dim, spacedim>::hessian_type
4274  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
4275  const unsigned int q_point) const
4276  {
4277  // this function works like in the case above
4278  Assert(shape_function < fe_values->fe->dofs_per_cell,
4279  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4280  Assert(fe_values->update_flags & update_hessians,
4282  "update_hessians")));
4283 
4284  // same as for the scalar case except that we have one more index
4285  const int snc =
4286  shape_function_data[shape_function].single_nonzero_component;
4287  if (snc == -2)
4288  return hessian_type();
4289  else if (snc != -1)
4290  {
4291  hessian_type return_value;
4292  return_value[shape_function_data[shape_function]
4293  .single_nonzero_component_index] =
4294  fe_values->finite_element_output.shape_hessians[snc][q_point];
4295  return return_value;
4296  }
4297  else
4298  {
4299  hessian_type return_value;
4300  for (unsigned int d = 0; d < dim; ++d)
4301  if (shape_function_data[shape_function]
4302  .is_nonzero_shape_function_component[d])
4303  return_value[d] =
4304  fe_values->finite_element_output.shape_hessians
4305  [shape_function_data[shape_function].row_index[d]][q_point];
4306 
4307  return return_value;
4308  }
4309  }
4310 
4311 
4312 
4313  template <int dim, int spacedim>
4314  inline typename Vector<dim, spacedim>::third_derivative_type
4315  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
4316  const unsigned int q_point) const
4317  {
4318  // this function works like in the case above
4319  Assert(shape_function < fe_values->fe->dofs_per_cell,
4320  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4321  Assert(fe_values->update_flags & update_3rd_derivatives,
4323  "update_3rd_derivatives")));
4324 
4325  // same as for the scalar case except that we have one more index
4326  const int snc =
4327  shape_function_data[shape_function].single_nonzero_component;
4328  if (snc == -2)
4329  return third_derivative_type();
4330  else if (snc != -1)
4331  {
4332  third_derivative_type return_value;
4333  return_value[shape_function_data[shape_function]
4334  .single_nonzero_component_index] =
4335  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
4336  return return_value;
4337  }
4338  else
4339  {
4340  third_derivative_type return_value;
4341  for (unsigned int d = 0; d < dim; ++d)
4342  if (shape_function_data[shape_function]
4343  .is_nonzero_shape_function_component[d])
4344  return_value[d] =
4345  fe_values->finite_element_output.shape_3rd_derivatives
4346  [shape_function_data[shape_function].row_index[d]][q_point];
4347 
4348  return return_value;
4349  }
4350  }
4351 
4352 
4353 
4354  namespace internal
4355  {
4360  inline ::SymmetricTensor<2, 1>
4361  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
4362  {
4363  Assert(n < 1, ExcIndexRange(n, 0, 1));
4364  (void)n;
4365 
4366  return {{t[0]}};
4367  }
4368 
4369 
4370 
4371  inline ::SymmetricTensor<2, 2>
4372  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
4373  {
4374  switch (n)
4375  {
4376  case 0:
4377  {
4378  return {{t[0], 0, t[1] / 2}};
4379  }
4380  case 1:
4381  {
4382  return {{0, t[1], t[0] / 2}};
4383  }
4384  default:
4385  {
4386  Assert(false, ExcIndexRange(n, 0, 2));
4387  return {};
4388  }
4389  }
4390  }
4391 
4392 
4393 
4394  inline ::SymmetricTensor<2, 3>
4395  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
4396  {
4397  switch (n)
4398  {
4399  case 0:
4400  {
4401  return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
4402  }
4403  case 1:
4404  {
4405  return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
4406  }
4407  case 2:
4408  {
4409  return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
4410  }
4411  default:
4412  {
4413  Assert(false, ExcIndexRange(n, 0, 3));
4414  return {};
4415  }
4416  }
4417  }
4418  } // namespace internal
4419 
4420 
4421 
4422  template <int dim, int spacedim>
4423  inline typename Vector<dim, spacedim>::symmetric_gradient_type
4424  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
4425  const unsigned int q_point) const
4426  {
4427  Assert(shape_function < fe_values->fe->dofs_per_cell,
4428  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4429  Assert(fe_values->update_flags & update_gradients,
4431  "update_gradients")));
4432 
4433  // same as for the scalar case except that we have one more index
4434  const int snc =
4435  shape_function_data[shape_function].single_nonzero_component;
4436  if (snc == -2)
4437  return symmetric_gradient_type();
4438  else if (snc != -1)
4439  return internal::symmetrize_single_row(
4440  shape_function_data[shape_function].single_nonzero_component_index,
4441  fe_values->finite_element_output.shape_gradients[snc][q_point]);
4442  else
4443  {
4444  gradient_type return_value;
4445  for (unsigned int d = 0; d < dim; ++d)
4446  if (shape_function_data[shape_function]
4447  .is_nonzero_shape_function_component[d])
4448  return_value[d] =
4449  fe_values->finite_element_output.shape_gradients
4450  [shape_function_data[shape_function].row_index[d]][q_point];
4451 
4452  return symmetrize(return_value);
4453  }
4454  }
4455 
4456 
4457 
4458  template <int dim, int spacedim>
4460  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
4461  const unsigned int q_point) const
4462  {
4463  Assert(shape_function < fe_values->fe->dofs_per_cell,
4464  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4465  Assert(fe_values->update_flags & update_values,
4467  "update_values")));
4468 
4469  // similar to the vector case where we have more then one index and we need
4470  // to convert between unrolled and component indexing for tensors
4471  const int snc =
4472  shape_function_data[shape_function].single_nonzero_component;
4473 
4474  if (snc == -2)
4475  {
4476  // shape function is zero for the selected components
4477  return value_type();
4478  }
4479  else if (snc != -1)
4480  {
4481  value_type return_value;
4482  const unsigned int comp =
4483  shape_function_data[shape_function].single_nonzero_component_index;
4484  return_value[value_type::unrolled_to_component_indices(comp)] =
4485  fe_values->finite_element_output.shape_values(snc, q_point);
4486  return return_value;
4487  }
4488  else
4489  {
4490  value_type return_value;
4491  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
4492  if (shape_function_data[shape_function]
4493  .is_nonzero_shape_function_component[d])
4494  return_value[value_type::unrolled_to_component_indices(d)] =
4495  fe_values->finite_element_output.shape_values(
4496  shape_function_data[shape_function].row_index[d], q_point);
4497  return return_value;
4498  }
4499  }
4500 
4501 
4502 
4503  template <int dim, int spacedim>
4506  const unsigned int shape_function,
4507  const unsigned int q_point) const
4508  {
4509  Assert(shape_function < fe_values->fe->dofs_per_cell,
4510  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4511  Assert(fe_values->update_flags & update_gradients,
4513  "update_gradients")));
4514 
4515  const int snc =
4516  shape_function_data[shape_function].single_nonzero_component;
4517 
4518  if (snc == -2)
4519  {
4520  // shape function is zero for the selected components
4521  return divergence_type();
4522  }
4523  else if (snc != -1)
4524  {
4525  // we have a single non-zero component when the symmetric tensor is
4526  // represented in unrolled form. this implies we potentially have
4527  // two non-zero components when represented in component form! we
4528  // will only have one non-zero entry if the non-zero component lies on
4529  // the diagonal of the tensor.
4530  //
4531  // the divergence of a second-order tensor is a first order tensor.
4532  //
4533  // assume the second-order tensor is A with components A_{ij}. then
4534  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
4535  // entries in the tensorial representation. define the
4536  // divergence as:
4537  // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
4538  // (which is incidentally also
4539  // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
4540  // In both cases, a sum is implied.
4541  //
4542  // Now, we know the nonzero component in unrolled form: it is indicated
4543  // by 'snc'. we can figure out which tensor components belong to this:
4544  const unsigned int comp =
4545  shape_function_data[shape_function].single_nonzero_component_index;
4546  const unsigned int ii =
4547  value_type::unrolled_to_component_indices(comp)[0];
4548  const unsigned int jj =
4549  value_type::unrolled_to_component_indices(comp)[1];
4550 
4551  // given the form of the divergence above, if ii=jj there is only a
4552  // single nonzero component of the full tensor and the gradient
4553  // equals
4554  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
4555  // all other entries of 'b' are zero
4556  //
4557  // on the other hand, if ii!=jj, then there are two nonzero entries in
4558  // the full tensor and
4559  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
4560  // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
4561  // again, all other entries of 'b' are zero
4562  const ::Tensor<1, spacedim> &phi_grad =
4563  fe_values->finite_element_output.shape_gradients[snc][q_point];
4564 
4565  divergence_type return_value;
4566  return_value[ii] = phi_grad[jj];
4567 
4568  if (ii != jj)
4569  return_value[jj] = phi_grad[ii];
4570 
4571  return return_value;
4572  }
4573  else
4574  {
4575  Assert(false, ExcNotImplemented());
4576  divergence_type return_value;
4577  return return_value;
4578  }
4579  }
4580 
4581 
4582 
4583  template <int dim, int spacedim>
4584  inline typename Tensor<2, dim, spacedim>::value_type
4585  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
4586  const unsigned int q_point) const
4587  {
4588  Assert(shape_function < fe_values->fe->dofs_per_cell,
4589  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4590  Assert(fe_values->update_flags & update_values,
4592  "update_values")));
4593 
4594  // similar to the vector case where we have more then one index and we need
4595  // to convert between unrolled and component indexing for tensors
4596  const int snc =
4597  shape_function_data[shape_function].single_nonzero_component;
4598 
4599  if (snc == -2)
4600  {
4601  // shape function is zero for the selected components
4602  return value_type();
4603  }
4604  else if (snc != -1)
4605  {
4606  value_type return_value;
4607  const unsigned int comp =
4608  shape_function_data[shape_function].single_nonzero_component_index;
4609  const TableIndices<2> indices =
4611  return_value[indices] =
4612  fe_values->finite_element_output.shape_values(snc, q_point);
4613  return return_value;
4614  }
4615  else
4616  {
4617  value_type return_value;
4618  for (unsigned int d = 0; d < dim * dim; ++d)
4619  if (shape_function_data[shape_function]
4620  .is_nonzero_shape_function_component[d])
4621  {
4622  const TableIndices<2> indices =
4624  return_value[indices] =
4625  fe_values->finite_element_output.shape_values(
4626  shape_function_data[shape_function].row_index[d], q_point);
4627  }
4628  return return_value;
4629  }
4630  }
4631 
4632 
4633 
4634  template <int dim, int spacedim>
4636  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
4637  const unsigned int q_point) const
4638  {
4639  Assert(shape_function < fe_values->fe->dofs_per_cell,
4640  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4641  Assert(fe_values->update_flags & update_gradients,
4643  "update_gradients")));
4644 
4645  const int snc =
4646  shape_function_data[shape_function].single_nonzero_component;
4647 
4648  if (snc == -2)
4649  {
4650  // shape function is zero for the selected components
4651  return divergence_type();
4652  }
4653  else if (snc != -1)
4654  {
4655  // we have a single non-zero component when the tensor is
4656  // represented in unrolled form.
4657  //
4658  // the divergence of a second-order tensor is a first order tensor.
4659  //
4660  // assume the second-order tensor is A with components A_{ij},
4661  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
4662  //
4663  // Now, we know the nonzero component in unrolled form: it is indicated
4664  // by 'snc'. we can figure out which tensor components belong to this:
4665  const unsigned int comp =
4666  shape_function_data[shape_function].single_nonzero_component_index;
4667  const TableIndices<2> indices =
4669  const unsigned int ii = indices[0];
4670  const unsigned int jj = indices[1];
4671 
4672  const ::Tensor<1, spacedim> &phi_grad =
4673  fe_values->finite_element_output.shape_gradients[snc][q_point];
4674 
4675  divergence_type return_value;
4676  // note that we contract \nabla from the right
4677  return_value[ii] = phi_grad[jj];
4678 
4679  return return_value;
4680  }
4681  else
4682  {
4683  Assert(false, ExcNotImplemented());
4684  divergence_type return_value;
4685  return return_value;
4686  }
4687  }
4688 
4689 
4690 
4691  template <int dim, int spacedim>
4693  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
4694  const unsigned int q_point) const
4695  {
4696  Assert(shape_function < fe_values->fe->dofs_per_cell,
4697  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4698  Assert(fe_values->update_flags & update_gradients,
4700  "update_gradients")));
4701 
4702  const int snc =
4703  shape_function_data[shape_function].single_nonzero_component;
4704 
4705  if (snc == -2)
4706  {
4707  // shape function is zero for the selected components
4708  return gradient_type();
4709  }
4710  else if (snc != -1)
4711  {
4712  // we have a single non-zero component when the tensor is
4713  // represented in unrolled form.
4714  //
4715  // the gradient of a second-order tensor is a third order tensor.
4716  //
4717  // assume the second-order tensor is A with components A_{ij},
4718  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
4719  //
4720  // Now, we know the nonzero component in unrolled form: it is indicated
4721  // by 'snc'. we can figure out which tensor components belong to this:
4722  const unsigned int comp =
4723  shape_function_data[shape_function].single_nonzero_component_index;
4724  const TableIndices<2> indices =
4726  const unsigned int ii = indices[0];
4727  const unsigned int jj = indices[1];
4728 
4729  const ::Tensor<1, spacedim> &phi_grad =
4730  fe_values->finite_element_output.shape_gradients[snc][q_point];
4731 
4732  gradient_type return_value;
4733  return_value[ii][jj] = phi_grad;
4734 
4735  return return_value;
4736  }
4737  else
4738  {
4739  Assert(false, ExcNotImplemented());
4740  gradient_type return_value;
4741  return return_value;
4742  }
4743  }
4744 
4745 } // namespace FEValuesViews
4746 
4747 
4748 
4749 /*---------------------- Inline functions: FEValuesBase ---------------------*/
4750 
4751 
4752 
4753 template <int dim, int spacedim>
4755  operator[](const FEValuesExtractors::Scalar &scalar) const
4756 {
4757  Assert(scalar.component < fe_values_views_cache.scalars.size(),
4758  ExcIndexRange(scalar.component,
4759  0,
4760  fe_values_views_cache.scalars.size()));
4761 
4762  return fe_values_views_cache.scalars[scalar.component];
4763 }
4764 
4765 
4766 
4767 template <int dim, int spacedim>
4769  operator[](const FEValuesExtractors::Vector &vector) const
4770 {
4771  Assert(vector.first_vector_component < fe_values_views_cache.vectors.size(),
4773  0,
4774  fe_values_views_cache.vectors.size()));
4775 
4776  return fe_values_views_cache.vectors[vector.first_vector_component];
4777 }
4778 
4779 
4780 
4781 template <int dim, int spacedim>
4785 {
4786  Assert(
4787  tensor.first_tensor_component <
4788  fe_values_views_cache.symmetric_second_order_tensors.size(),
4790  0,
4791  fe_values_views_cache.symmetric_second_order_tensors.size()));
4792 
4793  return fe_values_views_cache
4794  .symmetric_second_order_tensors[tensor.first_tensor_component];
4795 }
4796 
4797 
4798 
4799 template <int dim, int spacedim>
4802  operator[](const FEValuesExtractors::Tensor<2> &tensor) const
4803 {
4805  fe_values_views_cache.second_order_tensors.size(),
4807  0,
4808  fe_values_views_cache.second_order_tensors.size()));
4809 
4810  return fe_values_views_cache
4811  .second_order_tensors[tensor.first_tensor_component];
4812 }
4813 
4814 
4815 
4816 template <int dim, int spacedim>
4817 inline const double &
4818 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
4819  const unsigned int j) const
4820 {
4821  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4822  Assert(this->update_flags & update_values,
4823  ExcAccessToUninitializedField("update_values"));
4824  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4825  Assert(present_cell.get() != nullptr,
4826  ExcMessage("FEValues object is not reinit'ed to any cell"));
4827  // if the entire FE is primitive,
4828  // then we can take a short-cut:
4829  if (fe->is_primitive())
4830  return this->finite_element_output.shape_values(i, j);
4831  else
4832  {
4833  // otherwise, use the mapping
4834  // between shape function
4835  // numbers and rows. note that
4836  // by the assertions above, we
4837  // know that this particular
4838  // shape function is primitive,
4839  // so we can call
4840  // system_to_component_index
4841  const unsigned int row =
4842  this->finite_element_output
4843  .shape_function_to_row_table[i * fe->n_components() +
4844  fe->system_to_component_index(i).first];
4845  return this->finite_element_output.shape_values(row, j);
4846  }
4847 }
4848 
4849 
4850 
4851 template <int dim, int spacedim>
4852 inline double
4854  const unsigned int i,
4855  const unsigned int j,
4856  const unsigned int component) const
4857 {
4858  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4859  Assert(this->update_flags & update_values,
4860  ExcAccessToUninitializedField("update_values"));
4861  Assert(component < fe->n_components(),
4862  ExcIndexRange(component, 0, fe->n_components()));
4863  Assert(present_cell.get() != nullptr,
4864  ExcMessage("FEValues object is not reinit'ed to any cell"));
4865 
4866  // check whether the shape function
4867  // is non-zero at all within
4868  // this component:
4869  if (fe->get_nonzero_components(i)[component] == false)
4870  return 0;
4871 
4872  // look up the right row in the
4873  // table and take the data from
4874  // there
4875  const unsigned int row =
4876  this->finite_element_output
4877  .shape_function_to_row_table[i * fe->n_components() + component];
4878  return this->finite_element_output.shape_values(row, j);
4879 }
4880 
4881 
4882 
4883 template <int dim, int spacedim>
4884 inline const Tensor<1, spacedim> &
4885 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
4886  const unsigned int j) const
4887 {
4888  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4889  Assert(this->update_flags & update_gradients,
4890  ExcAccessToUninitializedField("update_gradients"));
4891  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4892  Assert(present_cell.get() != nullptr,
4893  ExcMessage("FEValues object is not reinit'ed to any cell"));
4894  // if the entire FE is primitive,
4895  // then we can take a short-cut:
4896  if (fe->is_primitive())
4897  return this->finite_element_output.shape_gradients[i][j];
4898  else
4899  {
4900  // otherwise, use the mapping
4901  // between shape function
4902  // numbers and rows. note that
4903  // by the assertions above, we
4904  // know that this particular
4905  // shape function is primitive,
4906  // so we can call
4907  // system_to_component_index
4908  const unsigned int row =
4909  this->finite_element_output
4910  .shape_function_to_row_table[i * fe->n_components() +
4911  fe->system_to_component_index(i).first];
4912  return this->finite_element_output.shape_gradients[row][j];
4913  }
4914 }
4915 
4916 
4917 
4918 template <int dim, int spacedim>
4919 inline Tensor<1, spacedim>
4921  const unsigned int i,
4922  const unsigned int j,
4923  const unsigned int component) const
4924 {
4925  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4926  Assert(this->update_flags & update_gradients,
4927  ExcAccessToUninitializedField("update_gradients"));
4928  Assert(component < fe->n_components(),
4929  ExcIndexRange(component, 0, fe->n_components()));
4930  Assert(present_cell.get() != nullptr,
4931  ExcMessage("FEValues object is not reinit'ed to any cell"));
4932  // check whether the shape function
4933  // is non-zero at all within
4934  // this component:
4935  if (fe->get_nonzero_components(i)[component] == false)
4936  return Tensor<1, spacedim>();
4937 
4938  // look up the right row in the
4939  // table and take the data from
4940  // there
4941  const unsigned int row =
4942  this->finite_element_output
4943  .shape_function_to_row_table[i * fe->n_components() + component];
4944  return this->finite_element_output.shape_gradients[row][j];
4945 }
4946 
4947 
4948 
4949 template <int dim, int spacedim>
4950 inline const Tensor<2, spacedim> &
4951 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
4952  const unsigned int j) const
4953 {
4954  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4955  Assert(this->update_flags & update_hessians,
4956  ExcAccessToUninitializedField("update_hessians"));
4957  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4958  Assert(present_cell.get() != nullptr,
4959  ExcMessage("FEValues object is not reinit'ed to any cell"));
4960  // if the entire FE is primitive,
4961  // then we can take a short-cut:
4962  if (fe->is_primitive())
4963  return this->finite_element_output.shape_hessians[i][j];
4964  else
4965  {
4966  // otherwise, use the mapping
4967  // between shape function
4968  // numbers and rows. note that
4969  // by the assertions above, we
4970  // know that this particular
4971  // shape function is primitive,
4972  // so we can call
4973  // system_to_component_index
4974  const unsigned int row =
4975  this->finite_element_output
4976  .shape_function_to_row_table[i * fe->n_components() +
4977  fe->system_to_component_index(i).first];
4978  return this->finite_element_output.shape_hessians[row][j];
4979  }
4980 }
4981 
4982 
4983 
4984 template <int dim, int spacedim>
4985 inline Tensor<2, spacedim>
4987  const unsigned int i,
4988  const unsigned int j,
4989  const unsigned int component) const
4990 {
4991  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4992  Assert(this->update_flags & update_hessians,
4993  ExcAccessToUninitializedField("update_hessians"));
4994  Assert(component < fe->n_components(),
4995  ExcIndexRange(component, 0, fe->n_components()));
4996  Assert(present_cell.get() != nullptr,
4997  ExcMessage("FEValues object is not reinit'ed to any cell"));
4998  // check whether the shape function
4999  // is non-zero at all within
5000  // this component:
5001  if (fe->get_nonzero_components(i)[component] == false)
5002  return Tensor<2, spacedim>();
5003 
5004  // look up the right row in the
5005  // table and take the data from
5006  // there
5007  const unsigned int row =
5008  this->finite_element_output
5009  .shape_function_to_row_table[i * fe->n_components() + component];
5010  return this->finite_element_output.shape_hessians[row][j];
5011 }
5012 
5013 
5014 
5015 template <int dim, int spacedim>
5016 inline const Tensor<3, spacedim> &
5018  const unsigned int j) const
5019 {
5020  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5021  Assert(this->update_flags & update_hessians,
5022  ExcAccessToUninitializedField("update_3rd_derivatives"));
5023  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5024  Assert(present_cell.get() != nullptr,
5025  ExcMessage("FEValues object is not reinit'ed to any cell"));
5026  // if the entire FE is primitive,
5027  // then we can take a short-cut:
5028  if (fe->is_primitive())
5029  return this->finite_element_output.shape_3rd_derivatives[i][j];
5030  else
5031  {
5032  // otherwise, use the mapping
5033  // between shape function
5034  // numbers and rows. note that
5035  // by the assertions above, we
5036  // know that this particular
5037  // shape function is primitive,
5038  // so we can call
5039  // system_to_component_index
5040  const unsigned int row =
5041  this->finite_element_output
5042  .shape_function_to_row_table[i * fe->n_components() +
5043  fe->system_to_component_index(i).first];
5044  return this->finite_element_output.shape_3rd_derivatives[row][j];
5045  }
5046 }
5047 
5048 
5049 
5050 template <int dim, int spacedim>
5051 inline Tensor<3, spacedim>
5053  const unsigned int i,
5054  const unsigned int j,
5055  const unsigned int component) const
5056 {
5057  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5058  Assert(this->update_flags & update_hessians,
5059  ExcAccessToUninitializedField("update_3rd_derivatives"));
5060  Assert(component < fe->n_components(),
5061  ExcIndexRange(component, 0, fe->n_components()));
5062  Assert(present_cell.get() != nullptr,
5063  ExcMessage("FEValues object is not reinit'ed to any cell"));
5064  // check whether the shape function
5065  // is non-zero at all within
5066  // this component:
5067  if (fe->get_nonzero_components(i)[component] == false)
5068  return Tensor<3, spacedim>();
5069 
5070  // look up the right row in the
5071  // table and take the data from
5072  // there
5073  const unsigned int row =
5074  this->finite_element_output
5075  .shape_function_to_row_table[i * fe->n_components() + component];
5076  return this->finite_element_output.shape_3rd_derivatives[row][j];
5077 }
5078 
5079 
5080 
5081 template <int dim, int spacedim>
5082 inline const FiniteElement<dim, spacedim> &
5084 {
5085  return *fe;
5086 }
5087 
5088 
5089 
5090 template <int dim, int spacedim>
5091 inline const Mapping<dim, spacedim> &
5093 {
5094  return *mapping;
5095 }
5096 
5097 
5098 
5099 template <int dim, int spacedim>
5100 inline UpdateFlags
5102 {
5103  return this->update_flags;
5104 }
5105 
5106 
5107 
5108 template <int dim, int spacedim>
5109 inline const std::vector<Point<spacedim>> &
5111 {
5112  Assert(this->update_flags & update_quadrature_points,
5113  ExcAccessToUninitializedField("update_quadrature_points"));
5114  Assert(present_cell.get() != nullptr,
5115  ExcMessage("FEValues object is not reinit'ed to any cell"));
5116  return this->mapping_output.quadrature_points;
5117 }
5118 
5119 
5120 
5121 template <int dim, int spacedim>
5122 inline const std::vector<double> &
5124 {
5125  Assert(this->update_flags & update_JxW_values,
5126  ExcAccessToUninitializedField("update_JxW_values"));
5127  Assert(present_cell.get() != nullptr,
5128  ExcMessage("FEValues object is not reinit'ed to any cell"));
5129  return this->mapping_output.JxW_values;
5130 }
5131 
5132 
5133 
5134 template <int dim, int spacedim>
5135 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5137 {
5138  Assert(this->update_flags & update_jacobians,
5139  ExcAccessToUninitializedField("update_jacobians"));
5140  Assert(present_cell.get() != nullptr,
5141  ExcMessage("FEValues object is not reinit'ed to any cell"));
5142  return this->mapping_output.jacobians;
5143 }
5144 
5145 
5146 
5147 template <int dim, int spacedim>
5148 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5150 {
5151  Assert(this->update_flags & update_jacobian_grads,
5152  ExcAccessToUninitializedField("update_jacobians_grads"));
5153  Assert(present_cell.get() != nullptr,
5154  ExcMessage("FEValues object is not reinit'ed to any cell"));
5155  return this->mapping_output.jacobian_grads;
5156 }
5157 
5158 
5159 
5160 template <int dim, int spacedim>
5161 inline const Tensor<3, spacedim> &
5163  const unsigned int i) const
5164 {
5165  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5166  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5167  Assert(present_cell.get() != nullptr,
5168  ExcMessage("FEValues object is not reinit'ed to any cell"));
5169  return this->mapping_output.jacobian_pushed_forward_grads[i];
5170 }
5171 
5172 
5173 
5174 template <int dim, int spacedim>
5175 inline const std::vector<Tensor<3, spacedim>> &
5177 {
5178  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5179  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5180  Assert(present_cell.get() != nullptr,
5181  ExcMessage("FEValues object is not reinit'ed to any cell"));
5182  return this->mapping_output.jacobian_pushed_forward_grads;
5183 }
5184 
5185 
5186 
5187 template <int dim, int spacedim>
5188 inline const DerivativeForm<3, dim, spacedim> &
5189 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5190 {
5191  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5192  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5193  Assert(present_cell.get() != nullptr,
5194  ExcMessage("FEValues object is not reinit'ed to any cell"));
5195  return this->mapping_output.jacobian_2nd_derivatives[i];
5196 }
5197 
5198 
5199 
5200 template <int dim, int spacedim>
5201 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5203 {
5204  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5205  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5206  Assert(present_cell.get() != nullptr,
5207  ExcMessage("FEValues object is not reinit'ed to any cell"));
5208  return this->mapping_output.jacobian_2nd_derivatives;
5209 }
5210 
5211 
5212 
5213 template <int dim, int spacedim>
5214 inline const Tensor<4, spacedim> &
5216  const unsigned int i) const
5217 {
5220  "update_jacobian_pushed_forward_2nd_derivatives"));
5221  Assert(present_cell.get() != nullptr,
5222  ExcMessage("FEValues object is not reinit'ed to any cell"));
5223  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5224 }
5225 
5226 
5227 
5228 template <int dim, int spacedim>
5229 inline const std::vector<Tensor<4, spacedim>> &
5231 {
5232  Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
5234  "update_jacobian_pushed_forward_2nd_derivatives"));
5235  Assert(present_cell.get() != nullptr,
5236  ExcMessage("FEValues object is not reinit'ed to any cell"));
5237  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5238 }
5239 
5240 
5241 
5242 template <int dim, int spacedim>
5243 inline const DerivativeForm<4, dim, spacedim> &
5244 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5245 {
5246  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5247  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5248  Assert(present_cell.get() != nullptr,
5249  ExcMessage("FEValues object is not reinit'ed to any cell"));
5250  return this->mapping_output.jacobian_3rd_derivatives[i];
5251 }
5252 
5253 
5254 
5255 template <int dim, int spacedim>
5256 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5258 {
5259  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5260  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5261  Assert(present_cell.get() != nullptr,
5262  ExcMessage("FEValues object is not reinit'ed to any cell"));
5263  return this->mapping_output.jacobian_3rd_derivatives;
5264 }
5265 
5266 
5267 
5268 template <int dim, int spacedim>
5269 inline const Tensor<5, spacedim> &
5271  const unsigned int i) const
5272 {
5275  "update_jacobian_pushed_forward_3rd_derivatives"));
5276  Assert(present_cell.get() != nullptr,
5277  ExcMessage("FEValues object is not reinit'ed to any cell"));
5278  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
5279 }
5280 
5281 
5282 
5283 template <int dim, int spacedim>
5284 inline const std::vector<Tensor<5, spacedim>> &
5286 {
5287  Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
5289  "update_jacobian_pushed_forward_3rd_derivatives"));
5290  Assert(present_cell.get() != nullptr,
5291  ExcMessage("FEValues object is not reinit'ed to any cell"));
5292  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
5293 }
5294 
5295 
5296 
5297 template <int dim, int spacedim>
5298 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
5300 {
5301  Assert(this->update_flags & update_inverse_jacobians,
5302  ExcAccessToUninitializedField("update_inverse_jacobians"));
5303  Assert(present_cell.get() != nullptr,
5304  ExcMessage("FEValues object is not reinit'ed to any cell"));
5305  return this->mapping_output.inverse_jacobians;
5306 }
5307 
5308 
5309 
5310 template <int dim, int spacedim>
5311 inline const Point<spacedim> &
5312 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
5313 {
5314  Assert(this->update_flags & update_quadrature_points,
5315  ExcAccessToUninitializedField("update_quadrature_points"));
5316  Assert(i < this->mapping_output.quadrature_points.size(),
5317  ExcIndexRange(i, 0, this->mapping_output.quadrature_points.size()));
5318  Assert(present_cell.get() != nullptr,
5319  ExcMessage("FEValues object is not reinit'ed to any cell"));
5320 
5321  return this->mapping_output.quadrature_points[i];
5322 }
5323 
5324 
5325 
5326 template <int dim, int spacedim>
5327 inline double
5328 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
5329 {
5330  Assert(this->update_flags & update_JxW_values,
5331  ExcAccessToUninitializedField("update_JxW_values"));
5332  Assert(i < this->mapping_output.JxW_values.size(),
5333  ExcIndexRange(i, 0, this->mapping_output.JxW_values.size()));
5334  Assert(present_cell.get() != nullptr,
5335  ExcMessage("FEValues object is not reinit'ed to any cell"));
5336 
5337  return this->mapping_output.JxW_values[i];
5338 }
5339 
5340 
5341 
5342 template <int dim, int spacedim>
5343 inline const DerivativeForm<1, dim, spacedim> &
5344 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
5345 {
5346  Assert(this->update_flags & update_jacobians,
5347  ExcAccessToUninitializedField("update_jacobians"));
5348  Assert(i < this->mapping_output.jacobians.size(),
5349  ExcIndexRange(i, 0, this->mapping_output.jacobians.size()));
5350  Assert(present_cell.get() != nullptr,
5351  ExcMessage("FEValues object is not reinit'ed to any cell"));
5352 
5353  return this->mapping_output.jacobians[i];
5354 }
5355 
5356 
5357 
5358 template <int dim, int spacedim>
5359 inline const DerivativeForm<2, dim, spacedim> &
5360 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
5361 {
5362  Assert(this->update_flags & update_jacobian_grads,
5363  ExcAccessToUninitializedField("update_jacobians_grads"));
5364  Assert(i < this->mapping_output.jacobian_grads.size(),
5365  ExcIndexRange(i, 0, this->mapping_output.jacobian_grads.size()));
5366  Assert(present_cell.get() != nullptr,
5367  ExcMessage("FEValues object is not reinit'ed to any cell"));
5368 
5369  return this->mapping_output.jacobian_grads[i];
5370 }
5371 
5372 
5373 
5374 template <int dim, int spacedim>
5375 inline const DerivativeForm<1, spacedim, dim> &
5376 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
5377 {
5378  Assert(this->update_flags & update_inverse_jacobians,
5379  ExcAccessToUninitializedField("update_inverse_jacobians"));
5380  Assert(i < this->mapping_output.inverse_jacobians.size(),
5381  ExcIndexRange(i, 0, this->mapping_output.inverse_jacobians.size()));
5382  Assert(present_cell.get() != nullptr,
5383  ExcMessage("FEValues object is not reinit'ed to any cell"));
5384 
5385  return this->mapping_output.inverse_jacobians[i];
5386 }
5387 
5388 
5389 
5390 template <int dim, int spacedim>
5391 inline const Tensor<1, spacedim> &
5392 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
5393 {
5394  Assert(this->update_flags & update_normal_vectors,
5396  "update_normal_vectors")));
5397  Assert(i < this->mapping_output.normal_vectors.size(),
5398  ExcIndexRange(i, 0, this->mapping_output.normal_vectors.size()));
5399  Assert(present_cell.get() != nullptr,
5400  ExcMessage("FEValues object is not reinit'ed to any cell"));
5401 
5402  return this->mapping_output.normal_vectors[i];
5403 }
5404 
5405 
5406 
5407 /*--------------------- Inline functions: FEValues --------------------------*/
5408 
5409 
5410 template <int dim, int spacedim>
5411 inline const Quadrature<dim> &
5413 {
5414  return quadrature;
5415 }
5416 
5417 
5418 
5419 template <int dim, int spacedim>
5420 inline const FEValues<dim, spacedim> &
5422 {
5423  return *this;
5424 }
5425 
5426 
5427 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
5428 
5429 
5430 template <int dim, int spacedim>
5431 inline unsigned int
5433 {
5434  return present_face_index;
5435 }
5436 
5437 
5438 /*----------------------- Inline functions: FE*FaceValues -------------------*/
5439 
5440 template <int dim, int spacedim>
5441 inline const Quadrature<dim - 1> &
5443 {
5444  return quadrature;
5445 }
5446 
5447 
5448 
5449 template <int dim, int spacedim>
5450 inline const FEFaceValues<dim, spacedim> &
5452 {
5453  return *this;
5454 }
5455 
5456 
5457 
5458 template <int dim, int spacedim>
5459 inline const FESubfaceValues<dim, spacedim> &
5461 {
5462  return *this;
5463 }
5464 
5465 
5466 
5467 template <int dim, int spacedim>
5468 inline const Tensor<1, spacedim> &
5469 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
5470 {
5471  Assert(i < this->mapping_output.boundary_forms.size(),
5472  ExcIndexRange(i, 0, this->mapping_output.boundary_forms.size()));
5473  Assert(this->update_flags & update_boundary_forms,
5475  "update_boundary_forms")));
5476 
5477  return this->mapping_output.boundary_forms[i];
5478 }
5479 
5480 #endif // DOXYGEN
5481 
5482 DEAL_II_NAMESPACE_CLOSE
5483 
5484 #endif
Vector & operator=(const Vector< dim, spacedim > &)=delete
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:2165
virtual ~FEValuesBase() override
Definition: fe_values.cc:3104
void get_function_values(const InputVector &fe_function, std::vector< typename InputVector::value_type > &values) const
Definition: fe_values.cc:3577
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:2078
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
typename ProductType< Number, typename Vector< dim, spacedim >::curl_type >::type curl_type
Definition: fe_values.h:694
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:1743
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:4939
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:1574
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:3356
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:3388
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:1595
::Tensor< 1, spacedim > divergence_type
Definition: fe_values.h:1558
static const unsigned int dimension
Definition: fe_values.h:3769
const FEFaceValues< dim, spacedim > & get_present_fe_values() const
Cache(const FEValuesBase< dim, spacedim > &fe_values)
Definition: fe_values.cc:2596
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1292
std::size_t memory_consumption() const
Definition: fe_values.cc:4658
unsigned int present_face_index
Definition: fe_values.h:3623
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:1687
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:4081
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1518
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:1830
static::ExceptionBase & ExcAccessToUninitializedField()
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:534
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:1997
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:2032
typename::internal::CurlType< spacedim >::type curl_type
Definition: fe_values.h:625
const unsigned int component
Definition: fe_values.h:540
const Quadrature< dim-1 > quadrature
Definition: fe_values.h:3628
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:3081
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:3718
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:3655
unsigned int get_face_index() const
static::ExceptionBase & ExcAccessToUninitializedField(std::string arg1)
static::ExceptionBase & ExcReinitCalledWithBoundaryFace()
Outer normal vector, not normalized.
typename ProductType< Number, typename Scalar< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:211
std::unique_ptr< const CellIteratorBase > present_cell
Definition: fe_values.h:3272
::Tensor< 3, spacedim > gradient_type
Definition: fe_values.h:1564
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()
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:2109
void do_reinit(const unsigned int face_no)
Definition: fe_values.cc:4813
const Triangulation< dim, spacedim >::cell_iterator get_cell() const
Definition: fe_values.cc:4195
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
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:187
unsigned int row_index[spacedim]
Definition: fe_values.h:738
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:1885
void check_cell_similarity(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
Definition: fe_values.cc:4333
static const unsigned int space_dimension
Definition: fe_values.h:2020
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:3324
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:1774
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:1799
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:1966
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:1910
const std::vector< double > & get_JxW_values() const
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition: fe_values.h:3403
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:5012
UpdateFlags get_update_flags() const
static const unsigned int integral_dimension
Definition: fe_values.h:3780
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:2226
const Quadrature< dim > & get_quadrature() const
UpdateFlags compute_update_flags(const UpdateFlags update_flags) const
Definition: fe_values.cc:4249
typename ProductType< Number, typename Vector< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:662
static const unsigned int integral_dimension
Definition: fe_values.h:3663
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:1941
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:3944
typename ProductType< Number, typename Scalar< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:195
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1579
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:2134
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
Number value_type
Definition: vector.h:113
typename ProductType< Number, typename Vector< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:710
const std::vector< Tensor< 1, spacedim > > & get_boundary_forms() const
Definition: fe_values.cc:4646
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:518
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:1631
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1865
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:1718
#define Assert(cond, exc)
Definition: exceptions.h:1407
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:2254
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
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
Scalar & operator=(const Scalar< dim, spacedim > &)=delete
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1225
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
const Quadrature< dim > quadrature
Definition: fe_values.h:3521
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
const unsigned int first_vector_component
Definition: fe_values.h:1220
std::size_t memory_consumption() const
Definition: fe_values.cc:4229
#define DeclException0(Exception0)
Definition: exceptions.h:473
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:1662
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
void invalidate_present_cell()
Definition: fe_values.cc:4268
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition: fe_values.h:3332
typename Tensor< rank_-1, dim, Number >::tensor_type value_type
Definition: tensor.h:427
static const unsigned int integral_dimension
Definition: fe_values.h:3439
SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
const std::vector< Tensor< 1, spacedim > > & get_normal_vectors() const
Definition: fe_values.cc:4216
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:2053
Tensor()=default
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:727
const FEValues< dim, spacedim > & get_present_fe_values() const
static const unsigned int integral_dimension
Definition: fe_values.h:3558
Second derivatives of shape functions.
static::ExceptionBase & ExcFENotPrimitive()
Gradient of volume element.
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1300
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:1390
std::vector<::FEValuesViews::Scalar< dim, spacedim > > scalars
Definition: fe_values.h:1889
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access >> &cell)
Definition: fe_values.cc:4554
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:1854
::SymmetricTensor< 2, spacedim > value_type
Definition: fe_values.h:1266
const unsigned int n_quadrature_points
Definition: fe_values.h:2025
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
const FESubfaceValues< dim, spacedim > & get_present_fe_values() const
static const unsigned int dimension
Definition: fe_values.h:2015
boost::signals2::connection tria_listener_mesh_transform
Definition: fe_values.h:3297
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:654
static const unsigned int space_dimension
Definition: fe_values.h:3774
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:88
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:4626
typename ProductType< Number, typename Vector< dim, spacedim >::symmetric_gradient_type >::type symmetric_gradient_type
Definition: fe_values.h:670
typename ProductType< Number, typename Vector< dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:678
void initialize(const UpdateFlags update_flags)
Definition: fe_values.cc:4445
typename ProductType< Number, typename Vector< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:702
Shape function gradients.
Normal vectors.
void maybe_invalidate_previous_present_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
Definition: fe_values.cc:4286
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1587
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:2021
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1854
Definition: fe.h:38
void initialize(const UpdateFlags update_flags)
Definition: fe_values.cc:4890
const std::vector< DerivativeForm< 1, spacedim, dim > > & get_inverse_jacobians() const
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1507
const std::vector< Tensor< 1, spacedim > > & get_all_normal_vectors() const
Definition: fe_values.cc:4204
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:2190
static::ExceptionBase & ExcNotImplemented()
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access >> &cell, const unsigned int face_no)
Definition: fe_values.cc:4761
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1214
boost::signals2::connection tria_listener_refinement
Definition: fe_values.h:3288
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3364
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int quadrature_point) const
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:686
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:4714
::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:3339
FEValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4411
const Mapping< dim, spacedim > & get_mapping() const
void do_reinit()
Definition: fe_values.cc:4581
FEValuesBase & operator=(const FEValuesBase &)=delete
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:4854
Tensor< 1, spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:203
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:1606
FEFaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim-1 > &quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4678
CellSimilarity::Similarity get_cell_similarity() const
Definition: fe_values.cc:4388
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:2285
typename ProductType< Number, typename Scalar< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:219
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:545
const DerivativeForm< 2, dim, spacedim > & jacobian_grad(const unsigned int quadrature_point) const
UpdateFlags update_flags
Definition: fe_values.h:3370
std::size_t memory_consumption() const
Definition: fe_values.cc:4615
static::ExceptionBase & ExcInternalError()
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
Definition: fe_values.h:3348
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:3831