Reference documentation for deal.II version Git 169fded 2018-12-12 16:55:28 +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 &
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 &
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> &);
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> &);
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 
2048 
2052  virtual ~FEValuesBase() override;
2053 
2054 
2057 
2058 
2079  const double &
2080  shape_value(const unsigned int function_no,
2081  const unsigned int point_no) const;
2082 
2103  double
2104  shape_value_component(const unsigned int function_no,
2105  const unsigned int point_no,
2106  const unsigned int component) const;
2107 
2133  const Tensor<1, spacedim> &
2134  shape_grad(const unsigned int function_no,
2135  const unsigned int quadrature_point) const;
2136 
2154  shape_grad_component(const unsigned int function_no,
2155  const unsigned int point_no,
2156  const unsigned int component) const;
2157 
2177  const Tensor<2, spacedim> &
2178  shape_hessian(const unsigned int function_no,
2179  const unsigned int point_no) const;
2180 
2198  shape_hessian_component(const unsigned int function_no,
2199  const unsigned int point_no,
2200  const unsigned int component) const;
2201 
2221  const Tensor<3, spacedim> &
2222  shape_3rd_derivative(const unsigned int function_no,
2223  const unsigned int point_no) const;
2224 
2242  shape_3rd_derivative_component(const unsigned int function_no,
2243  const unsigned int point_no,
2244  const unsigned int component) const;
2245 
2247 
2249 
2288  template <class InputVector>
2289  void
2291  const InputVector & fe_function,
2292  std::vector<typename InputVector::value_type> &values) const;
2293 
2307  template <class InputVector>
2308  void
2310  const InputVector & fe_function,
2311  std::vector<Vector<typename InputVector::value_type>> &values) const;
2312 
2331  template <class InputVector>
2332  void
2334  const InputVector & fe_function,
2335  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2336  std::vector<typename InputVector::value_type> &values) const;
2337 
2359  template <class InputVector>
2360  void
2362  const InputVector & fe_function,
2363  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2364  std::vector<Vector<typename InputVector::value_type>> &values) const;
2365 
2366 
2397  template <class InputVector>
2398  void
2400  const InputVector & fe_function,
2401  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2402  VectorSlice<std::vector<std::vector<typename InputVector::value_type>>>
2403  values,
2404  const bool quadrature_points_fastest) const;
2405 
2407 
2409 
2448  template <class InputVector>
2449  void
2451  const InputVector &fe_function,
2453  &gradients) const;
2454 
2471  template <class InputVector>
2472  void
2474  const InputVector &fe_function,
2475  std::vector<
2477  &gradients) const;
2478 
2485  template <class InputVector>
2486  void
2488  const InputVector & fe_function,
2489  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2491  &gradients) const;
2492 
2499  template <class InputVector>
2500  void
2502  const InputVector & fe_function,
2503  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2504  VectorSlice<std::vector<
2506  gradients,
2507  bool quadrature_points_fastest = false) const;
2508 
2510 
2513 
2553  template <class InputVector>
2554  void
2556  const InputVector &fe_function,
2558  &hessians) const;
2559 
2577  template <class InputVector>
2578  void
2580  const InputVector &fe_function,
2581  std::vector<
2583  & hessians,
2584  bool quadrature_points_fastest = false) const;
2585 
2590  template <class InputVector>
2591  void
2593  const InputVector & fe_function,
2594  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2596  &hessians) const;
2597 
2604  template <class InputVector>
2605  void
2607  const InputVector & fe_function,
2608  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2609  VectorSlice<std::vector<
2611  hessians,
2612  bool quadrature_points_fastest = false) const;
2613 
2656  template <class InputVector>
2657  void
2659  const InputVector & fe_function,
2660  std::vector<typename InputVector::value_type> &laplacians) const;
2661 
2681  template <class InputVector>
2682  void
2684  const InputVector & fe_function,
2685  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2686 
2693  template <class InputVector>
2694  void
2696  const InputVector & fe_function,
2697  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2698  std::vector<typename InputVector::value_type> &laplacians) const;
2699 
2706  template <class InputVector>
2707  void
2709  const InputVector & fe_function,
2710  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2711  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
2712 
2719  template <class InputVector>
2720  void
2722  const InputVector & fe_function,
2723  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2724  std::vector<std::vector<typename InputVector::value_type>> & laplacians,
2725  bool quadrature_points_fastest = false) const;
2726 
2728 
2730 
2771  template <class InputVector>
2772  void
2774  const InputVector &fe_function,
2776  &third_derivatives) const;
2777 
2796  template <class InputVector>
2797  void
2799  const InputVector &fe_function,
2800  std::vector<
2802  & third_derivatives,
2803  bool quadrature_points_fastest = false) const;
2804 
2809  template <class InputVector>
2810  void
2812  const InputVector & fe_function,
2813  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2815  &third_derivatives) const;
2816 
2823  template <class InputVector>
2824  void
2826  const InputVector & fe_function,
2827  const VectorSlice<const std::vector<types::global_dof_index>> &indices,
2828  VectorSlice<std::vector<
2830  third_derivatives,
2831  bool quadrature_points_fastest = false) const;
2833 
2835 
2836 
2842  const Point<spacedim> &
2843  quadrature_point(const unsigned int q) const;
2844 
2850  const std::vector<Point<spacedim>> &
2851  get_quadrature_points() const;
2852 
2868  double
2869  JxW(const unsigned int quadrature_point) const;
2870 
2874  const std::vector<double> &
2875  get_JxW_values() const;
2876 
2884  jacobian(const unsigned int quadrature_point) const;
2885 
2892  const std::vector<DerivativeForm<1, dim, spacedim>> &
2893  get_jacobians() const;
2894 
2903  jacobian_grad(const unsigned int quadrature_point) const;
2904 
2911  const std::vector<DerivativeForm<2, dim, spacedim>> &
2912  get_jacobian_grads() const;
2913 
2922  const Tensor<3, spacedim> &
2923  jacobian_pushed_forward_grad(const unsigned int quadrature_point) const;
2924 
2931  const std::vector<Tensor<3, spacedim>> &
2933 
2942  jacobian_2nd_derivative(const unsigned int quadrature_point) const;
2943 
2950  const std::vector<DerivativeForm<3, dim, spacedim>> &
2952 
2962  const Tensor<4, spacedim> &
2964  const unsigned int quadrature_point) const;
2965 
2972  const std::vector<Tensor<4, spacedim>> &
2974 
2984  jacobian_3rd_derivative(const unsigned int quadrature_point) const;
2985 
2992  const std::vector<DerivativeForm<4, dim, spacedim>> &
2994 
3004  const Tensor<5, spacedim> &
3006  const unsigned int quadrature_point) const;
3007 
3014  const std::vector<Tensor<5, spacedim>> &
3016 
3024  inverse_jacobian(const unsigned int quadrature_point) const;
3025 
3032  const std::vector<DerivativeForm<1, spacedim, dim>> &
3033  get_inverse_jacobians() const;
3034 
3048  const Tensor<1, spacedim> &
3049  normal_vector(const unsigned int i) const;
3050 
3061  DEAL_II_DEPRECATED
3062  const std::vector<Tensor<1, spacedim>> &
3063  get_all_normal_vectors() const;
3064 
3072  const std::vector<Tensor<1, spacedim>> &
3073  get_normal_vectors() const;
3074 
3076 
3078 
3079 
3089  operator[](const FEValuesExtractors::Scalar &scalar) const;
3090 
3100  operator[](const FEValuesExtractors::Vector &vector) const;
3101 
3113 
3114 
3124  operator[](const FEValuesExtractors::Tensor<2> &tensor) const;
3125 
3127 
3129 
3130 
3134  const Mapping<dim, spacedim> &
3135  get_mapping() const;
3136 
3141  get_fe() const;
3142 
3146  UpdateFlags
3147  get_update_flags() const;
3148 
3153  get_cell() const;
3154 
3161  get_cell_similarity() const;
3162 
3167  std::size_t
3168  memory_consumption() const;
3170 
3171 
3180  std::string,
3181  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3182  << "object for which this kind of information has not been computed. What "
3183  << "information these objects compute is determined by the update_* flags you "
3184  << "pass to the constructor. Here, the operation you are attempting requires "
3185  << "the <" << arg1
3186  << "> flag to be set, but it was apparently not specified "
3187  << "upon construction.");
3188 
3197  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3198  "to the DoFHandler that provided the cell iterator do not match.");
3205  int,
3206  << "The shape function with index " << arg1
3207  << " is not primitive, i.e. it is vector-valued and "
3208  << "has more than one non-zero vector component. This "
3209  << "function cannot be called for these shape functions. "
3210  << "Maybe you want to use the same function with the "
3211  << "_component suffix?");
3212 
3221  "The given FiniteElement is not a primitive element but the requested operation "
3222  "only works for those. See FiniteElement::is_primitive() for more information.");
3223 
3224 protected:
3253  class CellIteratorBase;
3254 
3259  template <typename CI>
3262 
3268  std::unique_ptr<const CellIteratorBase> present_cell;
3269 
3277  boost::signals2::connection tria_listener_refinement;
3278 
3286  boost::signals2::connection tria_listener_mesh_transform;
3287 
3293  void
3295 
3305  void
3307  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3308 
3314 
3320  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3322 
3329 
3330 
3338 
3344  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3346 
3352  spacedim>
3354 
3355 
3360 
3369  UpdateFlags
3370  compute_update_flags(const UpdateFlags update_flags) const;
3371 
3378 
3384  void
3386  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3387 
3388 private:
3393  FEValuesBase(const FEValuesBase &);
3394 
3399  FEValuesBase &
3400  operator=(const FEValuesBase &);
3401 
3406 
3411  template <int, int>
3413  template <int, int>
3414  friend class FEValuesViews::Vector;
3415  template <int, int, int>
3416  friend class FEValuesViews::SymmetricTensor;
3417  template <int, int, int>
3418  friend class FEValuesViews::Tensor;
3419 };
3420 
3421 
3422 
3433 template <int dim, int spacedim = dim>
3434 class FEValues : public FEValuesBase<dim, spacedim>
3435 {
3436 public:
3441  static const unsigned int integral_dimension = dim;
3442 
3449  const Quadrature<dim> & quadrature,
3450  const UpdateFlags update_flags);
3451 
3458  const Quadrature<dim> & quadrature,
3459  const UpdateFlags update_flags);
3460 
3467  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3468  void
3469  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3470  level_dof_access>> &cell);
3471 
3485  void
3487 
3492  const Quadrature<dim> &
3493  get_quadrature() const;
3494 
3499  std::size_t
3500  memory_consumption() const;
3501 
3516  const FEValues<dim, spacedim> &
3517  get_present_fe_values() const;
3518 
3519 private:
3524 
3528  void
3530 
3537  void
3538  do_reinit();
3539 };
3540 
3541 
3552 template <int dim, int spacedim = dim>
3553 class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
3554 {
3555 public:
3560  static const unsigned int integral_dimension = dim - 1;
3561 
3573  FEFaceValuesBase(const unsigned int n_q_points,
3574  const unsigned int dofs_per_cell,
3575  const UpdateFlags update_flags,
3579 
3587  const Tensor<1, spacedim> &
3588  boundary_form(const unsigned int i) const;
3589 
3596  const std::vector<Tensor<1, spacedim>> &
3597  get_boundary_forms() const;
3598 
3603  unsigned int
3604  get_face_index() const;
3605 
3610  const Quadrature<dim - 1> &
3611  get_quadrature() const;
3612 
3617  std::size_t
3618  memory_consumption() const;
3619 
3620 protected:
3625  unsigned int present_face_index;
3626 
3630  const Quadrature<dim - 1> quadrature;
3631 };
3632 
3633 
3634 
3649 template <int dim, int spacedim = dim>
3650 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
3651 {
3652 public:
3657  static const unsigned int dimension = dim;
3658 
3659  static const unsigned int space_dimension = spacedim;
3660 
3665  static const unsigned int integral_dimension = dim - 1;
3666 
3674  const UpdateFlags update_flags);
3675 
3683  const UpdateFlags update_flags);
3684 
3689  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3690  void
3691  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3692  level_dof_access>> &cell,
3693  const unsigned int face_no);
3694 
3708  void
3710  const unsigned int face_no);
3711 
3727  get_present_fe_values() const;
3728 
3729 private:
3733  void
3735 
3742  void
3743  do_reinit(const unsigned int face_no);
3744 };
3745 
3746 
3764 template <int dim, int spacedim = dim>
3765 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
3766 {
3767 public:
3771  static const unsigned int dimension = dim;
3772 
3776  static const unsigned int space_dimension = spacedim;
3777 
3782  static const unsigned int integral_dimension = dim - 1;
3783 
3790  const Quadrature<dim - 1> & face_quadrature,
3791  const UpdateFlags update_flags);
3792 
3799  const Quadrature<dim - 1> & face_quadrature,
3800  const UpdateFlags update_flags);
3801 
3808  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3809  void
3810  reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType<dim, spacedim>,
3811  level_dof_access>> &cell,
3812  const unsigned int face_no,
3813  const unsigned int subface_no);
3814 
3828  void
3830  const unsigned int face_no,
3831  const unsigned int subface_no);
3832 
3848  get_present_fe_values() const;
3849 
3856 
3863 
3864 private:
3868  void
3870 
3877  void
3878  do_reinit(const unsigned int face_no, const unsigned int subface_no);
3879 };
3880 
3881 
3882 #ifndef DOXYGEN
3883 
3884 
3885 /*------------------------ Inline functions: namespace FEValuesViews --------*/
3886 
3887 namespace FEValuesViews
3888 {
3889  template <int dim, int spacedim>
3890  inline typename Scalar<dim, spacedim>::value_type
3891  Scalar<dim, spacedim>::value(const unsigned int shape_function,
3892  const unsigned int q_point) const
3893  {
3894  Assert(shape_function < fe_values->fe->dofs_per_cell,
3895  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3896  Assert(
3897  fe_values->update_flags & update_values,
3899  "update_values"))));
3900 
3901  // an adaptation of the FEValuesBase::shape_value_component function
3902  // except that here we know the component as fixed and we have
3903  // pre-computed and cached a bunch of information. See the comments there.
3904  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3905  return fe_values->finite_element_output.shape_values(
3906  shape_function_data[shape_function].row_index, q_point);
3907  else
3908  return 0;
3909  }
3910 
3911 
3912 
3913  template <int dim, int spacedim>
3914  inline typename Scalar<dim, spacedim>::gradient_type
3915  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
3916  const unsigned int q_point) const
3917  {
3918  Assert(shape_function < fe_values->fe->dofs_per_cell,
3919  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3920  Assert(fe_values->update_flags & update_gradients,
3922  "update_gradients")));
3923 
3924  // an adaptation of the FEValuesBase::shape_grad_component
3925  // function except that here we know the component as fixed and we have
3926  // pre-computed and cached a bunch of information. See the comments there.
3927  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3928  return fe_values->finite_element_output
3929  .shape_gradients[shape_function_data[shape_function].row_index]
3930  [q_point];
3931  else
3932  return gradient_type();
3933  }
3934 
3935 
3936 
3937  template <int dim, int spacedim>
3938  inline typename Scalar<dim, spacedim>::hessian_type
3939  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
3940  const unsigned int q_point) const
3941  {
3942  Assert(shape_function < fe_values->fe->dofs_per_cell,
3943  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3944  Assert(fe_values->update_flags & update_hessians,
3946  "update_hessians")));
3947 
3948  // an adaptation of the FEValuesBase::shape_hessian_component
3949  // function except that here we know the component as fixed and we have
3950  // pre-computed and cached a bunch of information. See the comments there.
3951  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3952  return fe_values->finite_element_output
3953  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
3954  else
3955  return hessian_type();
3956  }
3957 
3958 
3959 
3960  template <int dim, int spacedim>
3961  inline typename Scalar<dim, spacedim>::third_derivative_type
3962  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
3963  const unsigned int q_point) const
3964  {
3965  Assert(shape_function < fe_values->fe->dofs_per_cell,
3966  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3967  Assert(fe_values->update_flags & update_3rd_derivatives,
3969  "update_3rd_derivatives")));
3970 
3971  // an adaptation of the FEValuesBase::shape_3rdderivative_component
3972  // function except that here we know the component as fixed and we have
3973  // pre-computed and cached a bunch of information. See the comments there.
3974  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3975  return fe_values->finite_element_output
3976  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
3977  [q_point];
3978  else
3979  return third_derivative_type();
3980  }
3981 
3982 
3983 
3984  template <int dim, int spacedim>
3985  inline typename Vector<dim, spacedim>::value_type
3986  Vector<dim, spacedim>::value(const unsigned int shape_function,
3987  const unsigned int q_point) const
3988  {
3989  Assert(shape_function < fe_values->fe->dofs_per_cell,
3990  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
3991  Assert(fe_values->update_flags & update_values,
3993  "update_values")));
3994 
3995  // same as for the scalar case except that we have one more index
3996  const int snc =
3997  shape_function_data[shape_function].single_nonzero_component;
3998  if (snc == -2)
3999  return value_type();
4000  else if (snc != -1)
4001  {
4002  value_type return_value;
4003  return_value[shape_function_data[shape_function]
4004  .single_nonzero_component_index] =
4005  fe_values->finite_element_output.shape_values(snc, q_point);
4006  return return_value;
4007  }
4008  else
4009  {
4010  value_type return_value;
4011  for (unsigned int d = 0; d < dim; ++d)
4012  if (shape_function_data[shape_function]
4013  .is_nonzero_shape_function_component[d])
4014  return_value[d] = fe_values->finite_element_output.shape_values(
4015  shape_function_data[shape_function].row_index[d], q_point);
4016 
4017  return return_value;
4018  }
4019  }
4020 
4021 
4022 
4023  template <int dim, int spacedim>
4024  inline typename Vector<dim, spacedim>::gradient_type
4025  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4026  const unsigned int q_point) const
4027  {
4028  Assert(shape_function < fe_values->fe->dofs_per_cell,
4029  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4030  Assert(fe_values->update_flags & update_gradients,
4032  "update_gradients")));
4033 
4034  // same as for the scalar case except that we have one more index
4035  const int snc =
4036  shape_function_data[shape_function].single_nonzero_component;
4037  if (snc == -2)
4038  return gradient_type();
4039  else if (snc != -1)
4040  {
4041  gradient_type return_value;
4042  return_value[shape_function_data[shape_function]
4043  .single_nonzero_component_index] =
4044  fe_values->finite_element_output.shape_gradients[snc][q_point];
4045  return return_value;
4046  }
4047  else
4048  {
4049  gradient_type return_value;
4050  for (unsigned int d = 0; d < dim; ++d)
4051  if (shape_function_data[shape_function]
4052  .is_nonzero_shape_function_component[d])
4053  return_value[d] =
4054  fe_values->finite_element_output.shape_gradients
4055  [shape_function_data[shape_function].row_index[d]][q_point];
4056 
4057  return return_value;
4058  }
4059  }
4060 
4061 
4062 
4063  template <int dim, int spacedim>
4064  inline typename Vector<dim, spacedim>::divergence_type
4065  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4066  const unsigned int q_point) const
4067  {
4068  // this function works like in the case above
4069  Assert(shape_function < fe_values->fe->dofs_per_cell,
4070  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4071  Assert(fe_values->update_flags & update_gradients,
4073  "update_gradients")));
4074 
4075  // same as for the scalar case except that we have one more index
4076  const int snc =
4077  shape_function_data[shape_function].single_nonzero_component;
4078  if (snc == -2)
4079  return divergence_type();
4080  else if (snc != -1)
4081  return fe_values->finite_element_output
4082  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4083  .single_nonzero_component_index];
4084  else
4085  {
4086  divergence_type return_value = 0;
4087  for (unsigned int d = 0; d < dim; ++d)
4088  if (shape_function_data[shape_function]
4089  .is_nonzero_shape_function_component[d])
4090  return_value +=
4091  fe_values->finite_element_output.shape_gradients
4092  [shape_function_data[shape_function].row_index[d]][q_point][d];
4093 
4094  return return_value;
4095  }
4096  }
4097 
4098 
4099 
4100  template <int dim, int spacedim>
4101  inline typename Vector<dim, spacedim>::curl_type
4102  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4103  const unsigned int q_point) const
4104  {
4105  // this function works like in the case above
4106 
4107  Assert(shape_function < fe_values->fe->dofs_per_cell,
4108  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4109  Assert(fe_values->update_flags & update_gradients,
4111  "update_gradients")));
4112  // same as for the scalar case except that we have one more index
4113  const int snc =
4114  shape_function_data[shape_function].single_nonzero_component;
4115 
4116  if (snc == -2)
4117  return curl_type();
4118 
4119  else
4120  switch (dim)
4121  {
4122  case 1:
4123  {
4124  Assert(false,
4125  ExcMessage(
4126  "Computing the curl in 1d is not a useful operation"));
4127  return curl_type();
4128  }
4129 
4130  case 2:
4131  {
4132  if (snc != -1)
4133  {
4134  curl_type return_value;
4135 
4136  // the single nonzero component can only be zero or one in 2d
4137  if (shape_function_data[shape_function]
4138  .single_nonzero_component_index == 0)
4139  return_value[0] =
4140  -1.0 * fe_values->finite_element_output
4141  .shape_gradients[snc][q_point][1];
4142  else
4143  return_value[0] = fe_values->finite_element_output
4144  .shape_gradients[snc][q_point][0];
4145 
4146  return return_value;
4147  }
4148 
4149  else
4150  {
4151  curl_type return_value;
4152 
4153  return_value[0] = 0.0;
4154 
4155  if (shape_function_data[shape_function]
4156  .is_nonzero_shape_function_component[0])
4157  return_value[0] -=
4158  fe_values->finite_element_output
4159  .shape_gradients[shape_function_data[shape_function]
4160  .row_index[0]][q_point][1];
4161 
4162  if (shape_function_data[shape_function]
4163  .is_nonzero_shape_function_component[1])
4164  return_value[0] +=
4165  fe_values->finite_element_output
4166  .shape_gradients[shape_function_data[shape_function]
4167  .row_index[1]][q_point][0];
4168 
4169  return return_value;
4170  }
4171  }
4172 
4173  case 3:
4174  {
4175  if (snc != -1)
4176  {
4177  curl_type return_value;
4178 
4179  switch (shape_function_data[shape_function]
4180  .single_nonzero_component_index)
4181  {
4182  case 0:
4183  {
4184  return_value[0] = 0;
4185  return_value[1] = fe_values->finite_element_output
4186  .shape_gradients[snc][q_point][2];
4187  return_value[2] =
4188  -1.0 * fe_values->finite_element_output
4189  .shape_gradients[snc][q_point][1];
4190  return return_value;
4191  }
4192 
4193  case 1:
4194  {
4195  return_value[0] =
4196  -1.0 * fe_values->finite_element_output
4197  .shape_gradients[snc][q_point][2];
4198  return_value[1] = 0;
4199  return_value[2] = fe_values->finite_element_output
4200  .shape_gradients[snc][q_point][0];
4201  return return_value;
4202  }
4203 
4204  default:
4205  {
4206  return_value[0] = fe_values->finite_element_output
4207  .shape_gradients[snc][q_point][1];
4208  return_value[1] =
4209  -1.0 * fe_values->finite_element_output
4210  .shape_gradients[snc][q_point][0];
4211  return_value[2] = 0;
4212  return return_value;
4213  }
4214  }
4215  }
4216 
4217  else
4218  {
4219  curl_type return_value;
4220 
4221  for (unsigned int i = 0; i < dim; ++i)
4222  return_value[i] = 0.0;
4223 
4224  if (shape_function_data[shape_function]
4225  .is_nonzero_shape_function_component[0])
4226  {
4227  return_value[1] +=
4228  fe_values->finite_element_output
4229  .shape_gradients[shape_function_data[shape_function]
4230  .row_index[0]][q_point][2];
4231  return_value[2] -=
4232  fe_values->finite_element_output
4233  .shape_gradients[shape_function_data[shape_function]
4234  .row_index[0]][q_point][1];
4235  }
4236 
4237  if (shape_function_data[shape_function]
4238  .is_nonzero_shape_function_component[1])
4239  {
4240  return_value[0] -=
4241  fe_values->finite_element_output
4242  .shape_gradients[shape_function_data[shape_function]
4243  .row_index[1]][q_point][2];
4244  return_value[2] +=
4245  fe_values->finite_element_output
4246  .shape_gradients[shape_function_data[shape_function]
4247  .row_index[1]][q_point][0];
4248  }
4249 
4250  if (shape_function_data[shape_function]
4251  .is_nonzero_shape_function_component[2])
4252  {
4253  return_value[0] +=
4254  fe_values->finite_element_output
4255  .shape_gradients[shape_function_data[shape_function]
4256  .row_index[2]][q_point][1];
4257  return_value[1] -=
4258  fe_values->finite_element_output
4259  .shape_gradients[shape_function_data[shape_function]
4260  .row_index[2]][q_point][0];
4261  }
4262 
4263  return return_value;
4264  }
4265  }
4266  }
4267  // should not end up here
4268  Assert(false, ExcInternalError());
4269  return curl_type();
4270  }
4271 
4272 
4273 
4274  template <int dim, int spacedim>
4275  inline typename Vector<dim, spacedim>::hessian_type
4276  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
4277  const unsigned int q_point) const
4278  {
4279  // this function works like in the case above
4280  Assert(shape_function < fe_values->fe->dofs_per_cell,
4281  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4282  Assert(fe_values->update_flags & update_hessians,
4284  "update_hessians")));
4285 
4286  // same as for the scalar case except that we have one more index
4287  const int snc =
4288  shape_function_data[shape_function].single_nonzero_component;
4289  if (snc == -2)
4290  return hessian_type();
4291  else if (snc != -1)
4292  {
4293  hessian_type return_value;
4294  return_value[shape_function_data[shape_function]
4295  .single_nonzero_component_index] =
4296  fe_values->finite_element_output.shape_hessians[snc][q_point];
4297  return return_value;
4298  }
4299  else
4300  {
4301  hessian_type return_value;
4302  for (unsigned int d = 0; d < dim; ++d)
4303  if (shape_function_data[shape_function]
4304  .is_nonzero_shape_function_component[d])
4305  return_value[d] =
4306  fe_values->finite_element_output.shape_hessians
4307  [shape_function_data[shape_function].row_index[d]][q_point];
4308 
4309  return return_value;
4310  }
4311  }
4312 
4313 
4314 
4315  template <int dim, int spacedim>
4316  inline typename Vector<dim, spacedim>::third_derivative_type
4317  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
4318  const unsigned int q_point) const
4319  {
4320  // this function works like in the case above
4321  Assert(shape_function < fe_values->fe->dofs_per_cell,
4322  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4323  Assert(fe_values->update_flags & update_3rd_derivatives,
4325  "update_3rd_derivatives")));
4326 
4327  // same as for the scalar case except that we have one more index
4328  const int snc =
4329  shape_function_data[shape_function].single_nonzero_component;
4330  if (snc == -2)
4331  return third_derivative_type();
4332  else if (snc != -1)
4333  {
4334  third_derivative_type return_value;
4335  return_value[shape_function_data[shape_function]
4336  .single_nonzero_component_index] =
4337  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
4338  return return_value;
4339  }
4340  else
4341  {
4342  third_derivative_type return_value;
4343  for (unsigned int d = 0; d < dim; ++d)
4344  if (shape_function_data[shape_function]
4345  .is_nonzero_shape_function_component[d])
4346  return_value[d] =
4347  fe_values->finite_element_output.shape_3rd_derivatives
4348  [shape_function_data[shape_function].row_index[d]][q_point];
4349 
4350  return return_value;
4351  }
4352  }
4353 
4354 
4355 
4356  namespace internal
4357  {
4362  inline ::SymmetricTensor<2, 1>
4363  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
4364  {
4365  Assert(n < 1, ExcIndexRange(n, 0, 1));
4366  (void)n;
4367 
4368  const double array[1] = {t[0]};
4369  return ::SymmetricTensor<2, 1>(array);
4370  }
4371 
4372 
4373 
4374  inline ::SymmetricTensor<2, 2>
4375  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
4376  {
4377  switch (n)
4378  {
4379  case 0:
4380  {
4381  const double array[3] = {t[0], 0, t[1] / 2};
4382  return ::SymmetricTensor<2, 2>(array);
4383  }
4384  case 1:
4385  {
4386  const double array[3] = {0, t[1], t[0] / 2};
4387  return ::SymmetricTensor<2, 2>(array);
4388  }
4389  default:
4390  {
4391  Assert(false, ExcIndexRange(n, 0, 2));
4392  return ::SymmetricTensor<2, 2>();
4393  }
4394  }
4395  }
4396 
4397 
4398 
4399  inline ::SymmetricTensor<2, 3>
4400  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
4401  {
4402  switch (n)
4403  {
4404  case 0:
4405  {
4406  const double array[6] = {t[0], 0, 0, t[1] / 2, t[2] / 2, 0};
4407  return ::SymmetricTensor<2, 3>(array);
4408  }
4409  case 1:
4410  {
4411  const double array[6] = {0, t[1], 0, t[0] / 2, 0, t[2] / 2};
4412  return ::SymmetricTensor<2, 3>(array);
4413  }
4414  case 2:
4415  {
4416  const double array[6] = {0, 0, t[2], 0, t[0] / 2, t[1] / 2};
4417  return ::SymmetricTensor<2, 3>(array);
4418  }
4419  default:
4420  {
4421  Assert(false, ExcIndexRange(n, 0, 3));
4422  return ::SymmetricTensor<2, 3>();
4423  }
4424  }
4425  }
4426  } // namespace internal
4427 
4428 
4429 
4430  template <int dim, int spacedim>
4431  inline typename Vector<dim, spacedim>::symmetric_gradient_type
4432  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
4433  const unsigned int q_point) const
4434  {
4435  Assert(shape_function < fe_values->fe->dofs_per_cell,
4436  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4437  Assert(fe_values->update_flags & update_gradients,
4439  "update_gradients")));
4440 
4441  // same as for the scalar case except that we have one more index
4442  const int snc =
4443  shape_function_data[shape_function].single_nonzero_component;
4444  if (snc == -2)
4445  return symmetric_gradient_type();
4446  else if (snc != -1)
4447  return internal::symmetrize_single_row(
4448  shape_function_data[shape_function].single_nonzero_component_index,
4449  fe_values->finite_element_output.shape_gradients[snc][q_point]);
4450  else
4451  {
4452  gradient_type return_value;
4453  for (unsigned int d = 0; d < dim; ++d)
4454  if (shape_function_data[shape_function]
4455  .is_nonzero_shape_function_component[d])
4456  return_value[d] =
4457  fe_values->finite_element_output.shape_gradients
4458  [shape_function_data[shape_function].row_index[d]][q_point];
4459 
4460  return symmetrize(return_value);
4461  }
4462  }
4463 
4464 
4465 
4466  template <int dim, int spacedim>
4468  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
4469  const unsigned int q_point) const
4470  {
4471  Assert(shape_function < fe_values->fe->dofs_per_cell,
4472  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4473  Assert(fe_values->update_flags & update_values,
4475  "update_values")));
4476 
4477  // similar to the vector case where we have more then one index and we need
4478  // to convert between unrolled and component indexing for tensors
4479  const int snc =
4480  shape_function_data[shape_function].single_nonzero_component;
4481 
4482  if (snc == -2)
4483  {
4484  // shape function is zero for the selected components
4485  return value_type();
4486  }
4487  else if (snc != -1)
4488  {
4489  value_type return_value;
4490  const unsigned int comp =
4491  shape_function_data[shape_function].single_nonzero_component_index;
4492  return_value[value_type::unrolled_to_component_indices(comp)] =
4493  fe_values->finite_element_output.shape_values(snc, q_point);
4494  return return_value;
4495  }
4496  else
4497  {
4498  value_type return_value;
4499  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
4500  if (shape_function_data[shape_function]
4501  .is_nonzero_shape_function_component[d])
4502  return_value[value_type::unrolled_to_component_indices(d)] =
4503  fe_values->finite_element_output.shape_values(
4504  shape_function_data[shape_function].row_index[d], q_point);
4505  return return_value;
4506  }
4507  }
4508 
4509 
4510 
4511  template <int dim, int spacedim>
4514  const unsigned int shape_function,
4515  const unsigned int q_point) const
4516  {
4517  Assert(shape_function < fe_values->fe->dofs_per_cell,
4518  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4519  Assert(fe_values->update_flags & update_gradients,
4521  "update_gradients")));
4522 
4523  const int snc =
4524  shape_function_data[shape_function].single_nonzero_component;
4525 
4526  if (snc == -2)
4527  {
4528  // shape function is zero for the selected components
4529  return divergence_type();
4530  }
4531  else if (snc != -1)
4532  {
4533  // we have a single non-zero component when the symmetric tensor is
4534  // represented in unrolled form. this implies we potentially have
4535  // two non-zero components when represented in component form! we
4536  // will only have one non-zero entry if the non-zero component lies on
4537  // the diagonal of the tensor.
4538  //
4539  // the divergence of a second-order tensor is a first order tensor.
4540  //
4541  // assume the second-order tensor is A with components A_{ij}. then
4542  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
4543  // entries in the tensorial representation. define the
4544  // divergence as:
4545  // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
4546  // (which is incidentally also
4547  // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
4548  // In both cases, a sum is implied.
4549  //
4550  // Now, we know the nonzero component in unrolled form: it is indicated
4551  // by 'snc'. we can figure out which tensor components belong to this:
4552  const unsigned int comp =
4553  shape_function_data[shape_function].single_nonzero_component_index;
4554  const unsigned int ii =
4555  value_type::unrolled_to_component_indices(comp)[0];
4556  const unsigned int jj =
4557  value_type::unrolled_to_component_indices(comp)[1];
4558 
4559  // given the form of the divergence above, if ii=jj there is only a
4560  // single nonzero component of the full tensor and the gradient
4561  // equals
4562  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
4563  // all other entries of 'b' are zero
4564  //
4565  // on the other hand, if ii!=jj, then there are two nonzero entries in
4566  // the full tensor and
4567  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
4568  // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
4569  // again, all other entries of 'b' are zero
4570  const ::Tensor<1, spacedim> &phi_grad =
4571  fe_values->finite_element_output.shape_gradients[snc][q_point];
4572 
4573  divergence_type return_value;
4574  return_value[ii] = phi_grad[jj];
4575 
4576  if (ii != jj)
4577  return_value[jj] = phi_grad[ii];
4578 
4579  return return_value;
4580  }
4581  else
4582  {
4583  Assert(false, ExcNotImplemented());
4584  divergence_type return_value;
4585  return return_value;
4586  }
4587  }
4588 
4589 
4590 
4591  template <int dim, int spacedim>
4592  inline typename Tensor<2, dim, spacedim>::value_type
4593  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
4594  const unsigned int q_point) const
4595  {
4596  Assert(shape_function < fe_values->fe->dofs_per_cell,
4597  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4598  Assert(fe_values->update_flags & update_values,
4600  "update_values")));
4601 
4602  // similar to the vector case where we have more then one index and we need
4603  // to convert between unrolled and component indexing for tensors
4604  const int snc =
4605  shape_function_data[shape_function].single_nonzero_component;
4606 
4607  if (snc == -2)
4608  {
4609  // shape function is zero for the selected components
4610  return value_type();
4611  }
4612  else if (snc != -1)
4613  {
4614  value_type return_value;
4615  const unsigned int comp =
4616  shape_function_data[shape_function].single_nonzero_component_index;
4617  const TableIndices<2> indices =
4619  return_value[indices] =
4620  fe_values->finite_element_output.shape_values(snc, q_point);
4621  return return_value;
4622  }
4623  else
4624  {
4625  value_type return_value;
4626  for (unsigned int d = 0; d < dim * dim; ++d)
4627  if (shape_function_data[shape_function]
4628  .is_nonzero_shape_function_component[d])
4629  {
4630  const TableIndices<2> indices =
4632  return_value[indices] =
4633  fe_values->finite_element_output.shape_values(
4634  shape_function_data[shape_function].row_index[d], q_point);
4635  }
4636  return return_value;
4637  }
4638  }
4639 
4640 
4641 
4642  template <int dim, int spacedim>
4644  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
4645  const unsigned int q_point) const
4646  {
4647  Assert(shape_function < fe_values->fe->dofs_per_cell,
4648  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4649  Assert(fe_values->update_flags & update_gradients,
4651  "update_gradients")));
4652 
4653  const int snc =
4654  shape_function_data[shape_function].single_nonzero_component;
4655 
4656  if (snc == -2)
4657  {
4658  // shape function is zero for the selected components
4659  return divergence_type();
4660  }
4661  else if (snc != -1)
4662  {
4663  // we have a single non-zero component when the tensor is
4664  // represented in unrolled form.
4665  //
4666  // the divergence of a second-order tensor is a first order tensor.
4667  //
4668  // assume the second-order tensor is A with components A_{ij},
4669  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
4670  //
4671  // Now, we know the nonzero component in unrolled form: it is indicated
4672  // by 'snc'. we can figure out which tensor components belong to this:
4673  const unsigned int comp =
4674  shape_function_data[shape_function].single_nonzero_component_index;
4675  const TableIndices<2> indices =
4677  const unsigned int ii = indices[0];
4678  const unsigned int jj = indices[1];
4679 
4680  const ::Tensor<1, spacedim> &phi_grad =
4681  fe_values->finite_element_output.shape_gradients[snc][q_point];
4682 
4683  divergence_type return_value;
4684  // note that we contract \nabla from the right
4685  return_value[ii] = phi_grad[jj];
4686 
4687  return return_value;
4688  }
4689  else
4690  {
4691  Assert(false, ExcNotImplemented());
4692  divergence_type return_value;
4693  return return_value;
4694  }
4695  }
4696 
4697 
4698 
4699  template <int dim, int spacedim>
4701  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
4702  const unsigned int q_point) const
4703  {
4704  Assert(shape_function < fe_values->fe->dofs_per_cell,
4705  ExcIndexRange(shape_function, 0, fe_values->fe->dofs_per_cell));
4706  Assert(fe_values->update_flags & update_gradients,
4708  "update_gradients")));
4709 
4710  const int snc =
4711  shape_function_data[shape_function].single_nonzero_component;
4712 
4713  if (snc == -2)
4714  {
4715  // shape function is zero for the selected components
4716  return gradient_type();
4717  }
4718  else if (snc != -1)
4719  {
4720  // we have a single non-zero component when the tensor is
4721  // represented in unrolled form.
4722  //
4723  // the gradient of a second-order tensor is a third order tensor.
4724  //
4725  // assume the second-order tensor is A with components A_{ij},
4726  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
4727  //
4728  // Now, we know the nonzero component in unrolled form: it is indicated
4729  // by 'snc'. we can figure out which tensor components belong to this:
4730  const unsigned int comp =
4731  shape_function_data[shape_function].single_nonzero_component_index;
4732  const TableIndices<2> indices =
4734  const unsigned int ii = indices[0];
4735  const unsigned int jj = indices[1];
4736 
4737  const ::Tensor<1, spacedim> &phi_grad =
4738  fe_values->finite_element_output.shape_gradients[snc][q_point];
4739 
4740  gradient_type return_value;
4741  return_value[ii][jj] = phi_grad;
4742 
4743  return return_value;
4744  }
4745  else
4746  {
4747  Assert(false, ExcNotImplemented());
4748  gradient_type return_value;
4749  return return_value;
4750  }
4751  }
4752 
4753 } // namespace FEValuesViews
4754 
4755 
4756 
4757 /*---------------------- Inline functions: FEValuesBase ---------------------*/
4758 
4759 
4760 
4761 template <int dim, int spacedim>
4763  operator[](const FEValuesExtractors::Scalar &scalar) const
4764 {
4765  Assert(scalar.component < fe_values_views_cache.scalars.size(),
4766  ExcIndexRange(scalar.component,
4767  0,
4768  fe_values_views_cache.scalars.size()));
4769 
4770  return fe_values_views_cache.scalars[scalar.component];
4771 }
4772 
4773 
4774 
4775 template <int dim, int spacedim>
4777  operator[](const FEValuesExtractors::Vector &vector) const
4778 {
4779  Assert(vector.first_vector_component < fe_values_views_cache.vectors.size(),
4781  0,
4782  fe_values_views_cache.vectors.size()));
4783 
4784  return fe_values_views_cache.vectors[vector.first_vector_component];
4785 }
4786 
4787 
4788 
4789 template <int dim, int spacedim>
4793 {
4794  Assert(
4795  tensor.first_tensor_component <
4796  fe_values_views_cache.symmetric_second_order_tensors.size(),
4798  0,
4799  fe_values_views_cache.symmetric_second_order_tensors.size()));
4800 
4801  return fe_values_views_cache
4802  .symmetric_second_order_tensors[tensor.first_tensor_component];
4803 }
4804 
4805 
4806 
4807 template <int dim, int spacedim>
4810  operator[](const FEValuesExtractors::Tensor<2> &tensor) const
4811 {
4813  fe_values_views_cache.second_order_tensors.size(),
4815  0,
4816  fe_values_views_cache.second_order_tensors.size()));
4817 
4818  return fe_values_views_cache
4819  .second_order_tensors[tensor.first_tensor_component];
4820 }
4821 
4822 
4823 
4824 template <int dim, int spacedim>
4825 inline const double &
4826 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
4827  const unsigned int j) const
4828 {
4829  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4830  Assert(this->update_flags & update_values,
4831  ExcAccessToUninitializedField("update_values"));
4832  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4833  Assert(present_cell.get() != nullptr,
4834  ExcMessage("FEValues object is not reinit'ed to any cell"));
4835  // if the entire FE is primitive,
4836  // then we can take a short-cut:
4837  if (fe->is_primitive())
4838  return this->finite_element_output.shape_values(i, j);
4839  else
4840  {
4841  // otherwise, use the mapping
4842  // between shape function
4843  // numbers and rows. note that
4844  // by the assertions above, we
4845  // know that this particular
4846  // shape function is primitive,
4847  // so we can call
4848  // system_to_component_index
4849  const unsigned int row =
4850  this->finite_element_output
4851  .shape_function_to_row_table[i * fe->n_components() +
4852  fe->system_to_component_index(i).first];
4853  return this->finite_element_output.shape_values(row, j);
4854  }
4855 }
4856 
4857 
4858 
4859 template <int dim, int spacedim>
4860 inline double
4862  const unsigned int i,
4863  const unsigned int j,
4864  const unsigned int component) const
4865 {
4866  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4867  Assert(this->update_flags & update_values,
4868  ExcAccessToUninitializedField("update_values"));
4869  Assert(component < fe->n_components(),
4870  ExcIndexRange(component, 0, fe->n_components()));
4871  Assert(present_cell.get() != nullptr,
4872  ExcMessage("FEValues object is not reinit'ed to any cell"));
4873 
4874  // check whether the shape function
4875  // is non-zero at all within
4876  // this component:
4877  if (fe->get_nonzero_components(i)[component] == false)
4878  return 0;
4879 
4880  // look up the right row in the
4881  // table and take the data from
4882  // there
4883  const unsigned int row =
4884  this->finite_element_output
4885  .shape_function_to_row_table[i * fe->n_components() + component];
4886  return this->finite_element_output.shape_values(row, j);
4887 }
4888 
4889 
4890 
4891 template <int dim, int spacedim>
4892 inline const Tensor<1, spacedim> &
4893 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
4894  const unsigned int j) const
4895 {
4896  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4897  Assert(this->update_flags & update_gradients,
4898  ExcAccessToUninitializedField("update_gradients"));
4899  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4900  Assert(present_cell.get() != nullptr,
4901  ExcMessage("FEValues object is not reinit'ed to any cell"));
4902  // if the entire FE is primitive,
4903  // then we can take a short-cut:
4904  if (fe->is_primitive())
4905  return this->finite_element_output.shape_gradients[i][j];
4906  else
4907  {
4908  // otherwise, use the mapping
4909  // between shape function
4910  // numbers and rows. note that
4911  // by the assertions above, we
4912  // know that this particular
4913  // shape function is primitive,
4914  // so we can call
4915  // system_to_component_index
4916  const unsigned int row =
4917  this->finite_element_output
4918  .shape_function_to_row_table[i * fe->n_components() +
4919  fe->system_to_component_index(i).first];
4920  return this->finite_element_output.shape_gradients[row][j];
4921  }
4922 }
4923 
4924 
4925 
4926 template <int dim, int spacedim>
4927 inline Tensor<1, spacedim>
4929  const unsigned int i,
4930  const unsigned int j,
4931  const unsigned int component) const
4932 {
4933  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4934  Assert(this->update_flags & update_gradients,
4935  ExcAccessToUninitializedField("update_gradients"));
4936  Assert(component < fe->n_components(),
4937  ExcIndexRange(component, 0, fe->n_components()));
4938  Assert(present_cell.get() != nullptr,
4939  ExcMessage("FEValues object is not reinit'ed to any cell"));
4940  // check whether the shape function
4941  // is non-zero at all within
4942  // this component:
4943  if (fe->get_nonzero_components(i)[component] == false)
4944  return Tensor<1, spacedim>();
4945 
4946  // look up the right row in the
4947  // table and take the data from
4948  // there
4949  const unsigned int row =
4950  this->finite_element_output
4951  .shape_function_to_row_table[i * fe->n_components() + component];
4952  return this->finite_element_output.shape_gradients[row][j];
4953 }
4954 
4955 
4956 
4957 template <int dim, int spacedim>
4958 inline const Tensor<2, spacedim> &
4959 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
4960  const unsigned int j) const
4961 {
4962  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
4963  Assert(this->update_flags & update_hessians,
4964  ExcAccessToUninitializedField("update_hessians"));
4965  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
4966  Assert(present_cell.get() != nullptr,
4967  ExcMessage("FEValues object is not reinit'ed to any cell"));
4968  // if the entire FE is primitive,
4969  // then we can take a short-cut:
4970  if (fe->is_primitive())
4971  return this->finite_element_output.shape_hessians[i][j];
4972  else
4973  {
4974  // otherwise, use the mapping
4975  // between shape function
4976  // numbers and rows. note that
4977  // by the assertions above, we
4978  // know that this particular
4979  // shape function is primitive,
4980  // so we can call
4981  // system_to_component_index
4982  const unsigned int row =
4983  this->finite_element_output
4984  .shape_function_to_row_table[i * fe->n_components() +
4985  fe->system_to_component_index(i).first];
4986  return this->finite_element_output.shape_hessians[row][j];
4987  }
4988 }
4989 
4990 
4991 
4992 template <int dim, int spacedim>
4993 inline Tensor<2, spacedim>
4995  const unsigned int i,
4996  const unsigned int j,
4997  const unsigned int component) const
4998 {
4999  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5000  Assert(this->update_flags & update_hessians,
5001  ExcAccessToUninitializedField("update_hessians"));
5002  Assert(component < fe->n_components(),
5003  ExcIndexRange(component, 0, fe->n_components()));
5004  Assert(present_cell.get() != nullptr,
5005  ExcMessage("FEValues object is not reinit'ed to any cell"));
5006  // check whether the shape function
5007  // is non-zero at all within
5008  // this component:
5009  if (fe->get_nonzero_components(i)[component] == false)
5010  return Tensor<2, spacedim>();
5011 
5012  // look up the right row in the
5013  // table and take the data from
5014  // there
5015  const unsigned int row =
5016  this->finite_element_output
5017  .shape_function_to_row_table[i * fe->n_components() + component];
5018  return this->finite_element_output.shape_hessians[row][j];
5019 }
5020 
5021 
5022 
5023 template <int dim, int spacedim>
5024 inline const Tensor<3, spacedim> &
5026  const unsigned int j) const
5027 {
5028  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5029  Assert(this->update_flags & update_hessians,
5030  ExcAccessToUninitializedField("update_3rd_derivatives"));
5031  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5032  Assert(present_cell.get() != nullptr,
5033  ExcMessage("FEValues object is not reinit'ed to any cell"));
5034  // if the entire FE is primitive,
5035  // then we can take a short-cut:
5036  if (fe->is_primitive())
5037  return this->finite_element_output.shape_3rd_derivatives[i][j];
5038  else
5039  {
5040  // otherwise, use the mapping
5041  // between shape function
5042  // numbers and rows. note that
5043  // by the assertions above, we
5044  // know that this particular
5045  // shape function is primitive,
5046  // so we can call
5047  // system_to_component_index
5048  const unsigned int row =
5049  this->finite_element_output
5050  .shape_function_to_row_table[i * fe->n_components() +
5051  fe->system_to_component_index(i).first];
5052  return this->finite_element_output.shape_3rd_derivatives[row][j];
5053  }
5054 }
5055 
5056 
5057 
5058 template <int dim, int spacedim>
5059 inline Tensor<3, spacedim>
5061  const unsigned int i,
5062  const unsigned int j,
5063  const unsigned int component) const
5064 {
5065  Assert(i < fe->dofs_per_cell, ExcIndexRange(i, 0, fe->dofs_per_cell));
5066  Assert(this->update_flags & update_hessians,
5067  ExcAccessToUninitializedField("update_3rd_derivatives"));
5068  Assert(component < fe->n_components(),
5069  ExcIndexRange(component, 0, fe->n_components()));
5070  Assert(present_cell.get() != nullptr,
5071  ExcMessage("FEValues object is not reinit'ed to any cell"));
5072  // check whether the shape function
5073  // is non-zero at all within
5074  // this component:
5075  if (fe->get_nonzero_components(i)[component] == false)
5076  return Tensor<3, spacedim>();
5077 
5078  // look up the right row in the
5079  // table and take the data from
5080  // there
5081  const unsigned int row =
5082  this->finite_element_output
5083  .shape_function_to_row_table[i * fe->n_components() + component];
5084  return this->finite_element_output.shape_3rd_derivatives[row][j];
5085 }
5086 
5087 
5088 
5089 template <int dim, int spacedim>
5090 inline const FiniteElement<dim, spacedim> &
5092 {
5093  return *fe;
5094 }
5095 
5096 
5097 
5098 template <int dim, int spacedim>
5099 inline const Mapping<dim, spacedim> &
5101 {
5102  return *mapping;
5103 }
5104 
5105 
5106 
5107 template <int dim, int spacedim>
5108 inline UpdateFlags
5110 {
5111  return this->update_flags;
5112 }
5113 
5114 
5115 
5116 template <int dim, int spacedim>
5117 inline const std::vector<Point<spacedim>> &
5119 {
5120  Assert(this->update_flags & update_quadrature_points,
5121  ExcAccessToUninitializedField("update_quadrature_points"));
5122  Assert(present_cell.get() != nullptr,
5123  ExcMessage("FEValues object is not reinit'ed to any cell"));
5124  return this->mapping_output.quadrature_points;
5125 }
5126 
5127 
5128 
5129 template <int dim, int spacedim>
5130 inline const std::vector<double> &
5132 {
5133  Assert(this->update_flags & update_JxW_values,
5134  ExcAccessToUninitializedField("update_JxW_values"));
5135  Assert(present_cell.get() != nullptr,
5136  ExcMessage("FEValues object is not reinit'ed to any cell"));
5137  return this->mapping_output.JxW_values;
5138 }
5139 
5140 
5141 
5142 template <int dim, int spacedim>
5143 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5145 {
5146  Assert(this->update_flags & update_jacobians,
5147  ExcAccessToUninitializedField("update_jacobians"));
5148  Assert(present_cell.get() != nullptr,
5149  ExcMessage("FEValues object is not reinit'ed to any cell"));
5150  return this->mapping_output.jacobians;
5151 }
5152 
5153 
5154 
5155 template <int dim, int spacedim>
5156 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5158 {
5159  Assert(this->update_flags & update_jacobian_grads,
5160  ExcAccessToUninitializedField("update_jacobians_grads"));
5161  Assert(present_cell.get() != nullptr,
5162  ExcMessage("FEValues object is not reinit'ed to any cell"));
5163  return this->mapping_output.jacobian_grads;
5164 }
5165 
5166 
5167 
5168 template <int dim, int spacedim>
5169 inline const Tensor<3, spacedim> &
5171  const unsigned int i) const
5172 {
5173  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5174  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5175  Assert(present_cell.get() != nullptr,
5176  ExcMessage("FEValues object is not reinit'ed to any cell"));
5177  return this->mapping_output.jacobian_pushed_forward_grads[i];
5178 }
5179 
5180 
5181 
5182 template <int dim, int spacedim>
5183 inline const std::vector<Tensor<3, spacedim>> &
5185 {
5186  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5187  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5188  Assert(present_cell.get() != nullptr,
5189  ExcMessage("FEValues object is not reinit'ed to any cell"));
5190  return this->mapping_output.jacobian_pushed_forward_grads;
5191 }
5192 
5193 
5194 
5195 template <int dim, int spacedim>
5196 inline const DerivativeForm<3, dim, spacedim> &
5197 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5198 {
5199  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5200  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5201  Assert(present_cell.get() != nullptr,
5202  ExcMessage("FEValues object is not reinit'ed to any cell"));
5203  return this->mapping_output.jacobian_2nd_derivatives[i];
5204 }
5205 
5206 
5207 
5208 template <int dim, int spacedim>
5209 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5211 {
5212  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5213  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5214  Assert(present_cell.get() != nullptr,
5215  ExcMessage("FEValues object is not reinit'ed to any cell"));
5216  return this->mapping_output.jacobian_2nd_derivatives;
5217 }
5218 
5219 
5220 
5221 template <int dim, int spacedim>
5222 inline const Tensor<4, spacedim> &
5224  const unsigned int i) const
5225 {
5228  "update_jacobian_pushed_forward_2nd_derivatives"));
5229  Assert(present_cell.get() != nullptr,
5230  ExcMessage("FEValues object is not reinit'ed to any cell"));
5231  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5232 }
5233 
5234 
5235 
5236 template <int dim, int spacedim>
5237 inline const std::vector<Tensor<4, spacedim>> &
5239 {
5240  Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
5242  "update_jacobian_pushed_forward_2nd_derivatives"));
5243  Assert(present_cell.get() != nullptr,
5244  ExcMessage("FEValues object is not reinit'ed to any cell"));
5245  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5246 }
5247 
5248 
5249 
5250 template <int dim, int spacedim>
5251 inline const DerivativeForm<4, dim, spacedim> &
5252 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5253 {
5254  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5255  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5256  Assert(present_cell.get() != nullptr,
5257  ExcMessage("FEValues object is not reinit'ed to any cell"));
5258  return this->mapping_output.jacobian_3rd_derivatives[i];
5259 }
5260 
5261 
5262 
5263 template <int dim, int spacedim>
5264 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5266 {
5267  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5268  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5269  Assert(present_cell.get() != nullptr,
5270  ExcMessage("FEValues object is not reinit'ed to any cell"));
5271  return this->mapping_output.jacobian_3rd_derivatives;
5272 }
5273 
5274 
5275 
5276 template <int dim, int spacedim>
5277 inline const Tensor<5, spacedim> &
5279  const unsigned int i) const
5280 {
5283  "update_jacobian_pushed_forward_3rd_derivatives"));
5284  Assert(present_cell.get() != nullptr,
5285  ExcMessage("FEValues object is not reinit'ed to any cell"));
5286  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
5287 }
5288 
5289 
5290 
5291 template <int dim, int spacedim>
5292 inline const std::vector<Tensor<5, spacedim>> &
5294 {
5295  Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
5297  "update_jacobian_pushed_forward_3rd_derivatives"));
5298  Assert(present_cell.get() != nullptr,
5299  ExcMessage("FEValues object is not reinit'ed to any cell"));
5300  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
5301 }
5302 
5303 
5304 
5305 template <int dim, int spacedim>
5306 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
5308 {
5309  Assert(this->update_flags & update_inverse_jacobians,
5310  ExcAccessToUninitializedField("update_inverse_jacobians"));
5311  Assert(present_cell.get() != nullptr,
5312  ExcMessage("FEValues object is not reinit'ed to any cell"));
5313  return this->mapping_output.inverse_jacobians;
5314 }
5315 
5316 
5317 
5318 template <int dim, int spacedim>
5319 inline const Point<spacedim> &
5320 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
5321 {
5322  Assert(this->update_flags & update_quadrature_points,
5323  ExcAccessToUninitializedField("update_quadrature_points"));
5324  Assert(i < this->mapping_output.quadrature_points.size(),
5325  ExcIndexRange(i, 0, this->mapping_output.quadrature_points.size()));
5326  Assert(present_cell.get() != nullptr,
5327  ExcMessage("FEValues object is not reinit'ed to any cell"));
5328 
5329  return this->mapping_output.quadrature_points[i];
5330 }
5331 
5332 
5333 
5334 template <int dim, int spacedim>
5335 inline double
5336 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
5337 {
5338  Assert(this->update_flags & update_JxW_values,
5339  ExcAccessToUninitializedField("update_JxW_values"));
5340  Assert(i < this->mapping_output.JxW_values.size(),
5341  ExcIndexRange(i, 0, this->mapping_output.JxW_values.size()));
5342  Assert(present_cell.get() != nullptr,
5343  ExcMessage("FEValues object is not reinit'ed to any cell"));
5344 
5345  return this->mapping_output.JxW_values[i];
5346 }
5347 
5348 
5349 
5350 template <int dim, int spacedim>
5351 inline const DerivativeForm<1, dim, spacedim> &
5352 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
5353 {
5354  Assert(this->update_flags & update_jacobians,
5355  ExcAccessToUninitializedField("update_jacobians"));
5356  Assert(i < this->mapping_output.jacobians.size(),
5357  ExcIndexRange(i, 0, this->mapping_output.jacobians.size()));
5358  Assert(present_cell.get() != nullptr,
5359  ExcMessage("FEValues object is not reinit'ed to any cell"));
5360 
5361  return this->mapping_output.jacobians[i];
5362 }
5363 
5364 
5365 
5366 template <int dim, int spacedim>
5367 inline const DerivativeForm<2, dim, spacedim> &
5368 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
5369 {
5370  Assert(this->update_flags & update_jacobian_grads,
5371  ExcAccessToUninitializedField("update_jacobians_grads"));
5372  Assert(i < this->mapping_output.jacobian_grads.size(),
5373  ExcIndexRange(i, 0, this->mapping_output.jacobian_grads.size()));
5374  Assert(present_cell.get() != nullptr,
5375  ExcMessage("FEValues object is not reinit'ed to any cell"));
5376 
5377  return this->mapping_output.jacobian_grads[i];
5378 }
5379 
5380 
5381 
5382 template <int dim, int spacedim>
5383 inline const DerivativeForm<1, spacedim, dim> &
5384 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
5385 {
5386  Assert(this->update_flags & update_inverse_jacobians,
5387  ExcAccessToUninitializedField("update_inverse_jacobians"));
5388  Assert(i < this->mapping_output.inverse_jacobians.size(),
5389  ExcIndexRange(i, 0, this->mapping_output.inverse_jacobians.size()));
5390  Assert(present_cell.get() != nullptr,
5391  ExcMessage("FEValues object is not reinit'ed to any cell"));
5392 
5393  return this->mapping_output.inverse_jacobians[i];
5394 }
5395 
5396 
5397 
5398 template <int dim, int spacedim>
5399 inline const Tensor<1, spacedim> &
5400 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
5401 {
5402  Assert(this->update_flags & update_normal_vectors,
5404  "update_normal_vectors")));
5405  Assert(i < this->mapping_output.normal_vectors.size(),
5406  ExcIndexRange(i, 0, this->mapping_output.normal_vectors.size()));
5407  Assert(present_cell.get() != nullptr,
5408  ExcMessage("FEValues object is not reinit'ed to any cell"));
5409 
5410  return this->mapping_output.normal_vectors[i];
5411 }
5412 
5413 
5414 
5415 /*--------------------- Inline functions: FEValues --------------------------*/
5416 
5417 
5418 template <int dim, int spacedim>
5419 inline const Quadrature<dim> &
5421 {
5422  return quadrature;
5423 }
5424 
5425 
5426 
5427 template <int dim, int spacedim>
5428 inline const FEValues<dim, spacedim> &
5430 {
5431  return *this;
5432 }
5433 
5434 
5435 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
5436 
5437 
5438 template <int dim, int spacedim>
5439 inline unsigned int
5441 {
5442  return present_face_index;
5443 }
5444 
5445 
5446 /*----------------------- Inline functions: FE*FaceValues -------------------*/
5447 
5448 template <int dim, int spacedim>
5449 inline const Quadrature<dim - 1> &
5451 {
5452  return quadrature;
5453 }
5454 
5455 
5456 
5457 template <int dim, int spacedim>
5458 inline const FEFaceValues<dim, spacedim> &
5460 {
5461  return *this;
5462 }
5463 
5464 
5465 
5466 template <int dim, int spacedim>
5467 inline const FESubfaceValues<dim, spacedim> &
5469 {
5470  return *this;
5471 }
5472 
5473 
5474 
5475 template <int dim, int spacedim>
5476 inline const Tensor<1, spacedim> &
5477 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
5478 {
5479  Assert(i < this->mapping_output.boundary_forms.size(),
5480  ExcIndexRange(i, 0, this->mapping_output.boundary_forms.size()));
5481  Assert(this->update_flags & update_boundary_forms,
5483  "update_boundary_forms")));
5484 
5485  return this->mapping_output.boundary_forms[i];
5486 }
5487 
5488 #endif // DOXYGEN
5489 
5490 DEAL_II_NAMESPACE_CLOSE
5491 
5492 #endif
Transformed quadrature weights.
const std::vector< Tensor< 4, spacedim > > & get_jacobian_pushed_forward_2nd_derivatives() const
void get_function_hessians_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::hessian_type > &hessians) const
Definition: fe_values.cc:2210
virtual ~FEValuesBase() override
Definition: fe_values.cc:3184
void get_function_values(const InputVector &fe_function, std::vector< typename InputVector::value_type > &values) const
Definition: fe_values.cc:3657
const Tensor< 3, spacedim > & shape_3rd_derivative(const unsigned int function_no, const unsigned int point_no) const
void get_function_curls(const InputVector &fe_function, std::vector< typename ProductType< curl_type, typename InputVector::value_type >::type > &curls) const
Definition: fe_values.cc:2123
Shape function values.
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int quadrature_point) const
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:1788
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access >> &cell, const unsigned int face_no, const unsigned int subface_no)
Definition: fe_values.cc:5023
void get_function_values(const InputVector &fe_function, std::vector< typename ProductType< value_type, typename InputVector::value_type >::type > &values) const
Definition: fe_values.cc:1619
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > fe_data
Definition: fe_values.h:3345
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:3377
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:3771
const FEFaceValues< dim, spacedim > & get_present_fe_values() const
Cache(const FEValuesBase< dim, spacedim > &fe_values)
Definition: fe_values.cc:2641
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:4742
unsigned int present_face_index
Definition: fe_values.h:3625
void get_function_hessians(const InputVector &fe_function, std::vector< typename ProductType< hessian_type, typename InputVector::value_type >::type > &hessians) const
Definition: fe_values.cc:1732
void get_function_third_derivatives(const InputVector &fe_function, std::vector< Tensor< 3, spacedim, typename InputVector::value_type >> &third_derivatives) const
Definition: fe_values.cc:4164
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h: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:1875
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:2042
static::ExceptionBase & ExcFaceHasNoSubfaces()
const DerivativeForm< 3, dim, spacedim > & jacobian_2nd_derivative(const unsigned int quadrature_point) const
const unsigned int dofs_per_cell
Definition: fe_values.h: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:3630
FEValuesBase(const unsigned int n_q_points, const unsigned int dofs_per_cell, const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe)
Definition: fe_values.cc:3161
const Tensor< 1, spacedim > & normal_vector(const unsigned int i) const
void get_function_gradients(const InputVector &fe_function, std::vector< Tensor< 1, spacedim, typename InputVector::value_type >> &gradients) const
Definition: fe_values.cc:3799
Volume element.
const Tensor< 3, spacedim > & jacobian_pushed_forward_grad(const unsigned int quadrature_point) const
static const unsigned int dimension
Definition: fe_values.h:3657
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:3261
::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()
Scalar & operator=(const Scalar< dim, spacedim > &)
Definition: fe_values.cc:188
Transformed quadrature points.
void get_function_curls_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::curl_type > &curls) const
Definition: fe_values.cc:2154
void do_reinit(const unsigned int face_no)
Definition: fe_values.cc:4897
const Triangulation< dim, spacedim >::cell_iterator get_cell() const
Definition: fe_values.cc:4279
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:1930
void check_cell_similarity(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
Definition: fe_values.cc:4417
static const unsigned int space_dimension
Definition: fe_values.h:2020
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:3313
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
void get_function_laplacians_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::laplacian_type > &laplacians) const
Definition: fe_values.cc:1819
void get_function_third_derivatives(const InputVector &fe_function, std::vector< typename ProductType< third_derivative_type, typename InputVector::value_type >::type > &third_derivatives) const
Definition: fe_values.cc:1844
static::ExceptionBase & ExcShapeFunctionNotPrimitive(int arg1)
void get_function_symmetric_gradients(const InputVector &fe_function, std::vector< typename ProductType< symmetric_gradient_type, typename InputVector::value_type >::type > &symmetric_gradients) const
Definition: fe_values.cc:2011
void get_function_gradients(const InputVector &fe_function, std::vector< typename ProductType< gradient_type, typename InputVector::value_type >::type > &gradients) const
Definition: fe_values.cc:1955
const std::vector< double > & get_JxW_values() const
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition: fe_values.h:3405
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
void do_reinit(const unsigned int face_no, const unsigned int subface_no)
Definition: fe_values.cc:5096
UpdateFlags get_update_flags() const
static const unsigned int integral_dimension
Definition: fe_values.h:3782
void get_function_laplacians_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::laplacian_type > &laplacians) const
Definition: fe_values.cc:2271
const Quadrature< dim > & get_quadrature() const
UpdateFlags compute_update_flags(const UpdateFlags update_flags) const
Definition: fe_values.cc:4333
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:3665
void get_function_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::gradient_type > &gradients) const
Definition: fe_values.cc:1986
const std::vector< DerivativeForm< 3, dim, spacedim > > & get_jacobian_2nd_derivatives() const
void get_function_laplacians(const InputVector &fe_function, std::vector< typename InputVector::value_type > &laplacians) const
Definition: fe_values.cc:4027
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:2179
symmetric_gradient_type symmetric_gradient(const unsigned int shape_function, const unsigned int q_point) const
const std::vector< Tensor< 5, spacedim > > & get_jacobian_pushed_forward_3rd_derivatives() const
const Quadrature< dim-1 > & get_quadrature() const
divergence_type divergence(const unsigned int shape_function, const unsigned int q_point) const
Number value_type
Definition: vector.h:123
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:4730
#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:1676
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:1763
#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:2299
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
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:3523
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:4313
#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:1707
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
void invalidate_present_cell()
Definition: fe_values.cc:4352
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition: fe_values.h:3321
typename Tensor< rank_-1, dim, Number >::tensor_type value_type
Definition: tensor.h:427
static const unsigned int integral_dimension
Definition: fe_values.h:3441
Vector & operator=(const Vector< dim, spacedim > &)
Definition: fe_values.cc:280
SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
Tensor()
Definition: tensor.h:1015
const std::vector< Tensor< 1, spacedim > > & get_normal_vectors() const
Definition: fe_values.cc:4300
void get_function_divergences_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::divergence_type > &divergences) const
Definition: fe_values.cc:2098
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
bool is_nonzero_shape_function_component[spacedim]
Definition: fe_values.h:727
const FEValues< dim, spacedim > & get_present_fe_values() const
static const unsigned int integral_dimension
Definition: fe_values.h:3560
Second derivatives of shape functions.
static::ExceptionBase & ExcFENotPrimitive()
Gradient of volume element.
FEValuesBase & operator=(const FEValuesBase &)
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:1399
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:4638
void get_function_values(const InputVector &fe_function, std::vector< typename ProductType< value_type, typename InputVector::value_type >::type > &values) const
Definition: fe_values.cc:1899
::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:3286
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:3776
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:4710
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:4529
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:4370
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:2066
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:4974
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:4288
value_type value(const unsigned int shape_function, const unsigned int q_point) const
void get_function_laplacians(const InputVector &fe_function, std::vector< typename ProductType< value_type, typename InputVector::value_type >::type > &laplacians) const
Definition: fe_values.cc:2235
static::ExceptionBase & ExcNotImplemented()
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access >> &cell, const unsigned int face_no)
Definition: fe_values.cc:4845
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1214
boost::signals2::connection tria_listener_refinement
Definition: fe_values.h:3277
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3353
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:4798
::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:3328
FEValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4495
const Mapping< dim, spacedim > & get_mapping() const
void do_reinit()
Definition: fe_values.cc:4665
FESubfaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim-1 > &face_quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4938
Tensor< 1, spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
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:1651
FEFaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim-1 > &quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4762
CellSimilarity::Similarity get_cell_similarity() const
Definition: fe_values.cc:4472
void get_function_third_derivatives_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::third_derivative_type > &third_derivatives) const
Definition: fe_values.cc:2330
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:3359
std::size_t memory_consumption() const
Definition: fe_values.cc:4699
static::ExceptionBase & ExcInternalError()
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
Definition: fe_values.h:3337
const Tensor< 1, spacedim > & boundary_form(const unsigned int i) const
Tensor< 2, spacedim > shape_hessian_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
void get_function_hessians(const InputVector &fe_function, std::vector< Tensor< 2, spacedim, typename InputVector::value_type >> &hessians) const
Definition: fe_values.cc:3913