Reference documentation for deal.II version Git 39b7d3efb0 2021-05-07 15:49:09 -0400
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
fe_values.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2020 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 
24 #include <deal.II/base/point.h>
29 
32 
33 #include <deal.II/fe/fe.h>
36 #include <deal.II/fe/mapping.h>
37 
38 #include <deal.II/grid/tria.h>
40 
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 
56 
57 // Forward declaration
58 #ifndef DOXYGEN
59 template <int dim, int spacedim = dim>
60 class FEValuesBase;
61 #endif
62 
63 namespace internal
64 {
69  template <int dim, class NumberType = double>
70  struct CurlType;
71 
78  template <class NumberType>
79  struct CurlType<1, NumberType>
80  {
82  };
83 
90  template <class NumberType>
91  struct CurlType<2, NumberType>
92  {
94  };
95 
102  template <class NumberType>
103  struct CurlType<3, NumberType>
104  {
106  };
107 } // namespace internal
108 
109 
110 
132 namespace FEValuesViews
133 {
145  template <int dim, int spacedim = dim>
146  class Scalar
147  {
148  public:
155 
162 
169 
176 
183  template <typename Number>
185 
192  template <typename Number>
193  using solution_gradient_type =
195 
202  template <typename Number>
205 
212  template <typename Number>
213  using solution_hessian_type =
215 
222  template <typename Number>
225 
232  template <typename Number>
234  {
239  using value_type =
240  typename ProductType<Number,
242 
247  using gradient_type = typename ProductType<
248  Number,
250 
255  using laplacian_type =
256  typename ProductType<Number,
258 
263  using hessian_type = typename ProductType<
264  Number,
266 
271  using third_derivative_type = typename ProductType<
272  Number,
274  };
275 
281  {
291 
300  unsigned int row_index;
301  };
302 
306  Scalar();
307 
313  Scalar(const FEValuesBase<dim, spacedim> &fe_values_base,
314  const unsigned int component);
315 
320  Scalar &
321  operator=(const Scalar<dim, spacedim> &) = delete;
322 
336  value_type
337  value(const unsigned int shape_function, const unsigned int q_point) const;
338 
350  gradient(const unsigned int shape_function,
351  const unsigned int q_point) const;
352 
364  hessian(const unsigned int shape_function,
365  const unsigned int q_point) const;
366 
378  third_derivative(const unsigned int shape_function,
379  const unsigned int q_point) const;
380 
398  template <class InputVector>
399  void
400  get_function_values(
401  const InputVector &fe_function,
403  &values) const;
404 
439  template <class InputVector>
440  void
441  get_function_values_from_local_dof_values(
442  const InputVector &dof_values,
444  &values) const;
445 
463  template <class InputVector>
464  void
465  get_function_gradients(
466  const InputVector &fe_function,
468  &gradients) const;
469 
476  template <class InputVector>
477  void
478  get_function_gradients_from_local_dof_values(
479  const InputVector &dof_values,
481  &gradients) const;
482 
500  template <class InputVector>
501  void
502  get_function_hessians(
503  const InputVector &fe_function,
505  &hessians) const;
506 
513  template <class InputVector>
514  void
515  get_function_hessians_from_local_dof_values(
516  const InputVector &dof_values,
518  &hessians) const;
519 
520 
539  template <class InputVector>
540  void
541  get_function_laplacians(
542  const InputVector &fe_function,
544  &laplacians) const;
545 
552  template <class InputVector>
553  void
554  get_function_laplacians_from_local_dof_values(
555  const InputVector &dof_values,
557  &laplacians) const;
558 
559 
578  template <class InputVector>
579  void
580  get_function_third_derivatives(
581  const InputVector &fe_function,
582  std::vector<
584  &third_derivatives) const;
585 
592  template <class InputVector>
593  void
594  get_function_third_derivatives_from_local_dof_values(
595  const InputVector &dof_values,
596  std::vector<
598  &third_derivatives) const;
599 
600 
601  private:
606 
611  const unsigned int component;
612 
616  std::vector<ShapeFunctionData> shape_function_data;
617  };
618 
619 
620 
650  template <int dim, int spacedim = dim>
651  class Vector
652  {
653  public:
660 
670 
682 
689 
696  using curl_type = typename ::internal::CurlType<spacedim>::type;
697 
704 
711 
718  template <typename Number>
720 
727  template <typename Number>
728  using solution_gradient_type =
730 
737  template <typename Number>
740 
747  template <typename Number>
750 
757  template <typename Number>
760 
767  template <typename Number>
769 
776  template <typename Number>
777  using solution_hessian_type =
779 
786  template <typename Number>
789 
796  template <typename Number>
798  {
803  using value_type =
804  typename ProductType<Number,
806 
811  using gradient_type = typename ProductType<
812  Number,
814 
819  using symmetric_gradient_type = typename ProductType<
820  Number,
822 
827  using divergence_type = typename ProductType<
828  Number,
830 
835  using laplacian_type =
836  typename ProductType<Number,
838 
843  using curl_type =
844  typename ProductType<Number,
846 
851  using hessian_type = typename ProductType<
852  Number,
854 
859  using third_derivative_type = typename ProductType<
860  Number,
862  };
863 
869  {
878  bool is_nonzero_shape_function_component[spacedim];
879 
889  unsigned int row_index[spacedim];
890 
901  };
902 
906  Vector();
907 
916  Vector(const FEValuesBase<dim, spacedim> &fe_values_base,
917  const unsigned int first_vector_component);
918 
923  Vector &
924  operator=(const Vector<dim, spacedim> &) = delete;
925 
942  value_type
943  value(const unsigned int shape_function, const unsigned int q_point) const;
944 
959  gradient(const unsigned int shape_function,
960  const unsigned int q_point) const;
961 
978  symmetric_gradient(const unsigned int shape_function,
979  const unsigned int q_point) const;
980 
992  divergence(const unsigned int shape_function,
993  const unsigned int q_point) const;
994 
1015  curl_type
1016  curl(const unsigned int shape_function, const unsigned int q_point) const;
1017 
1028  hessian_type
1029  hessian(const unsigned int shape_function,
1030  const unsigned int q_point) const;
1031 
1043  third_derivative(const unsigned int shape_function,
1044  const unsigned int q_point) const;
1045 
1063  template <class InputVector>
1064  void
1065  get_function_values(
1066  const InputVector &fe_function,
1068  &values) const;
1069 
1104  template <class InputVector>
1105  void
1106  get_function_values_from_local_dof_values(
1107  const InputVector &dof_values,
1109  &values) const;
1110 
1128  template <class InputVector>
1129  void
1130  get_function_gradients(
1131  const InputVector &fe_function,
1133  &gradients) const;
1134 
1141  template <class InputVector>
1142  void
1143  get_function_gradients_from_local_dof_values(
1144  const InputVector &dof_values,
1146  &gradients) const;
1147 
1171  template <class InputVector>
1172  void
1173  get_function_symmetric_gradients(
1174  const InputVector &fe_function,
1175  std::vector<
1177  &symmetric_gradients) const;
1178 
1185  template <class InputVector>
1186  void
1187  get_function_symmetric_gradients_from_local_dof_values(
1188  const InputVector &dof_values,
1189  std::vector<
1191  &symmetric_gradients) const;
1192 
1211  template <class InputVector>
1212  void
1213  get_function_divergences(
1214  const InputVector &fe_function,
1216  &divergences) const;
1217 
1224  template <class InputVector>
1225  void
1226  get_function_divergences_from_local_dof_values(
1227  const InputVector &dof_values,
1229  &divergences) const;
1230 
1249  template <class InputVector>
1250  void
1251  get_function_curls(
1252  const InputVector &fe_function,
1254  const;
1255 
1262  template <class InputVector>
1263  void
1264  get_function_curls_from_local_dof_values(
1265  const InputVector &dof_values,
1267  const;
1268 
1286  template <class InputVector>
1287  void
1288  get_function_hessians(
1289  const InputVector &fe_function,
1291  &hessians) const;
1292 
1299  template <class InputVector>
1300  void
1301  get_function_hessians_from_local_dof_values(
1302  const InputVector &dof_values,
1304  &hessians) const;
1305 
1324  template <class InputVector>
1325  void
1326  get_function_laplacians(
1327  const InputVector &fe_function,
1329  &laplacians) const;
1330 
1337  template <class InputVector>
1338  void
1339  get_function_laplacians_from_local_dof_values(
1340  const InputVector &dof_values,
1342  &laplacians) const;
1343 
1362  template <class InputVector>
1363  void
1364  get_function_third_derivatives(
1365  const InputVector &fe_function,
1366  std::vector<
1368  &third_derivatives) const;
1369 
1376  template <class InputVector>
1377  void
1378  get_function_third_derivatives_from_local_dof_values(
1379  const InputVector &dof_values,
1380  std::vector<
1382  &third_derivatives) const;
1383 
1384  private:
1389 
1394  const unsigned int first_vector_component;
1395 
1399  std::vector<ShapeFunctionData> shape_function_data;
1400  };
1401 
1402 
1403  template <int rank, int dim, int spacedim = dim>
1405 
1428  template <int dim, int spacedim>
1429  class SymmetricTensor<2, dim, spacedim>
1430  {
1431  public:
1439 
1450 
1457  template <typename Number>
1459 
1466  template <typename Number>
1467  using solution_divergence_type =
1469 
1470 
1477  template <typename Number>
1478  struct DEAL_II_DEPRECATED_EARLY OutputType
1479  {
1484  using value_type = typename ProductType<
1485  Number,
1487 
1492  using divergence_type = typename ProductType<
1493  Number,
1495  };
1496 
1501  struct ShapeFunctionData
1502  {
1511  bool is_nonzero_shape_function_component
1512  [value_type::n_independent_components];
1513 
1523  unsigned int row_index[value_type::n_independent_components];
1524 
1534 
1539  };
1540 
1544  SymmetricTensor();
1545 
1555  SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1556  const unsigned int first_tensor_component);
1557 
1562  SymmetricTensor &
1563  operator=(const SymmetricTensor<2, dim, spacedim> &) = delete;
1564 
1582  value_type
1583  value(const unsigned int shape_function, const unsigned int q_point) const;
1584 
1599  divergence(const unsigned int shape_function,
1600  const unsigned int q_point) const;
1601 
1619  template <class InputVector>
1620  void
1621  get_function_values(
1622  const InputVector &fe_function,
1624  &values) const;
1625 
1660  template <class InputVector>
1661  void
1662  get_function_values_from_local_dof_values(
1663  const InputVector &dof_values,
1665  &values) const;
1666 
1688  template <class InputVector>
1689  void
1690  get_function_divergences(
1691  const InputVector &fe_function,
1693  &divergences) const;
1694 
1701  template <class InputVector>
1702  void
1703  get_function_divergences_from_local_dof_values(
1704  const InputVector &dof_values,
1706  &divergences) const;
1707 
1708  private:
1713 
1718  const unsigned int first_tensor_component;
1719 
1723  std::vector<ShapeFunctionData> shape_function_data;
1724  };
1725 
1726 
1727  template <int rank, int dim, int spacedim = dim>
1728  class Tensor;
1729 
1748  template <int dim, int spacedim>
1749  class Tensor<2, dim, spacedim>
1750  {
1751  public:
1757 
1762 
1768 
1775  template <typename Number>
1777 
1784  template <typename Number>
1785  using solution_divergence_type =
1787 
1794  template <typename Number>
1795  using solution_gradient_type =
1797 
1798 
1805  template <typename Number>
1806  struct DEAL_II_DEPRECATED_EARLY OutputType
1807  {
1812  using value_type = typename ProductType<
1813  Number,
1815 
1820  using divergence_type = typename ProductType<
1821  Number,
1823 
1828  using gradient_type = typename ProductType<
1829  Number,
1831  };
1832 
1837  struct ShapeFunctionData
1838  {
1847  bool is_nonzero_shape_function_component
1848  [value_type::n_independent_components];
1849 
1859  unsigned int row_index[value_type::n_independent_components];
1860 
1870 
1875  };
1876 
1880  Tensor();
1881 
1891  Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1892  const unsigned int first_tensor_component);
1893 
1894 
1899  Tensor &
1900  operator=(const Tensor<2, dim, spacedim> &) = delete;
1901 
1918  value_type
1919  value(const unsigned int shape_function, const unsigned int q_point) const;
1920 
1935  divergence(const unsigned int shape_function,
1936  const unsigned int q_point) const;
1937 
1952  gradient(const unsigned int shape_function,
1953  const unsigned int q_point) const;
1954 
1972  template <class InputVector>
1973  void
1974  get_function_values(
1975  const InputVector &fe_function,
1977  &values) const;
1978 
2013  template <class InputVector>
2014  void
2015  get_function_values_from_local_dof_values(
2016  const InputVector &dof_values,
2018  &values) const;
2019 
2041  template <class InputVector>
2042  void
2043  get_function_divergences(
2044  const InputVector &fe_function,
2046  &divergences) const;
2047 
2054  template <class InputVector>
2055  void
2056  get_function_divergences_from_local_dof_values(
2057  const InputVector &dof_values,
2059  &divergences) const;
2060 
2077  template <class InputVector>
2078  void
2079  get_function_gradients(
2080  const InputVector &fe_function,
2082  &gradients) const;
2083 
2090  template <class InputVector>
2091  void
2092  get_function_gradients_from_local_dof_values(
2093  const InputVector &dof_values,
2095  &gradients) const;
2096 
2097  private:
2102 
2107  const unsigned int first_tensor_component;
2108 
2112  std::vector<ShapeFunctionData> shape_function_data;
2113  };
2114 
2115 } // namespace FEValuesViews
2116 
2117 
2118 namespace internal
2119 {
2120  namespace FEValuesViews
2121  {
2126  template <int dim, int spacedim, typename Extractor>
2127  struct ViewType
2128  {};
2129 
2137  template <int dim, int spacedim>
2138  struct ViewType<dim, spacedim, FEValuesExtractors::Scalar>
2139  {
2140  using type = typename ::FEValuesViews::Scalar<dim, spacedim>;
2141  };
2142 
2150  template <int dim, int spacedim>
2151  struct ViewType<dim, spacedim, FEValuesExtractors::Vector>
2152  {
2153  using type = typename ::FEValuesViews::Vector<dim, spacedim>;
2154  };
2155 
2163  template <int dim, int spacedim, int rank>
2164  struct ViewType<dim, spacedim, FEValuesExtractors::Tensor<rank>>
2165  {
2166  using type = typename ::FEValuesViews::Tensor<rank, dim, spacedim>;
2167  };
2168 
2176  template <int dim, int spacedim, int rank>
2177  struct ViewType<dim, spacedim, FEValuesExtractors::SymmetricTensor<rank>>
2178  {
2179  using type =
2180  typename ::FEValuesViews::SymmetricTensor<rank, dim, spacedim>;
2181  };
2182 
2190  template <int dim, int spacedim>
2191  struct Cache
2192  {
2197  std::vector<::FEValuesViews::Scalar<dim, spacedim>> scalars;
2198  std::vector<::FEValuesViews::Vector<dim, spacedim>> vectors;
2199  std::vector<::FEValuesViews::SymmetricTensor<2, dim, spacedim>>
2201  std::vector<::FEValuesViews::Tensor<2, dim, spacedim>>
2203 
2207  Cache(const FEValuesBase<dim, spacedim> &fe_values);
2208  };
2209  } // namespace FEValuesViews
2210 } // namespace internal
2211 
2212 namespace FEValuesViews
2213 {
2218  template <int dim, int spacedim, typename Extractor>
2219  using View = typename ::internal::FEValuesViews::
2220  ViewType<dim, spacedim, Extractor>::type;
2221 } // namespace FEValuesViews
2222 
2223 
2323 template <int dim, int spacedim>
2324 class FEValuesBase : public Subscriptor
2325 {
2326 public:
2330  static const unsigned int dimension = dim;
2331 
2335  static const unsigned int space_dimension = spacedim;
2336 
2344  const unsigned int n_quadrature_points;
2345 
2355  const unsigned int max_n_quadrature_points;
2356 
2362  const unsigned int dofs_per_cell;
2363 
2364 
2372  FEValuesBase(const unsigned int n_q_points,
2373  const unsigned int dofs_per_cell,
2374  const UpdateFlags update_flags,
2375  const Mapping<dim, spacedim> & mapping,
2376  const FiniteElement<dim, spacedim> &fe);
2377 
2382  FEValuesBase &
2383  operator=(const FEValuesBase &) = delete;
2384 
2389  FEValuesBase(const FEValuesBase &) = delete;
2390 
2394  virtual ~FEValuesBase() override;
2395 
2396 
2400 
2401 
2422  const double &
2423  shape_value(const unsigned int function_no,
2424  const unsigned int point_no) const;
2425 
2446  double
2447  shape_value_component(const unsigned int function_no,
2448  const unsigned int point_no,
2449  const unsigned int component) const;
2450 
2476  const Tensor<1, spacedim> &
2477  shape_grad(const unsigned int function_no,
2478  const unsigned int quadrature_point) const;
2479 
2497  shape_grad_component(const unsigned int function_no,
2498  const unsigned int point_no,
2499  const unsigned int component) const;
2500 
2520  const Tensor<2, spacedim> &
2521  shape_hessian(const unsigned int function_no,
2522  const unsigned int point_no) const;
2523 
2541  shape_hessian_component(const unsigned int function_no,
2542  const unsigned int point_no,
2543  const unsigned int component) const;
2544 
2564  const Tensor<3, spacedim> &
2565  shape_3rd_derivative(const unsigned int function_no,
2566  const unsigned int point_no) const;
2567 
2585  shape_3rd_derivative_component(const unsigned int function_no,
2586  const unsigned int point_no,
2587  const unsigned int component) const;
2588 
2590 
2592 
2629  template <class InputVector>
2630  void
2631  get_function_values(
2632  const InputVector & fe_function,
2633  std::vector<typename InputVector::value_type> &values) const;
2634 
2648  template <class InputVector>
2649  void
2650  get_function_values(
2651  const InputVector & fe_function,
2652  std::vector<Vector<typename InputVector::value_type>> &values) const;
2653 
2710  template <class InputVector>
2711  void
2712  get_function_values(
2713  const InputVector & fe_function,
2715  std::vector<typename InputVector::value_type> & values) const;
2716 
2725  template <class InputVector>
2726  void
2727  get_function_values(
2728  const InputVector & fe_function,
2730  std::vector<Vector<typename InputVector::value_type>> &values) const;
2731 
2732 
2754  template <class InputVector>
2755  void
2756  get_function_values(
2757  const InputVector & fe_function,
2759  ArrayView<std::vector<typename InputVector::value_type>> values,
2760  const bool quadrature_points_fastest) const;
2761 
2763 
2765 
2802  template <class InputVector>
2803  void
2804  get_function_gradients(
2805  const InputVector &fe_function,
2807  &gradients) const;
2808 
2825  template <class InputVector>
2826  void
2827  get_function_gradients(
2828  const InputVector &fe_function,
2829  std::vector<
2831  &gradients) const;
2832 
2841  template <class InputVector>
2842  void
2843  get_function_gradients(
2844  const InputVector & fe_function,
2847  &gradients) const;
2848 
2857  template <class InputVector>
2858  void
2859  get_function_gradients(
2860  const InputVector & fe_function,
2862  ArrayView<
2864  gradients,
2865  const bool quadrature_points_fastest = false) const;
2866 
2868 
2872 
2910  template <class InputVector>
2911  void
2912  get_function_hessians(
2913  const InputVector &fe_function,
2915  &hessians) const;
2916 
2934  template <class InputVector>
2935  void
2936  get_function_hessians(
2937  const InputVector &fe_function,
2938  std::vector<
2940  & hessians,
2941  const bool quadrature_points_fastest = false) const;
2942 
2951  template <class InputVector>
2952  void
2953  get_function_hessians(
2954  const InputVector & fe_function,
2957  &hessians) const;
2958 
2967  template <class InputVector>
2968  void
2969  get_function_hessians(
2970  const InputVector & fe_function,
2972  ArrayView<
2974  hessians,
2975  const bool quadrature_points_fastest = false) const;
2976 
3017  template <class InputVector>
3018  void
3019  get_function_laplacians(
3020  const InputVector & fe_function,
3021  std::vector<typename InputVector::value_type> &laplacians) const;
3022 
3042  template <class InputVector>
3043  void
3044  get_function_laplacians(
3045  const InputVector & fe_function,
3046  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
3047 
3056  template <class InputVector>
3057  void
3058  get_function_laplacians(
3059  const InputVector & fe_function,
3061  std::vector<typename InputVector::value_type> & laplacians) const;
3062 
3071  template <class InputVector>
3072  void
3073  get_function_laplacians(
3074  const InputVector & fe_function,
3076  std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
3077 
3086  template <class InputVector>
3087  void
3088  get_function_laplacians(
3089  const InputVector & fe_function,
3091  std::vector<std::vector<typename InputVector::value_type>> &laplacians,
3092  const bool quadrature_points_fastest = false) const;
3093 
3095 
3097 
3136  template <class InputVector>
3137  void
3138  get_function_third_derivatives(
3139  const InputVector &fe_function,
3141  &third_derivatives) const;
3142 
3161  template <class InputVector>
3162  void
3163  get_function_third_derivatives(
3164  const InputVector &fe_function,
3165  std::vector<
3167  & third_derivatives,
3168  const bool quadrature_points_fastest = false) const;
3169 
3178  template <class InputVector>
3179  void
3180  get_function_third_derivatives(
3181  const InputVector & fe_function,
3184  &third_derivatives) const;
3185 
3194  template <class InputVector>
3195  void
3196  get_function_third_derivatives(
3197  const InputVector & fe_function,
3199  ArrayView<
3201  third_derivatives,
3202  const bool quadrature_points_fastest = false) const;
3204 
3206 
3207 
3232  dof_indices() const;
3233 
3267  dof_indices_starting_at(const unsigned int start_dof_index) const;
3268 
3300  dof_indices_ending_at(const unsigned int end_dof_index) const;
3301 
3303 
3305 
3306 
3328  quadrature_point_indices() const;
3329 
3335  const Point<spacedim> &
3336  quadrature_point(const unsigned int q) const;
3337 
3343  const std::vector<Point<spacedim>> &
3344  get_quadrature_points() const;
3345 
3361  double
3362  JxW(const unsigned int quadrature_point) const;
3363 
3367  const std::vector<double> &
3368  get_JxW_values() const;
3369 
3377  jacobian(const unsigned int quadrature_point) const;
3378 
3385  const std::vector<DerivativeForm<1, dim, spacedim>> &
3386  get_jacobians() const;
3387 
3396  jacobian_grad(const unsigned int quadrature_point) const;
3397 
3404  const std::vector<DerivativeForm<2, dim, spacedim>> &
3405  get_jacobian_grads() const;
3406 
3415  const Tensor<3, spacedim> &
3416  jacobian_pushed_forward_grad(const unsigned int quadrature_point) const;
3417 
3424  const std::vector<Tensor<3, spacedim>> &
3425  get_jacobian_pushed_forward_grads() const;
3426 
3435  jacobian_2nd_derivative(const unsigned int quadrature_point) const;
3436 
3443  const std::vector<DerivativeForm<3, dim, spacedim>> &
3444  get_jacobian_2nd_derivatives() const;
3445 
3455  const Tensor<4, spacedim> &
3456  jacobian_pushed_forward_2nd_derivative(
3457  const unsigned int quadrature_point) const;
3458 
3465  const std::vector<Tensor<4, spacedim>> &
3466  get_jacobian_pushed_forward_2nd_derivatives() const;
3467 
3477  jacobian_3rd_derivative(const unsigned int quadrature_point) const;
3478 
3485  const std::vector<DerivativeForm<4, dim, spacedim>> &
3486  get_jacobian_3rd_derivatives() const;
3487 
3497  const Tensor<5, spacedim> &
3498  jacobian_pushed_forward_3rd_derivative(
3499  const unsigned int quadrature_point) const;
3500 
3507  const std::vector<Tensor<5, spacedim>> &
3508  get_jacobian_pushed_forward_3rd_derivatives() const;
3509 
3517  inverse_jacobian(const unsigned int quadrature_point) const;
3518 
3525  const std::vector<DerivativeForm<1, spacedim, dim>> &
3526  get_inverse_jacobians() const;
3527 
3547  const Tensor<1, spacedim> &
3548  normal_vector(const unsigned int i) const;
3549 
3557  const std::vector<Tensor<1, spacedim>> &
3558  get_normal_vectors() const;
3559 
3561 
3563 
3564 
3574  operator[](const FEValuesExtractors::Scalar &scalar) const;
3575 
3585  operator[](const FEValuesExtractors::Vector &vector) const;
3586 
3597  operator[](const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
3598 
3599 
3609  operator[](const FEValuesExtractors::Tensor<2> &tensor) const;
3610 
3612 
3614 
3615 
3619  const Mapping<dim, spacedim> &
3620  get_mapping() const;
3621 
3626  get_fe() const;
3627 
3631  UpdateFlags
3632  get_update_flags() const;
3633 
3638  get_cell() const;
3639 
3646  get_cell_similarity() const;
3647 
3652  std::size_t
3653  memory_consumption() const;
3655 
3656 
3665  std::string,
3666  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3667  << "object for which this kind of information has not been computed. What "
3668  << "information these objects compute is determined by the update_* flags you "
3669  << "pass to the constructor. Here, the operation you are attempting requires "
3670  << "the <" << arg1
3671  << "> flag to be set, but it was apparently not specified "
3672  << "upon construction.");
3673 
3681  ExcFEDontMatch,
3682  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3683  "to the DoFHandler that provided the cell iterator do not match.");
3689  DeclException1(ExcShapeFunctionNotPrimitive,
3690  int,
3691  << "The shape function with index " << arg1
3692  << " is not primitive, i.e. it is vector-valued and "
3693  << "has more than one non-zero vector component. This "
3694  << "function cannot be called for these shape functions. "
3695  << "Maybe you want to use the same function with the "
3696  << "_component suffix?");
3697 
3706  "The given FiniteElement is not a primitive element but the requested operation "
3707  "only works for those. See FiniteElement::is_primitive() for more information.");
3708 
3709 protected:
3740  class CellIteratorBase;
3741 
3746  template <typename CI>
3749 
3755  std::unique_ptr<const CellIteratorBase> present_cell;
3756 
3764  boost::signals2::connection tria_listener_refinement;
3765 
3773  boost::signals2::connection tria_listener_mesh_transform;
3774 
3780  void
3781  invalidate_present_cell();
3782 
3792  void
3793  maybe_invalidate_previous_present_cell(
3794  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3795 
3801 
3807  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3809 
3816 
3817 
3825 
3831  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3833 
3839  spacedim>
3841 
3842 
3847 
3856  UpdateFlags
3857  compute_update_flags(const UpdateFlags update_flags) const;
3858 
3865 
3871  void
3872  check_cell_similarity(
3873  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3874 
3875 private:
3880 
3881  // Make the view classes friends of this class, since they access internal
3882  // data.
3883  template <int, int>
3885  template <int, int>
3887  template <int, int, int>
3889  template <int, int, int>
3891 };
3892 
3893 
3894 
3904 template <int dim, int spacedim = dim>
3905 class FEValues : public FEValuesBase<dim, spacedim>
3906 {
3907 public:
3912  static const unsigned int integral_dimension = dim;
3913 
3918  FEValues(const Mapping<dim, spacedim> & mapping,
3919  const FiniteElement<dim, spacedim> &fe,
3920  const Quadrature<dim> & quadrature,
3921  const UpdateFlags update_flags);
3922 
3929  FEValues(const Mapping<dim, spacedim> & mapping,
3930  const FiniteElement<dim, spacedim> &fe,
3931  const hp::QCollection<dim> & quadrature,
3932  const UpdateFlags update_flags);
3933 
3940  const Quadrature<dim> & quadrature,
3941  const UpdateFlags update_flags);
3942 
3950  const hp::QCollection<dim> & quadrature,
3951  const UpdateFlags update_flags);
3952 
3959  template <bool level_dof_access>
3960  void
3961  reinit(
3963 
3977  void
3979 
3984  const Quadrature<dim> &
3985  get_quadrature() const;
3986 
3991  std::size_t
3992  memory_consumption() const;
3993 
4008  const FEValues<dim, spacedim> &
4009  get_present_fe_values() const;
4010 
4011 private:
4016 
4020  void
4021  initialize(const UpdateFlags update_flags);
4022 
4029  void
4030  do_reinit();
4031 };
4032 
4033 
4043 template <int dim, int spacedim = dim>
4044 class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
4045 {
4046 public:
4051  static const unsigned int integral_dimension = dim - 1;
4052 
4064  FEFaceValuesBase(const unsigned int dofs_per_cell,
4065  const UpdateFlags update_flags,
4066  const Mapping<dim, spacedim> & mapping,
4067  const FiniteElement<dim, spacedim> &fe,
4068  const Quadrature<dim - 1> & quadrature);
4069 
4076  FEFaceValuesBase(const unsigned int dofs_per_cell,
4077  const UpdateFlags update_flags,
4078  const Mapping<dim, spacedim> & mapping,
4079  const FiniteElement<dim, spacedim> &fe,
4080  const hp::QCollection<dim - 1> & quadrature);
4081 
4089  const Tensor<1, spacedim> &
4090  boundary_form(const unsigned int i) const;
4091 
4098  const std::vector<Tensor<1, spacedim>> &
4099  get_boundary_forms() const;
4100 
4105  unsigned int
4106  get_face_index() const;
4107 
4112  const Quadrature<dim - 1> &
4113  get_quadrature() const;
4114 
4119  std::size_t
4120  memory_consumption() const;
4121 
4122 protected:
4127  unsigned int present_face_no;
4128 
4133  unsigned int present_face_index;
4134 
4138  const hp::QCollection<dim - 1> quadrature;
4139 };
4140 
4141 
4142 
4156 template <int dim, int spacedim = dim>
4157 class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
4158 {
4159 public:
4164  static const unsigned int dimension = dim;
4165 
4166  static const unsigned int space_dimension = spacedim;
4167 
4172  static const unsigned int integral_dimension = dim - 1;
4173 
4178  FEFaceValues(const Mapping<dim, spacedim> & mapping,
4179  const FiniteElement<dim, spacedim> &fe,
4180  const Quadrature<dim - 1> & quadrature,
4181  const UpdateFlags update_flags);
4182 
4189  FEFaceValues(const Mapping<dim, spacedim> & mapping,
4190  const FiniteElement<dim, spacedim> &fe,
4191  const hp::QCollection<dim - 1> & quadrature,
4192  const UpdateFlags update_flags);
4193 
4200  const Quadrature<dim - 1> & quadrature,
4201  const UpdateFlags update_flags);
4202 
4210  const hp::QCollection<dim - 1> & quadrature,
4211  const UpdateFlags update_flags);
4212 
4217  template <bool level_dof_access>
4218  void
4219  reinit(
4221  const unsigned int face_no);
4222 
4229  template <bool level_dof_access>
4230  void
4231  reinit(
4233  const typename Triangulation<dim, spacedim>::face_iterator & face);
4234 
4248  void
4250  const unsigned int face_no);
4251 
4252  /*
4253  * Reinitialize the gradients, Jacobi determinants, etc for the given face
4254  * on a given cell of type "iterator into a Triangulation object", and the
4255  * given finite element. Since iterators into a triangulation alone only
4256  * convey information about the geometry of a cell, but not about degrees of
4257  * freedom possibly associated with this cell, you will not be able to call
4258  * some functions of this class if they need information about degrees of
4259  * freedom. These functions are, above all, the
4260  * <tt>get_function_value/gradients/hessians/third_derivatives</tt>
4261  * functions. If you want to call these functions, you have to call the @p
4262  * reinit variants that take iterators into DoFHandler or other DoF handler
4263  * type objects.
4264  *
4265  * @note @p face must be one of @p cell's face iterators.
4266  */
4267  void
4269  const typename Triangulation<dim, spacedim>::face_iterator &face);
4270 
4286  get_present_fe_values() const;
4287 
4288 private:
4292  void
4293  initialize(const UpdateFlags update_flags);
4294 
4301  void
4302  do_reinit(const unsigned int face_no);
4303 };
4304 
4305 
4322 template <int dim, int spacedim = dim>
4323 class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
4324 {
4325 public:
4329  static const unsigned int dimension = dim;
4330 
4334  static const unsigned int space_dimension = spacedim;
4335 
4340  static const unsigned int integral_dimension = dim - 1;
4341 
4346  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
4347  const FiniteElement<dim, spacedim> &fe,
4348  const Quadrature<dim - 1> & face_quadrature,
4349  const UpdateFlags update_flags);
4350 
4357  FESubfaceValues(const Mapping<dim, spacedim> & mapping,
4358  const FiniteElement<dim, spacedim> &fe,
4359  const hp::QCollection<dim - 1> & face_quadrature,
4360  const UpdateFlags update_flags);
4361 
4368  const Quadrature<dim - 1> & face_quadrature,
4369  const UpdateFlags update_flags);
4370 
4378  const hp::QCollection<dim - 1> & face_quadrature,
4379  const UpdateFlags update_flags);
4380 
4387  template <bool level_dof_access>
4388  void
4389  reinit(
4391  const unsigned int face_no,
4392  const unsigned int subface_no);
4393 
4398  template <bool level_dof_access>
4399  void
4400  reinit(
4402  const typename Triangulation<dim, spacedim>::face_iterator & face,
4403  const typename Triangulation<dim, spacedim>::face_iterator &subface);
4404 
4418  void
4420  const unsigned int face_no,
4421  const unsigned int subface_no);
4422 
4442  void
4444  const typename Triangulation<dim, spacedim>::face_iterator &face,
4445  const typename Triangulation<dim, spacedim>::face_iterator &subface);
4446 
4462  get_present_fe_values() const;
4463 
4469  DeclException0(ExcReinitCalledWithBoundaryFace);
4470 
4476  DeclException0(ExcFaceHasNoSubfaces);
4477 
4478 private:
4482  void
4483  initialize(const UpdateFlags update_flags);
4484 
4491  void
4492  do_reinit(const unsigned int face_no, const unsigned int subface_no);
4493 };
4494 
4495 
4496 #ifndef DOXYGEN
4497 
4498 
4499 /*------------------------ Inline functions: namespace FEValuesViews --------*/
4500 
4501 namespace FEValuesViews
4502 {
4503  template <int dim, int spacedim>
4504  inline typename Scalar<dim, spacedim>::value_type
4505  Scalar<dim, spacedim>::value(const unsigned int shape_function,
4506  const unsigned int q_point) const
4507  {
4508  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4509  Assert(
4510  fe_values->update_flags & update_values,
4512  "update_values"))));
4513 
4514  // an adaptation of the FEValuesBase::shape_value_component function
4515  // except that here we know the component as fixed and we have
4516  // pre-computed and cached a bunch of information. See the comments there.
4517  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4518  return fe_values->finite_element_output.shape_values(
4519  shape_function_data[shape_function].row_index, q_point);
4520  else
4521  return 0;
4522  }
4523 
4524 
4525 
4526  template <int dim, int spacedim>
4527  inline typename Scalar<dim, spacedim>::gradient_type
4528  Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
4529  const unsigned int q_point) const
4530  {
4531  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4532  Assert(fe_values->update_flags & update_gradients,
4534  "update_gradients")));
4535 
4536  // an adaptation of the FEValuesBase::shape_grad_component
4537  // function except that here we know the component as fixed and we have
4538  // pre-computed and cached a bunch of information. See the comments there.
4539  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4540  return fe_values->finite_element_output
4541  .shape_gradients[shape_function_data[shape_function].row_index]
4542  [q_point];
4543  else
4544  return gradient_type();
4545  }
4546 
4547 
4548 
4549  template <int dim, int spacedim>
4550  inline typename Scalar<dim, spacedim>::hessian_type
4551  Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
4552  const unsigned int q_point) const
4553  {
4554  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4555  Assert(fe_values->update_flags & update_hessians,
4557  "update_hessians")));
4558 
4559  // an adaptation of the FEValuesBase::shape_hessian_component
4560  // function except that here we know the component as fixed and we have
4561  // pre-computed and cached a bunch of information. See the comments there.
4562  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4563  return fe_values->finite_element_output
4564  .shape_hessians[shape_function_data[shape_function].row_index][q_point];
4565  else
4566  return hessian_type();
4567  }
4568 
4569 
4570 
4571  template <int dim, int spacedim>
4572  inline typename Scalar<dim, spacedim>::third_derivative_type
4573  Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
4574  const unsigned int q_point) const
4575  {
4576  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4577  Assert(fe_values->update_flags & update_3rd_derivatives,
4579  "update_3rd_derivatives")));
4580 
4581  // an adaptation of the FEValuesBase::shape_3rdderivative_component
4582  // function except that here we know the component as fixed and we have
4583  // pre-computed and cached a bunch of information. See the comments there.
4584  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4585  return fe_values->finite_element_output
4586  .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
4587  [q_point];
4588  else
4589  return third_derivative_type();
4590  }
4591 
4592 
4593 
4594  template <int dim, int spacedim>
4595  inline typename Vector<dim, spacedim>::value_type
4596  Vector<dim, spacedim>::value(const unsigned int shape_function,
4597  const unsigned int q_point) const
4598  {
4599  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4600  Assert(fe_values->update_flags & update_values,
4602  "update_values")));
4603 
4604  // same as for the scalar case except that we have one more index
4605  const int snc =
4606  shape_function_data[shape_function].single_nonzero_component;
4607  if (snc == -2)
4608  return value_type();
4609  else if (snc != -1)
4610  {
4611  value_type return_value;
4612  return_value[shape_function_data[shape_function]
4613  .single_nonzero_component_index] =
4614  fe_values->finite_element_output.shape_values(snc, q_point);
4615  return return_value;
4616  }
4617  else
4618  {
4619  value_type return_value;
4620  for (unsigned int d = 0; d < dim; ++d)
4621  if (shape_function_data[shape_function]
4622  .is_nonzero_shape_function_component[d])
4623  return_value[d] = fe_values->finite_element_output.shape_values(
4624  shape_function_data[shape_function].row_index[d], q_point);
4625 
4626  return return_value;
4627  }
4628  }
4629 
4630 
4631 
4632  template <int dim, int spacedim>
4633  inline typename Vector<dim, spacedim>::gradient_type
4634  Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4635  const unsigned int q_point) const
4636  {
4637  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4638  Assert(fe_values->update_flags & update_gradients,
4640  "update_gradients")));
4641 
4642  // same as for the scalar case except that we have one more index
4643  const int snc =
4644  shape_function_data[shape_function].single_nonzero_component;
4645  if (snc == -2)
4646  return gradient_type();
4647  else if (snc != -1)
4648  {
4649  gradient_type return_value;
4650  return_value[shape_function_data[shape_function]
4651  .single_nonzero_component_index] =
4652  fe_values->finite_element_output.shape_gradients[snc][q_point];
4653  return return_value;
4654  }
4655  else
4656  {
4657  gradient_type return_value;
4658  for (unsigned int d = 0; d < dim; ++d)
4659  if (shape_function_data[shape_function]
4660  .is_nonzero_shape_function_component[d])
4661  return_value[d] =
4662  fe_values->finite_element_output.shape_gradients
4663  [shape_function_data[shape_function].row_index[d]][q_point];
4664 
4665  return return_value;
4666  }
4667  }
4668 
4669 
4670 
4671  template <int dim, int spacedim>
4672  inline typename Vector<dim, spacedim>::divergence_type
4673  Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4674  const unsigned int q_point) const
4675  {
4676  // this function works like in the case above
4677  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4678  Assert(fe_values->update_flags & update_gradients,
4680  "update_gradients")));
4681 
4682  // same as for the scalar case except that we have one more index
4683  const int snc =
4684  shape_function_data[shape_function].single_nonzero_component;
4685  if (snc == -2)
4686  return divergence_type();
4687  else if (snc != -1)
4688  return fe_values->finite_element_output
4689  .shape_gradients[snc][q_point][shape_function_data[shape_function]
4690  .single_nonzero_component_index];
4691  else
4692  {
4693  divergence_type return_value = 0;
4694  for (unsigned int d = 0; d < dim; ++d)
4695  if (shape_function_data[shape_function]
4696  .is_nonzero_shape_function_component[d])
4697  return_value +=
4698  fe_values->finite_element_output.shape_gradients
4699  [shape_function_data[shape_function].row_index[d]][q_point][d];
4700 
4701  return return_value;
4702  }
4703  }
4704 
4705 
4706 
4707  template <int dim, int spacedim>
4708  inline typename Vector<dim, spacedim>::curl_type
4709  Vector<dim, spacedim>::curl(const unsigned int shape_function,
4710  const unsigned int q_point) const
4711  {
4712  // this function works like in the case above
4713 
4714  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4715  Assert(fe_values->update_flags & update_gradients,
4717  "update_gradients")));
4718  // same as for the scalar case except that we have one more index
4719  const int snc =
4720  shape_function_data[shape_function].single_nonzero_component;
4721 
4722  if (snc == -2)
4723  return curl_type();
4724 
4725  else
4726  switch (dim)
4727  {
4728  case 1:
4729  {
4730  Assert(false,
4731  ExcMessage(
4732  "Computing the curl in 1d is not a useful operation"));
4733  return curl_type();
4734  }
4735 
4736  case 2:
4737  {
4738  if (snc != -1)
4739  {
4740  curl_type return_value;
4741 
4742  // the single nonzero component can only be zero or one in 2d
4743  if (shape_function_data[shape_function]
4744  .single_nonzero_component_index == 0)
4745  return_value[0] =
4746  -1.0 * fe_values->finite_element_output
4747  .shape_gradients[snc][q_point][1];
4748  else
4749  return_value[0] = fe_values->finite_element_output
4750  .shape_gradients[snc][q_point][0];
4751 
4752  return return_value;
4753  }
4754 
4755  else
4756  {
4757  curl_type return_value;
4758 
4759  return_value[0] = 0.0;
4760 
4761  if (shape_function_data[shape_function]
4762  .is_nonzero_shape_function_component[0])
4763  return_value[0] -=
4764  fe_values->finite_element_output
4765  .shape_gradients[shape_function_data[shape_function]
4766  .row_index[0]][q_point][1];
4767 
4768  if (shape_function_data[shape_function]
4769  .is_nonzero_shape_function_component[1])
4770  return_value[0] +=
4771  fe_values->finite_element_output
4772  .shape_gradients[shape_function_data[shape_function]
4773  .row_index[1]][q_point][0];
4774 
4775  return return_value;
4776  }
4777  }
4778 
4779  case 3:
4780  {
4781  if (snc != -1)
4782  {
4783  curl_type return_value;
4784 
4785  switch (shape_function_data[shape_function]
4786  .single_nonzero_component_index)
4787  {
4788  case 0:
4789  {
4790  return_value[0] = 0;
4791  return_value[1] = fe_values->finite_element_output
4792  .shape_gradients[snc][q_point][2];
4793  return_value[2] =
4794  -1.0 * fe_values->finite_element_output
4795  .shape_gradients[snc][q_point][1];
4796  return return_value;
4797  }
4798 
4799  case 1:
4800  {
4801  return_value[0] =
4802  -1.0 * fe_values->finite_element_output
4803  .shape_gradients[snc][q_point][2];
4804  return_value[1] = 0;
4805  return_value[2] = fe_values->finite_element_output
4806  .shape_gradients[snc][q_point][0];
4807  return return_value;
4808  }
4809 
4810  default:
4811  {
4812  return_value[0] = fe_values->finite_element_output
4813  .shape_gradients[snc][q_point][1];
4814  return_value[1] =
4815  -1.0 * fe_values->finite_element_output
4816  .shape_gradients[snc][q_point][0];
4817  return_value[2] = 0;
4818  return return_value;
4819  }
4820  }
4821  }
4822 
4823  else
4824  {
4825  curl_type return_value;
4826 
4827  for (unsigned int i = 0; i < dim; ++i)
4828  return_value[i] = 0.0;
4829 
4830  if (shape_function_data[shape_function]
4831  .is_nonzero_shape_function_component[0])
4832  {
4833  return_value[1] +=
4834  fe_values->finite_element_output
4835  .shape_gradients[shape_function_data[shape_function]
4836  .row_index[0]][q_point][2];
4837  return_value[2] -=
4838  fe_values->finite_element_output
4839  .shape_gradients[shape_function_data[shape_function]
4840  .row_index[0]][q_point][1];
4841  }
4842 
4843  if (shape_function_data[shape_function]
4844  .is_nonzero_shape_function_component[1])
4845  {
4846  return_value[0] -=
4847  fe_values->finite_element_output
4848  .shape_gradients[shape_function_data[shape_function]
4849  .row_index[1]][q_point][2];
4850  return_value[2] +=
4851  fe_values->finite_element_output
4852  .shape_gradients[shape_function_data[shape_function]
4853  .row_index[1]][q_point][0];
4854  }
4855 
4856  if (shape_function_data[shape_function]
4857  .is_nonzero_shape_function_component[2])
4858  {
4859  return_value[0] +=
4860  fe_values->finite_element_output
4861  .shape_gradients[shape_function_data[shape_function]
4862  .row_index[2]][q_point][1];
4863  return_value[1] -=
4864  fe_values->finite_element_output
4865  .shape_gradients[shape_function_data[shape_function]
4866  .row_index[2]][q_point][0];
4867  }
4868 
4869  return return_value;
4870  }
4871  }
4872  }
4873  // should not end up here
4874  Assert(false, ExcInternalError());
4875  return curl_type();
4876  }
4877 
4878 
4879 
4880  template <int dim, int spacedim>
4881  inline typename Vector<dim, spacedim>::hessian_type
4882  Vector<dim, spacedim>::hessian(const unsigned int shape_function,
4883  const unsigned int q_point) const
4884  {
4885  // this function works like in the case above
4886  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4887  Assert(fe_values->update_flags & update_hessians,
4889  "update_hessians")));
4890 
4891  // same as for the scalar case except that we have one more index
4892  const int snc =
4893  shape_function_data[shape_function].single_nonzero_component;
4894  if (snc == -2)
4895  return hessian_type();
4896  else if (snc != -1)
4897  {
4898  hessian_type return_value;
4899  return_value[shape_function_data[shape_function]
4900  .single_nonzero_component_index] =
4901  fe_values->finite_element_output.shape_hessians[snc][q_point];
4902  return return_value;
4903  }
4904  else
4905  {
4906  hessian_type return_value;
4907  for (unsigned int d = 0; d < dim; ++d)
4908  if (shape_function_data[shape_function]
4909  .is_nonzero_shape_function_component[d])
4910  return_value[d] =
4911  fe_values->finite_element_output.shape_hessians
4912  [shape_function_data[shape_function].row_index[d]][q_point];
4913 
4914  return return_value;
4915  }
4916  }
4917 
4918 
4919 
4920  template <int dim, int spacedim>
4921  inline typename Vector<dim, spacedim>::third_derivative_type
4922  Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
4923  const unsigned int q_point) const
4924  {
4925  // this function works like in the case above
4926  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4927  Assert(fe_values->update_flags & update_3rd_derivatives,
4929  "update_3rd_derivatives")));
4930 
4931  // same as for the scalar case except that we have one more index
4932  const int snc =
4933  shape_function_data[shape_function].single_nonzero_component;
4934  if (snc == -2)
4935  return third_derivative_type();
4936  else if (snc != -1)
4937  {
4938  third_derivative_type return_value;
4939  return_value[shape_function_data[shape_function]
4940  .single_nonzero_component_index] =
4941  fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
4942  return return_value;
4943  }
4944  else
4945  {
4946  third_derivative_type return_value;
4947  for (unsigned int d = 0; d < dim; ++d)
4948  if (shape_function_data[shape_function]
4949  .is_nonzero_shape_function_component[d])
4950  return_value[d] =
4951  fe_values->finite_element_output.shape_3rd_derivatives
4952  [shape_function_data[shape_function].row_index[d]][q_point];
4953 
4954  return return_value;
4955  }
4956  }
4957 
4958 
4959 
4960  namespace internal
4961  {
4966  inline ::SymmetricTensor<2, 1>
4967  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
4968  {
4969  AssertIndexRange(n, 1);
4970  (void)n;
4971 
4972  return {{t[0]}};
4973  }
4974 
4975 
4976 
4977  inline ::SymmetricTensor<2, 2>
4978  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
4979  {
4980  switch (n)
4981  {
4982  case 0:
4983  {
4984  return {{t[0], 0, t[1] / 2}};
4985  }
4986  case 1:
4987  {
4988  return {{0, t[1], t[0] / 2}};
4989  }
4990  default:
4991  {
4992  AssertIndexRange(n, 2);
4993  return {};
4994  }
4995  }
4996  }
4997 
4998 
4999 
5000  inline ::SymmetricTensor<2, 3>
5001  symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
5002  {
5003  switch (n)
5004  {
5005  case 0:
5006  {
5007  return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
5008  }
5009  case 1:
5010  {
5011  return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
5012  }
5013  case 2:
5014  {
5015  return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
5016  }
5017  default:
5018  {
5019  AssertIndexRange(n, 3);
5020  return {};
5021  }
5022  }
5023  }
5024  } // namespace internal
5025 
5026 
5027 
5028  template <int dim, int spacedim>
5029  inline typename Vector<dim, spacedim>::symmetric_gradient_type
5030  Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
5031  const unsigned int q_point) const
5032  {
5033  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5034  Assert(fe_values->update_flags & update_gradients,
5036  "update_gradients")));
5037 
5038  // same as for the scalar case except that we have one more index
5039  const int snc =
5040  shape_function_data[shape_function].single_nonzero_component;
5041  if (snc == -2)
5042  return symmetric_gradient_type();
5043  else if (snc != -1)
5044  return internal::symmetrize_single_row(
5045  shape_function_data[shape_function].single_nonzero_component_index,
5046  fe_values->finite_element_output.shape_gradients[snc][q_point]);
5047  else
5048  {
5049  gradient_type return_value;
5050  for (unsigned int d = 0; d < dim; ++d)
5051  if (shape_function_data[shape_function]
5052  .is_nonzero_shape_function_component[d])
5053  return_value[d] =
5054  fe_values->finite_element_output.shape_gradients
5055  [shape_function_data[shape_function].row_index[d]][q_point];
5056 
5057  return symmetrize(return_value);
5058  }
5059  }
5060 
5061 
5062 
5063  template <int dim, int spacedim>
5065  SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
5066  const unsigned int q_point) const
5067  {
5068  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5069  Assert(fe_values->update_flags & update_values,
5071  "update_values")));
5072 
5073  // similar to the vector case where we have more then one index and we need
5074  // to convert between unrolled and component indexing for tensors
5075  const int snc =
5076  shape_function_data[shape_function].single_nonzero_component;
5077 
5078  if (snc == -2)
5079  {
5080  // shape function is zero for the selected components
5081  return value_type();
5082  }
5083  else if (snc != -1)
5084  {
5085  value_type return_value;
5086  const unsigned int comp =
5087  shape_function_data[shape_function].single_nonzero_component_index;
5088  return_value[value_type::unrolled_to_component_indices(comp)] =
5089  fe_values->finite_element_output.shape_values(snc, q_point);
5090  return return_value;
5091  }
5092  else
5093  {
5094  value_type return_value;
5095  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
5096  if (shape_function_data[shape_function]
5097  .is_nonzero_shape_function_component[d])
5098  return_value[value_type::unrolled_to_component_indices(d)] =
5099  fe_values->finite_element_output.shape_values(
5100  shape_function_data[shape_function].row_index[d], q_point);
5101  return return_value;
5102  }
5103  }
5104 
5105 
5106 
5107  template <int dim, int spacedim>
5110  const unsigned int shape_function,
5111  const unsigned int q_point) const
5112  {
5113  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5114  Assert(fe_values->update_flags & update_gradients,
5116  "update_gradients")));
5117 
5118  const int snc =
5119  shape_function_data[shape_function].single_nonzero_component;
5120 
5121  if (snc == -2)
5122  {
5123  // shape function is zero for the selected components
5124  return divergence_type();
5125  }
5126  else if (snc != -1)
5127  {
5128  // we have a single non-zero component when the symmetric tensor is
5129  // represented in unrolled form. this implies we potentially have
5130  // two non-zero components when represented in component form! we
5131  // will only have one non-zero entry if the non-zero component lies on
5132  // the diagonal of the tensor.
5133  //
5134  // the divergence of a second-order tensor is a first order tensor.
5135  //
5136  // assume the second-order tensor is A with components A_{ij}. then
5137  // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
5138  // entries in the tensorial representation. define the
5139  // divergence as:
5140  // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
5141  // (which is incidentally also
5142  // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
5143  // In both cases, a sum is implied.
5144  //
5145  // Now, we know the nonzero component in unrolled form: it is indicated
5146  // by 'snc'. we can figure out which tensor components belong to this:
5147  const unsigned int comp =
5148  shape_function_data[shape_function].single_nonzero_component_index;
5149  const unsigned int ii =
5150  value_type::unrolled_to_component_indices(comp)[0];
5151  const unsigned int jj =
5152  value_type::unrolled_to_component_indices(comp)[1];
5153 
5154  // given the form of the divergence above, if ii=jj there is only a
5155  // single nonzero component of the full tensor and the gradient
5156  // equals
5157  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
5158  // all other entries of 'b' are zero
5159  //
5160  // on the other hand, if ii!=jj, then there are two nonzero entries in
5161  // the full tensor and
5162  // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
5163  // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
5164  // again, all other entries of 'b' are zero
5165  const ::Tensor<1, spacedim> &phi_grad =
5166  fe_values->finite_element_output.shape_gradients[snc][q_point];
5167 
5168  divergence_type return_value;
5169  return_value[ii] = phi_grad[jj];
5170 
5171  if (ii != jj)
5172  return_value[jj] = phi_grad[ii];
5173 
5174  return return_value;
5175  }
5176  else
5177  {
5178  Assert(false, ExcNotImplemented());
5179  divergence_type return_value;
5180  return return_value;
5181  }
5182  }
5183 
5184 
5185 
5186  template <int dim, int spacedim>
5187  inline typename Tensor<2, dim, spacedim>::value_type
5188  Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
5189  const unsigned int q_point) const
5190  {
5191  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5192  Assert(fe_values->update_flags & update_values,
5194  "update_values")));
5195 
5196  // similar to the vector case where we have more then one index and we need
5197  // to convert between unrolled and component indexing for tensors
5198  const int snc =
5199  shape_function_data[shape_function].single_nonzero_component;
5200 
5201  if (snc == -2)
5202  {
5203  // shape function is zero for the selected components
5204  return value_type();
5205  }
5206  else if (snc != -1)
5207  {
5208  value_type return_value;
5209  const unsigned int comp =
5210  shape_function_data[shape_function].single_nonzero_component_index;
5211  const TableIndices<2> indices =
5213  return_value[indices] =
5214  fe_values->finite_element_output.shape_values(snc, q_point);
5215  return return_value;
5216  }
5217  else
5218  {
5219  value_type return_value;
5220  for (unsigned int d = 0; d < dim * dim; ++d)
5221  if (shape_function_data[shape_function]
5222  .is_nonzero_shape_function_component[d])
5223  {
5224  const TableIndices<2> indices =
5226  return_value[indices] =
5227  fe_values->finite_element_output.shape_values(
5228  shape_function_data[shape_function].row_index[d], q_point);
5229  }
5230  return return_value;
5231  }
5232  }
5233 
5234 
5235 
5236  template <int dim, int spacedim>
5238  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
5239  const unsigned int q_point) const
5240  {
5241  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5242  Assert(fe_values->update_flags & update_gradients,
5244  "update_gradients")));
5245 
5246  const int snc =
5247  shape_function_data[shape_function].single_nonzero_component;
5248 
5249  if (snc == -2)
5250  {
5251  // shape function is zero for the selected components
5252  return divergence_type();
5253  }
5254  else if (snc != -1)
5255  {
5256  // we have a single non-zero component when the tensor is
5257  // represented in unrolled form.
5258  //
5259  // the divergence of a second-order tensor is a first order tensor.
5260  //
5261  // assume the second-order tensor is A with components A_{ij},
5262  // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
5263  //
5264  // Now, we know the nonzero component in unrolled form: it is indicated
5265  // by 'snc'. we can figure out which tensor components belong to this:
5266  const unsigned int comp =
5267  shape_function_data[shape_function].single_nonzero_component_index;
5268  const TableIndices<2> indices =
5270  const unsigned int ii = indices[0];
5271  const unsigned int jj = indices[1];
5272 
5273  const ::Tensor<1, spacedim> &phi_grad =
5274  fe_values->finite_element_output.shape_gradients[snc][q_point];
5275 
5276  divergence_type return_value;
5277  // note that we contract \nabla from the right
5278  return_value[ii] = phi_grad[jj];
5279 
5280  return return_value;
5281  }
5282  else
5283  {
5284  Assert(false, ExcNotImplemented());
5285  divergence_type return_value;
5286  return return_value;
5287  }
5288  }
5289 
5290 
5291 
5292  template <int dim, int spacedim>
5294  Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
5295  const unsigned int q_point) const
5296  {
5297  AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5298  Assert(fe_values->update_flags & update_gradients,
5300  "update_gradients")));
5301 
5302  const int snc =
5303  shape_function_data[shape_function].single_nonzero_component;
5304 
5305  if (snc == -2)
5306  {
5307  // shape function is zero for the selected components
5308  return gradient_type();
5309  }
5310  else if (snc != -1)
5311  {
5312  // we have a single non-zero component when the tensor is
5313  // represented in unrolled form.
5314  //
5315  // the gradient of a second-order tensor is a third order tensor.
5316  //
5317  // assume the second-order tensor is A with components A_{ij},
5318  // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
5319  //
5320  // Now, we know the nonzero component in unrolled form: it is indicated
5321  // by 'snc'. we can figure out which tensor components belong to this:
5322  const unsigned int comp =
5323  shape_function_data[shape_function].single_nonzero_component_index;
5324  const TableIndices<2> indices =
5326  const unsigned int ii = indices[0];
5327  const unsigned int jj = indices[1];
5328 
5329  const ::Tensor<1, spacedim> &phi_grad =
5330  fe_values->finite_element_output.shape_gradients[snc][q_point];
5331 
5332  gradient_type return_value;
5333  return_value[ii][jj] = phi_grad;
5334 
5335  return return_value;
5336  }
5337  else
5338  {
5339  Assert(false, ExcNotImplemented());
5340  gradient_type return_value;
5341  return return_value;
5342  }
5343  }
5344 
5345 } // namespace FEValuesViews
5346 
5347 
5348 
5349 /*---------------------- Inline functions: FEValuesBase ---------------------*/
5350 
5351 
5352 
5353 template <int dim, int spacedim>
5355  operator[](const FEValuesExtractors::Scalar &scalar) const
5356 {
5357  AssertIndexRange(scalar.component, fe_values_views_cache.scalars.size());
5358 
5359  return fe_values_views_cache.scalars[scalar.component];
5360 }
5361 
5362 
5363 
5364 template <int dim, int spacedim>
5366  operator[](const FEValuesExtractors::Vector &vector) const
5367 {
5369  fe_values_views_cache.vectors.size());
5370 
5371  return fe_values_views_cache.vectors[vector.first_vector_component];
5372 }
5373 
5374 
5375 
5376 template <int dim, int spacedim>
5380 {
5381  Assert(
5382  tensor.first_tensor_component <
5383  fe_values_views_cache.symmetric_second_order_tensors.size(),
5385  0,
5386  fe_values_views_cache.symmetric_second_order_tensors.size()));
5387 
5388  return fe_values_views_cache
5389  .symmetric_second_order_tensors[tensor.first_tensor_component];
5390 }
5391 
5392 
5393 
5394 template <int dim, int spacedim>
5397  operator[](const FEValuesExtractors::Tensor<2> &tensor) const
5398 {
5400  fe_values_views_cache.second_order_tensors.size());
5401 
5402  return fe_values_views_cache
5403  .second_order_tensors[tensor.first_tensor_component];
5404 }
5405 
5406 
5407 
5408 template <int dim, int spacedim>
5409 inline const double &
5410 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
5411  const unsigned int j) const
5412 {
5413  AssertIndexRange(i, fe->n_dofs_per_cell());
5414  Assert(this->update_flags & update_values,
5415  ExcAccessToUninitializedField("update_values"));
5416  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5417  Assert(present_cell.get() != nullptr,
5418  ExcMessage("FEValues object is not reinit'ed to any cell"));
5419  // if the entire FE is primitive,
5420  // then we can take a short-cut:
5421  if (fe->is_primitive())
5422  return this->finite_element_output.shape_values(i, j);
5423  else
5424  {
5425  // otherwise, use the mapping
5426  // between shape function
5427  // numbers and rows. note that
5428  // by the assertions above, we
5429  // know that this particular
5430  // shape function is primitive,
5431  // so we can call
5432  // system_to_component_index
5433  const unsigned int row =
5434  this->finite_element_output
5435  .shape_function_to_row_table[i * fe->n_components() +
5436  fe->system_to_component_index(i).first];
5437  return this->finite_element_output.shape_values(row, j);
5438  }
5439 }
5440 
5441 
5442 
5443 template <int dim, int spacedim>
5444 inline double
5446  const unsigned int i,
5447  const unsigned int j,
5448  const unsigned int component) const
5449 {
5450  AssertIndexRange(i, fe->n_dofs_per_cell());
5451  Assert(this->update_flags & update_values,
5452  ExcAccessToUninitializedField("update_values"));
5453  AssertIndexRange(component, fe->n_components());
5454  Assert(present_cell.get() != nullptr,
5455  ExcMessage("FEValues object is not reinit'ed to any cell"));
5456 
5457  // check whether the shape function
5458  // is non-zero at all within
5459  // this component:
5460  if (fe->get_nonzero_components(i)[component] == false)
5461  return 0;
5462 
5463  // look up the right row in the
5464  // table and take the data from
5465  // there
5466  const unsigned int row =
5467  this->finite_element_output
5468  .shape_function_to_row_table[i * fe->n_components() + component];
5469  return this->finite_element_output.shape_values(row, j);
5470 }
5471 
5472 
5473 
5474 template <int dim, int spacedim>
5475 inline const Tensor<1, spacedim> &
5476 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
5477  const unsigned int j) const
5478 {
5479  AssertIndexRange(i, fe->n_dofs_per_cell());
5480  Assert(this->update_flags & update_gradients,
5481  ExcAccessToUninitializedField("update_gradients"));
5482  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5483  Assert(present_cell.get() != nullptr,
5484  ExcMessage("FEValues object is not reinit'ed to any cell"));
5485  // if the entire FE is primitive,
5486  // then we can take a short-cut:
5487  if (fe->is_primitive())
5488  return this->finite_element_output.shape_gradients[i][j];
5489  else
5490  {
5491  // otherwise, use the mapping
5492  // between shape function
5493  // numbers and rows. note that
5494  // by the assertions above, we
5495  // know that this particular
5496  // shape function is primitive,
5497  // so we can call
5498  // system_to_component_index
5499  const unsigned int row =
5500  this->finite_element_output
5501  .shape_function_to_row_table[i * fe->n_components() +
5502  fe->system_to_component_index(i).first];
5503  return this->finite_element_output.shape_gradients[row][j];
5504  }
5505 }
5506 
5507 
5508 
5509 template <int dim, int spacedim>
5510 inline Tensor<1, spacedim>
5512  const unsigned int i,
5513  const unsigned int j,
5514  const unsigned int component) const
5515 {
5516  AssertIndexRange(i, fe->n_dofs_per_cell());
5517  Assert(this->update_flags & update_gradients,
5518  ExcAccessToUninitializedField("update_gradients"));
5519  AssertIndexRange(component, fe->n_components());
5520  Assert(present_cell.get() != nullptr,
5521  ExcMessage("FEValues object is not reinit'ed to any cell"));
5522  // check whether the shape function
5523  // is non-zero at all within
5524  // this component:
5525  if (fe->get_nonzero_components(i)[component] == false)
5526  return Tensor<1, spacedim>();
5527 
5528  // look up the right row in the
5529  // table and take the data from
5530  // there
5531  const unsigned int row =
5532  this->finite_element_output
5533  .shape_function_to_row_table[i * fe->n_components() + component];
5534  return this->finite_element_output.shape_gradients[row][j];
5535 }
5536 
5537 
5538 
5539 template <int dim, int spacedim>
5540 inline const Tensor<2, spacedim> &
5541 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
5542  const unsigned int j) const
5543 {
5544  AssertIndexRange(i, fe->n_dofs_per_cell());
5545  Assert(this->update_flags & update_hessians,
5546  ExcAccessToUninitializedField("update_hessians"));
5547  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5548  Assert(present_cell.get() != nullptr,
5549  ExcMessage("FEValues object is not reinit'ed to any cell"));
5550  // if the entire FE is primitive,
5551  // then we can take a short-cut:
5552  if (fe->is_primitive())
5553  return this->finite_element_output.shape_hessians[i][j];
5554  else
5555  {
5556  // otherwise, use the mapping
5557  // between shape function
5558  // numbers and rows. note that
5559  // by the assertions above, we
5560  // know that this particular
5561  // shape function is primitive,
5562  // so we can call
5563  // system_to_component_index
5564  const unsigned int row =
5565  this->finite_element_output
5566  .shape_function_to_row_table[i * fe->n_components() +
5567  fe->system_to_component_index(i).first];
5568  return this->finite_element_output.shape_hessians[row][j];
5569  }
5570 }
5571 
5572 
5573 
5574 template <int dim, int spacedim>
5575 inline Tensor<2, spacedim>
5577  const unsigned int i,
5578  const unsigned int j,
5579  const unsigned int component) const
5580 {
5581  AssertIndexRange(i, fe->n_dofs_per_cell());
5582  Assert(this->update_flags & update_hessians,
5583  ExcAccessToUninitializedField("update_hessians"));
5584  AssertIndexRange(component, fe->n_components());
5585  Assert(present_cell.get() != nullptr,
5586  ExcMessage("FEValues object is not reinit'ed to any cell"));
5587  // check whether the shape function
5588  // is non-zero at all within
5589  // this component:
5590  if (fe->get_nonzero_components(i)[component] == false)
5591  return Tensor<2, spacedim>();
5592 
5593  // look up the right row in the
5594  // table and take the data from
5595  // there
5596  const unsigned int row =
5597  this->finite_element_output
5598  .shape_function_to_row_table[i * fe->n_components() + component];
5599  return this->finite_element_output.shape_hessians[row][j];
5600 }
5601 
5602 
5603 
5604 template <int dim, int spacedim>
5605 inline const Tensor<3, spacedim> &
5607  const unsigned int j) const
5608 {
5609  AssertIndexRange(i, fe->n_dofs_per_cell());
5610  Assert(this->update_flags & update_3rd_derivatives,
5611  ExcAccessToUninitializedField("update_3rd_derivatives"));
5612  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5613  Assert(present_cell.get() != nullptr,
5614  ExcMessage("FEValues object is not reinit'ed to any cell"));
5615  // if the entire FE is primitive,
5616  // then we can take a short-cut:
5617  if (fe->is_primitive())
5618  return this->finite_element_output.shape_3rd_derivatives[i][j];
5619  else
5620  {
5621  // otherwise, use the mapping
5622  // between shape function
5623  // numbers and rows. note that
5624  // by the assertions above, we
5625  // know that this particular
5626  // shape function is primitive,
5627  // so we can call
5628  // system_to_component_index
5629  const unsigned int row =
5630  this->finite_element_output
5631  .shape_function_to_row_table[i * fe->n_components() +
5632  fe->system_to_component_index(i).first];
5633  return this->finite_element_output.shape_3rd_derivatives[row][j];
5634  }
5635 }
5636 
5637 
5638 
5639 template <int dim, int spacedim>
5640 inline Tensor<3, spacedim>
5642  const unsigned int i,
5643  const unsigned int j,
5644  const unsigned int component) const
5645 {
5646  AssertIndexRange(i, fe->n_dofs_per_cell());
5647  Assert(this->update_flags & update_3rd_derivatives,
5648  ExcAccessToUninitializedField("update_3rd_derivatives"));
5649  AssertIndexRange(component, fe->n_components());
5650  Assert(present_cell.get() != nullptr,
5651  ExcMessage("FEValues object is not reinit'ed to any cell"));
5652  // check whether the shape function
5653  // is non-zero at all within
5654  // this component:
5655  if (fe->get_nonzero_components(i)[component] == false)
5656  return Tensor<3, spacedim>();
5657 
5658  // look up the right row in the
5659  // table and take the data from
5660  // there
5661  const unsigned int row =
5662  this->finite_element_output
5663  .shape_function_to_row_table[i * fe->n_components() + component];
5664  return this->finite_element_output.shape_3rd_derivatives[row][j];
5665 }
5666 
5667 
5668 
5669 template <int dim, int spacedim>
5670 inline const FiniteElement<dim, spacedim> &
5672 {
5673  return *fe;
5674 }
5675 
5676 
5677 
5678 template <int dim, int spacedim>
5679 inline const Mapping<dim, spacedim> &
5681 {
5682  return *mapping;
5683 }
5684 
5685 
5686 
5687 template <int dim, int spacedim>
5688 inline UpdateFlags
5690 {
5691  return this->update_flags;
5692 }
5693 
5694 
5695 
5696 template <int dim, int spacedim>
5697 inline const std::vector<Point<spacedim>> &
5699 {
5700  Assert(this->update_flags & update_quadrature_points,
5701  ExcAccessToUninitializedField("update_quadrature_points"));
5702  Assert(present_cell.get() != nullptr,
5703  ExcMessage("FEValues object is not reinit'ed to any cell"));
5704  return this->mapping_output.quadrature_points;
5705 }
5706 
5707 
5708 
5709 template <int dim, int spacedim>
5710 inline const std::vector<double> &
5712 {
5713  Assert(this->update_flags & update_JxW_values,
5714  ExcAccessToUninitializedField("update_JxW_values"));
5715  Assert(present_cell.get() != nullptr,
5716  ExcMessage("FEValues object is not reinit'ed to any cell"));
5717  return this->mapping_output.JxW_values;
5718 }
5719 
5720 
5721 
5722 template <int dim, int spacedim>
5723 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5725 {
5726  Assert(this->update_flags & update_jacobians,
5727  ExcAccessToUninitializedField("update_jacobians"));
5728  Assert(present_cell.get() != nullptr,
5729  ExcMessage("FEValues object is not reinit'ed to any cell"));
5730  return this->mapping_output.jacobians;
5731 }
5732 
5733 
5734 
5735 template <int dim, int spacedim>
5736 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5738 {
5739  Assert(this->update_flags & update_jacobian_grads,
5740  ExcAccessToUninitializedField("update_jacobians_grads"));
5741  Assert(present_cell.get() != nullptr,
5742  ExcMessage("FEValues object is not reinit'ed to any cell"));
5743  return this->mapping_output.jacobian_grads;
5744 }
5745 
5746 
5747 
5748 template <int dim, int spacedim>
5749 inline const Tensor<3, spacedim> &
5751  const unsigned int i) const
5752 {
5753  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5754  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5755  Assert(present_cell.get() != nullptr,
5756  ExcMessage("FEValues object is not reinit'ed to any cell"));
5757  return this->mapping_output.jacobian_pushed_forward_grads[i];
5758 }
5759 
5760 
5761 
5762 template <int dim, int spacedim>
5763 inline const std::vector<Tensor<3, spacedim>> &
5765 {
5766  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5767  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5768  Assert(present_cell.get() != nullptr,
5769  ExcMessage("FEValues object is not reinit'ed to any cell"));
5770  return this->mapping_output.jacobian_pushed_forward_grads;
5771 }
5772 
5773 
5774 
5775 template <int dim, int spacedim>
5776 inline const DerivativeForm<3, dim, spacedim> &
5777 FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
5778 {
5779  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5780  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5781  Assert(present_cell.get() != nullptr,
5782  ExcMessage("FEValues object is not reinit'ed to any cell"));
5783  return this->mapping_output.jacobian_2nd_derivatives[i];
5784 }
5785 
5786 
5787 
5788 template <int dim, int spacedim>
5789 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5791 {
5792  Assert(this->update_flags & update_jacobian_2nd_derivatives,
5793  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5794  Assert(present_cell.get() != nullptr,
5795  ExcMessage("FEValues object is not reinit'ed to any cell"));
5796  return this->mapping_output.jacobian_2nd_derivatives;
5797 }
5798 
5799 
5800 
5801 template <int dim, int spacedim>
5802 inline const Tensor<4, spacedim> &
5804  const unsigned int i) const
5805 {
5808  "update_jacobian_pushed_forward_2nd_derivatives"));
5809  Assert(present_cell.get() != nullptr,
5810  ExcMessage("FEValues object is not reinit'ed to any cell"));
5811  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5812 }
5813 
5814 
5815 
5816 template <int dim, int spacedim>
5817 inline const std::vector<Tensor<4, spacedim>> &
5819 {
5820  Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
5822  "update_jacobian_pushed_forward_2nd_derivatives"));
5823  Assert(present_cell.get() != nullptr,
5824  ExcMessage("FEValues object is not reinit'ed to any cell"));
5825  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5826 }
5827 
5828 
5829 
5830 template <int dim, int spacedim>
5831 inline const DerivativeForm<4, dim, spacedim> &
5832 FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
5833 {
5834  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5835  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5836  Assert(present_cell.get() != nullptr,
5837  ExcMessage("FEValues object is not reinit'ed to any cell"));
5838  return this->mapping_output.jacobian_3rd_derivatives[i];
5839 }
5840 
5841 
5842 
5843 template <int dim, int spacedim>
5844 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5846 {
5847  Assert(this->update_flags & update_jacobian_3rd_derivatives,
5848  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5849  Assert(present_cell.get() != nullptr,
5850  ExcMessage("FEValues object is not reinit'ed to any cell"));
5851  return this->mapping_output.jacobian_3rd_derivatives;
5852 }
5853 
5854 
5855 
5856 template <int dim, int spacedim>
5857 inline const Tensor<5, spacedim> &
5859  const unsigned int i) const
5860 {
5863  "update_jacobian_pushed_forward_3rd_derivatives"));
5864  Assert(present_cell.get() != nullptr,
5865  ExcMessage("FEValues object is not reinit'ed to any cell"));
5866  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
5867 }
5868 
5869 
5870 
5871 template <int dim, int spacedim>
5872 inline const std::vector<Tensor<5, spacedim>> &
5874 {
5875  Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
5877  "update_jacobian_pushed_forward_3rd_derivatives"));
5878  Assert(present_cell.get() != nullptr,
5879  ExcMessage("FEValues object is not reinit'ed to any cell"));
5880  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
5881 }
5882 
5883 
5884 
5885 template <int dim, int spacedim>
5886 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
5888 {
5889  Assert(this->update_flags & update_inverse_jacobians,
5890  ExcAccessToUninitializedField("update_inverse_jacobians"));
5891  Assert(present_cell.get() != nullptr,
5892  ExcMessage("FEValues object is not reinit'ed to any cell"));
5893  return this->mapping_output.inverse_jacobians;
5894 }
5895 
5896 
5897 
5898 template <int dim, int spacedim>
5901 {
5902  return {0U, dofs_per_cell};
5903 }
5904 
5905 
5906 
5907 template <int dim, int spacedim>
5910  const unsigned int start_dof_index) const
5911 {
5912  Assert(start_dof_index <= dofs_per_cell,
5913  ExcIndexRange(start_dof_index, 0, dofs_per_cell + 1));
5914  return {start_dof_index, dofs_per_cell};
5915 }
5916 
5917 
5918 
5919 template <int dim, int spacedim>
5922  const unsigned int end_dof_index) const
5923 {
5924  Assert(end_dof_index < dofs_per_cell,
5925  ExcIndexRange(end_dof_index, 0, dofs_per_cell));
5926  return {0U, end_dof_index + 1};
5927 }
5928 
5929 
5930 
5931 template <int dim, int spacedim>
5934 {
5935  return {0U, n_quadrature_points};
5936 }
5937 
5938 
5939 
5940 template <int dim, int spacedim>
5941 inline const Point<spacedim> &
5942 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
5943 {
5944  Assert(this->update_flags & update_quadrature_points,
5945  ExcAccessToUninitializedField("update_quadrature_points"));
5946  AssertIndexRange(i, this->mapping_output.quadrature_points.size());
5947  Assert(present_cell.get() != nullptr,
5948  ExcMessage("FEValues object is not reinit'ed to any cell"));
5949 
5950  return this->mapping_output.quadrature_points[i];
5951 }
5952 
5953 
5954 
5955 template <int dim, int spacedim>
5956 inline double
5957 FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
5958 {
5959  Assert(this->update_flags & update_JxW_values,
5960  ExcAccessToUninitializedField("update_JxW_values"));
5961  AssertIndexRange(i, this->mapping_output.JxW_values.size());
5962  Assert(present_cell.get() != nullptr,
5963  ExcMessage("FEValues object is not reinit'ed to any cell"));
5964 
5965  return this->mapping_output.JxW_values[i];
5966 }
5967 
5968 
5969 
5970 template <int dim, int spacedim>
5971 inline const DerivativeForm<1, dim, spacedim> &
5972 FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
5973 {
5974  Assert(this->update_flags & update_jacobians,
5975  ExcAccessToUninitializedField("update_jacobians"));
5976  AssertIndexRange(i, this->mapping_output.jacobians.size());
5977  Assert(present_cell.get() != nullptr,
5978  ExcMessage("FEValues object is not reinit'ed to any cell"));
5979 
5980  return this->mapping_output.jacobians[i];
5981 }
5982 
5983 
5984 
5985 template <int dim, int spacedim>
5986 inline const DerivativeForm<2, dim, spacedim> &
5987 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
5988 {
5989  Assert(this->update_flags & update_jacobian_grads,
5990  ExcAccessToUninitializedField("update_jacobians_grads"));
5991  AssertIndexRange(i, this->mapping_output.jacobian_grads.size());
5992  Assert(present_cell.get() != nullptr,
5993  ExcMessage("FEValues object is not reinit'ed to any cell"));
5994 
5995  return this->mapping_output.jacobian_grads[i];
5996 }
5997 
5998 
5999 
6000 template <int dim, int spacedim>
6001 inline const DerivativeForm<1, spacedim, dim> &
6002 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
6003 {
6004  Assert(this->update_flags & update_inverse_jacobians,
6005  ExcAccessToUninitializedField("update_inverse_jacobians"));
6006  AssertIndexRange(i, this->mapping_output.inverse_jacobians.size());
6007  Assert(present_cell.get() != nullptr,
6008  ExcMessage("FEValues object is not reinit'ed to any cell"));
6009 
6010  return this->mapping_output.inverse_jacobians[i];
6011 }
6012 
6013 
6014 
6015 template <int dim, int spacedim>
6016 inline const Tensor<1, spacedim> &
6017 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
6018 {
6019  Assert(this->update_flags & update_normal_vectors,
6021  "update_normal_vectors")));
6022  AssertIndexRange(i, this->mapping_output.normal_vectors.size());
6023  Assert(present_cell.get() != nullptr,
6024  ExcMessage("FEValues object is not reinit'ed to any cell"));
6025 
6026  return this->mapping_output.normal_vectors[i];
6027 }
6028 
6029 
6030 
6031 /*--------------------- Inline functions: FEValues --------------------------*/
6032 
6033 
6034 template <int dim, int spacedim>
6035 inline const Quadrature<dim> &
6037 {
6038  return quadrature;
6039 }
6040 
6041 
6042 
6043 template <int dim, int spacedim>
6044 inline const FEValues<dim, spacedim> &
6046 {
6047  return *this;
6048 }
6049 
6050 
6051 /*---------------------- Inline functions: FEFaceValuesBase -----------------*/
6052 
6053 
6054 template <int dim, int spacedim>
6055 inline unsigned int
6057 {
6058  return present_face_index;
6059 }
6060 
6061 
6062 /*----------------------- Inline functions: FE*FaceValues -------------------*/
6063 
6064 template <int dim, int spacedim>
6065 inline const Quadrature<dim - 1> &
6067 {
6068  return quadrature[quadrature.size() == 1 ? 0 : present_face_no];
6069 }
6070 
6071 
6072 
6073 template <int dim, int spacedim>
6074 inline const FEFaceValues<dim, spacedim> &
6076 {
6077  return *this;
6078 }
6079 
6080 
6081 
6082 template <int dim, int spacedim>
6083 inline const FESubfaceValues<dim, spacedim> &
6085 {
6086  return *this;
6087 }
6088 
6089 
6090 
6091 template <int dim, int spacedim>
6092 inline const Tensor<1, spacedim> &
6093 FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int i) const
6094 {
6095  AssertIndexRange(i, this->mapping_output.boundary_forms.size());
6096  Assert(this->update_flags & update_boundary_forms,
6098  "update_boundary_forms")));
6099 
6100  return this->mapping_output.boundary_forms[i];
6101 }
6102 
6103 #endif // DOXYGEN
6104 
6106 
6107 #endif
Transformed quadrature weights.
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
Definition: fe_values.h:788
Shape function values.
typename ProductType< Number, typename Vector< dim, spacedim >::curl_type >::type curl_type
Definition: fe_values.h:845
typename ProductType< Number, divergence_type >::type solution_divergence_type
Definition: fe_values.h:1468
const FEFaceValues< dim, spacedim > & get_present_fe_values() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_starting_at(const unsigned int start_dof_index) const
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
const Tensor< 3, spacedim > & jacobian_pushed_forward_grad(const unsigned int quadrature_point) const
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > fe_data
Definition: fe_values.h:3832
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:3864
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:1830
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1486
unsigned int present_face_no
Definition: fe_values.h:4127
unsigned int present_face_index
Definition: fe_values.h:4133
std::vector<::FEValuesViews::Vector< dim, spacedim > > vectors
Definition: fe_values.h:2198
typename ::internal::FEValuesViews::ViewType< dim, spacedim, Extractor >::type View
Definition: fe_values.h:2220
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1723
#define DEAL_II_DEPRECATED_EARLY
Definition: config.h:165
constexpr SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
static ::ExceptionBase & ExcAccessToUninitializedField()
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:605
typename internal::ProductTypeImpl< typename std::decay< T >::type, typename std::decay< U >::type >::type type
const unsigned int dofs_per_cell
Definition: fe_values.h:2362
typename ::internal::CurlType< spacedim >::type curl_type
Definition: fe_values.h:696
const unsigned int component
Definition: fe_values.h:611
typename ProductType< Number, value_type >::type solution_value_type
Definition: fe_values.h:1458
Volume element.
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
Definition: fe_values.h:224
typename ::FEValuesViews::SymmetricTensor< rank, dim, spacedim > type
Definition: fe_values.h:2180
#define AssertIndexRange(index, range)
Definition: exceptions.h:1690
const Mapping< dim, spacedim > & get_mapping() const
Outer normal vector, not normalized.
typename ProductType< Number, typename Scalar< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:265
typename ProductType< Number, divergence_type >::type solution_divergence_type
Definition: fe_values.h:1786
const DerivativeForm< 2, dim, spacedim > & jacobian_grad(const unsigned int quadrature_point) const
const FiniteElement< dim, spacedim > & get_fe() const
std::unique_ptr< const CellIteratorBase > present_cell
Definition: fe_values.h:3748
UpdateFlags get_update_flags() const
static const char U
Transformed quadrature points.
std::vector< double > get_quadrature_points(const unsigned int n)
typename ProductType< Number, value_type >::type solution_value_type
Definition: fe_values.h:1776
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:241
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:3800
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int quadrature_point) const
const DerivativeForm< 3, dim, spacedim > & jacobian_2nd_derivative(const unsigned int quadrature_point) const
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition: fe_values.h:3879
std_cxx20::ranges::iota_view< unsigned int, unsigned int > quadrature_point_indices() const
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
static ::ExceptionBase & ExcFENotPrimitive()
typename ProductType< Number, value_type >::type solution_value_type
Definition: fe_values.h:719
Tensor< 2, spacedim > shape_hessian_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
std::vector<::FEValuesViews::SymmetricTensor< 2, dim, spacedim > > symmetric_second_order_tensors
Definition: fe_values.h:2200
typename ProductType< Number, typename Vector< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:813
const Quadrature< dim - 1 > & get_quadrature() const
const std::vector< Point< spacedim > > & get_quadrature_points() const
const Point< spacedim > & quadrature_point(const unsigned int q) const
const std::vector< Tensor< 4, spacedim > > & get_jacobian_pushed_forward_2nd_derivatives() const
typename ProductType< Number, typename Scalar< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:249
const hp::QCollection< dim - 1 > quadrature
Definition: fe_values.h:4138
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1814
static ::ExceptionBase & ExcMessage(std::string arg1)
typename ProductType< Number, curl_type >::type solution_curl_type
Definition: fe_values.h:768
Number value_type
Definition: vector.h:122
typename ProductType< Number, typename Vector< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:861
static constexpr TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:515
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_ending_at(const unsigned int end_dof_index) const
Tensor< 3, spacedim > shape_3rd_derivative_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
typename ProductType< Number, hessian_type >::type solution_hessian_type
Definition: fe_values.h:778
Third derivatives of shape functions.
const FEValues< dim, spacedim > & get_present_fe_values() const
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:2112
unsigned int get_face_index() const
#define Assert(cond, exc)
Definition: exceptions.h:1465
UpdateFlags
const Tensor< 2, spacedim > & shape_hessian(const unsigned int function_no, const unsigned int point_no) const
Abstract base class for mapping classes.
Definition: mapping.h:303
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1399
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
const Quadrature< dim > quadrature
Definition: fe_values.h:4015
const unsigned int first_vector_component
Definition: fe_values.h:1394
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:395
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition: fe_values.h:3808
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int quadrature_point) const
typename ProductType< Number, value_type >::type solution_laplacian_type
Definition: fe_values.h:204
typename Tensor< rank_ - 1, dim, Number >::tensor_type value_type
Definition: tensor.h:483
typename ProductType< Number, gradient_type >::type solution_gradient_type
Definition: fe_values.h:194
const std::vector< DerivativeForm< 4, dim, spacedim > > & get_jacobian_3rd_derivatives() const
Second derivatives of shape functions.
Gradient of volume element.
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1494
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
std::vector<::FEValuesViews::Scalar< dim, spacedim > > scalars
Definition: fe_values.h:2197
const Tensor< 1, spacedim > & boundary_form(const unsigned int i) const
const Quadrature< dim > & get_quadrature() const
const unsigned int n_quadrature_points
Definition: fe_values.h:2344
typename ProductType< Number, hessian_type >::type solution_hessian_type
Definition: fe_values.h:214
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices() const
Tensor< 1, spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
typename ProductType< Number, value_type >::type solution_value_type
Definition: fe_values.h:184
typename ProductType< Number, symmetric_gradient_type >::type solution_symmetric_gradient_type
Definition: fe_values.h:739
boost::signals2::connection tria_listener_mesh_transform
Definition: fe_values.h:3773
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:805
typename ::FEValuesViews::Tensor< rank, dim, spacedim > type
Definition: fe_values.h:2166
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
Definition: tensor.h:449
const std::vector< double > & get_JxW_values() const
double JxW(const unsigned int quadrature_point) const
typename ProductType< Number, typename Vector< dim, spacedim >::symmetric_gradient_type >::type symmetric_gradient_type
Definition: fe_values.h:821
typename ProductType< Number, typename Vector< dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:829
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:394
typename ProductType< Number, typename Vector< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:853
Shape function gradients.
Normal vectors.
const std::vector< DerivativeForm< 2, dim, spacedim > > & get_jacobian_grads() const
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1822
const Tensor< 4, spacedim > & jacobian_pushed_forward_2nd_derivative(const unsigned int quadrature_point) const
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:2101
Definition: fe.h:38
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1712
typename ProductType< Number, value_type >::type solution_laplacian_type
Definition: fe_values.h:759
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int quadrature_point) const
static ::ExceptionBase & ExcNotImplemented()
const std::vector< DerivativeForm< 3, dim, spacedim > > & get_jacobian_2nd_derivatives() const
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1388
const Tensor< 1, spacedim > & shape_grad(const unsigned int function_no, const unsigned int quadrature_point) const
boost::signals2::connection tria_listener_refinement
Definition: fe_values.h:3764
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3840
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:837
typename ProductType< Number, gradient_type >::type solution_gradient_type
Definition: fe_values.h:729
::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:3815
const std::vector< Tensor< 5, spacedim > > & get_jacobian_pushed_forward_3rd_derivatives() const
const FESubfaceValues< dim, spacedim > & get_present_fe_values() const
const unsigned int max_n_quadrature_points
Definition: fe_values.h:2355
std::vector<::FEValuesViews::Tensor< 2, dim, spacedim > > second_order_tensors
Definition: fe_values.h:2202
const std::vector< DerivativeForm< 1, spacedim, dim > > & get_inverse_jacobians() const
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:257
typename ProductType< Number, typename Scalar< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:273
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:616
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
UpdateFlags update_flags
Definition: fe_values.h:3846
const Tensor< 3, spacedim > & shape_3rd_derivative(const unsigned int function_no, const unsigned int point_no) const
static ::ExceptionBase & ExcInternalError()
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
Definition: fe_values.h:3824
const Tensor< 1, spacedim > & normal_vector(const unsigned int i) const
const Tensor< 5, spacedim > & jacobian_pushed_forward_3rd_derivative(const unsigned int quadrature_point) const
typename ProductType< Number, divergence_type >::type solution_divergence_type
Definition: fe_values.h:749
typename ProductType< Number, gradient_type >::type solution_gradient_type
Definition: fe_values.h:1796