Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
fe_values.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 1998 - 2023 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>
38
39#include <deal.II/grid/tria.h>
41
43
44#include <algorithm>
45#include <memory>
46#include <type_traits>
47
48
49// dummy include in order to have the
50// definition of PetscScalar available
51// without including other PETSc stuff
52#ifdef DEAL_II_WITH_PETSC
53# include <petsc.h>
54#endif
55
57
58// Forward declaration
59#ifndef DOXYGEN
60template <int dim, int spacedim = dim>
61class FEValuesBase;
62#endif
63
64namespace internal
65{
70 template <int dim, class NumberType = double>
71 struct CurlType;
72
79 template <class NumberType>
81 {
83 };
84
91 template <class NumberType>
93 {
95 };
96
103 template <class NumberType>
105 {
107 };
108} // namespace internal
109
110
111
134{
146 template <int dim, int spacedim = dim>
147 class Scalar
148 {
149 public:
155 using value_type = double;
156
163
170
177
184 template <typename Number>
186
193 template <typename Number>
196
203 template <typename Number>
206
213 template <typename Number>
216
223 template <typename Number>
226
233 template <typename Number>
235 {
241 typename ProductType<Number,
242 typename Scalar<dim, spacedim>::value_type>::type;
243
248 using gradient_type = typename ProductType<
249 Number,
251
257 typename ProductType<Number,
258 typename Scalar<dim, spacedim>::value_type>::type;
259
264 using hessian_type = typename ProductType<
265 Number,
267
273 Number,
275 };
276
282 {
292
301 unsigned int row_index;
302 };
303
307 Scalar();
308
314 Scalar(const FEValuesBase<dim, spacedim> &fe_values_base,
315 const unsigned int component);
316
321 Scalar(const Scalar<dim, spacedim> &) = delete;
322
326 // NOLINTNEXTLINE OSX does not compile with noexcept
328
332 ~Scalar() = default;
333
338 Scalar &
340
344 Scalar &
345 operator=(Scalar<dim, spacedim> &&) noexcept = default;
346
361 value(const unsigned int shape_function, const unsigned int q_point) const;
362
374 gradient(const unsigned int shape_function,
375 const unsigned int q_point) const;
376
388 hessian(const unsigned int shape_function,
389 const unsigned int q_point) const;
390
402 third_derivative(const unsigned int shape_function,
403 const unsigned int q_point) const;
404
422 template <class InputVector>
423 void
425 const InputVector &fe_function,
426 std::vector<solution_value_type<typename InputVector::value_type>>
427 &values) const;
428
463 template <class InputVector>
464 void
466 const InputVector &dof_values,
467 std::vector<solution_value_type<typename InputVector::value_type>>
468 &values) const;
469
487 template <class InputVector>
488 void
490 const InputVector &fe_function,
491 std::vector<solution_gradient_type<typename InputVector::value_type>>
492 &gradients) const;
493
500 template <class InputVector>
501 void
503 const InputVector &dof_values,
504 std::vector<solution_gradient_type<typename InputVector::value_type>>
505 &gradients) const;
506
524 template <class InputVector>
525 void
527 const InputVector &fe_function,
528 std::vector<solution_hessian_type<typename InputVector::value_type>>
529 &hessians) const;
530
537 template <class InputVector>
538 void
540 const InputVector &dof_values,
541 std::vector<solution_hessian_type<typename InputVector::value_type>>
542 &hessians) const;
543
544
563 template <class InputVector>
564 void
566 const InputVector &fe_function,
567 std::vector<solution_laplacian_type<typename InputVector::value_type>>
568 &laplacians) const;
569
576 template <class InputVector>
577 void
579 const InputVector &dof_values,
580 std::vector<solution_laplacian_type<typename InputVector::value_type>>
581 &laplacians) const;
582
583
602 template <class InputVector>
603 void
605 const InputVector &fe_function,
606 std::vector<
607 solution_third_derivative_type<typename InputVector::value_type>>
608 &third_derivatives) const;
609
616 template <class InputVector>
617 void
619 const InputVector &dof_values,
620 std::vector<
621 solution_third_derivative_type<typename InputVector::value_type>>
622 &third_derivatives) const;
623
624
625 private:
629 const SmartPointer<const FEValuesBase<dim, spacedim>> fe_values;
630
635 const unsigned int component;
636
641 };
642
643
644
674 template <int dim, int spacedim = dim>
675 class Vector
676 {
677 public:
684
694
706
712 using divergence_type = double;
713
720 using curl_type = typename ::internal::CurlType<spacedim>::type;
721
728
735
742 template <typename Number>
744
751 template <typename Number>
754
761 template <typename Number>
764
771 template <typename Number>
774
781 template <typename Number>
784
791 template <typename Number>
793
800 template <typename Number>
803
810 template <typename Number>
813
820 template <typename Number>
822 {
828 typename ProductType<Number,
829 typename Vector<dim, spacedim>::value_type>::type;
830
835 using gradient_type = typename ProductType<
836 Number,
838
844 Number,
846
851 using divergence_type = typename ProductType<
852 Number,
854
860 typename ProductType<Number,
861 typename Vector<dim, spacedim>::value_type>::type;
862
867 using curl_type =
868 typename ProductType<Number,
869 typename Vector<dim, spacedim>::curl_type>::type;
870
875 using hessian_type = typename ProductType<
876 Number,
878
884 Number,
886 };
887
893 {
903
913 unsigned int row_index[spacedim];
914
925 };
926
930 Vector();
931
940 Vector(const FEValuesBase<dim, spacedim> &fe_values_base,
941 const unsigned int first_vector_component);
942
947 Vector(const Vector<dim, spacedim> &) = delete;
948
952 // NOLINTNEXTLINE OSX does not compile with noexcept
954
958 ~Vector() = default;
959
964 Vector &
966
970 // NOLINTNEXTLINE OSX does not compile with noexcept
971 Vector &
972 operator=(Vector<dim, spacedim> &&) = default; // NOLINT
973
991 value(const unsigned int shape_function, const unsigned int q_point) const;
992
1007 gradient(const unsigned int shape_function,
1008 const unsigned int q_point) const;
1009
1026 symmetric_gradient(const unsigned int shape_function,
1027 const unsigned int q_point) const;
1028
1040 divergence(const unsigned int shape_function,
1041 const unsigned int q_point) const;
1042
1063 curl_type
1064 curl(const unsigned int shape_function, const unsigned int q_point) const;
1065
1077 hessian(const unsigned int shape_function,
1078 const unsigned int q_point) const;
1079
1091 third_derivative(const unsigned int shape_function,
1092 const unsigned int q_point) const;
1093
1111 template <class InputVector>
1112 void
1114 const InputVector &fe_function,
1116 &values) const;
1117
1152 template <class InputVector>
1153 void
1155 const InputVector &dof_values,
1157 &values) const;
1158
1176 template <class InputVector>
1177 void
1179 const InputVector &fe_function,
1181 &gradients) const;
1182
1189 template <class InputVector>
1190 void
1192 const InputVector &dof_values,
1194 &gradients) const;
1195
1219 template <class InputVector>
1220 void
1221 get_function_symmetric_gradients(
1222 const InputVector &fe_function,
1223 std::vector<
1225 &symmetric_gradients) const;
1226
1233 template <class InputVector>
1234 void
1235 get_function_symmetric_gradients_from_local_dof_values(
1236 const InputVector &dof_values,
1237 std::vector<
1239 &symmetric_gradients) const;
1240
1259 template <class InputVector>
1260 void
1261 get_function_divergences(
1262 const InputVector &fe_function,
1264 &divergences) const;
1265
1272 template <class InputVector>
1273 void
1274 get_function_divergences_from_local_dof_values(
1275 const InputVector &dof_values,
1277 &divergences) const;
1278
1297 template <class InputVector>
1298 void
1299 get_function_curls(
1300 const InputVector &fe_function,
1302 const;
1303
1310 template <class InputVector>
1311 void
1312 get_function_curls_from_local_dof_values(
1313 const InputVector &dof_values,
1315 const;
1316
1334 template <class InputVector>
1335 void
1337 const InputVector &fe_function,
1339 &hessians) const;
1340
1347 template <class InputVector>
1348 void
1350 const InputVector &dof_values,
1352 &hessians) const;
1353
1372 template <class InputVector>
1373 void
1375 const InputVector &fe_function,
1377 &laplacians) const;
1378
1385 template <class InputVector>
1386 void
1388 const InputVector &dof_values,
1390 &laplacians) const;
1391
1410 template <class InputVector>
1411 void
1413 const InputVector &fe_function,
1414 std::vector<
1416 &third_derivatives) const;
1417
1424 template <class InputVector>
1425 void
1427 const InputVector &dof_values,
1428 std::vector<
1430 &third_derivatives) const;
1431
1432 private:
1437
1442 const unsigned int first_vector_component;
1443
1447 std::vector<ShapeFunctionData> shape_function_data;
1448 };
1449
1450
1451 template <int rank, int dim, int spacedim = dim>
1453
1476 template <int dim, int spacedim>
1477 class SymmetricTensor<2, dim, spacedim>
1478 {
1479 public:
1487
1498
1505 template <typename Number>
1507
1514 template <typename Number>
1517
1518
1525 template <typename Number>
1526 struct DEAL_II_DEPRECATED OutputType
1527 {
1532 using value_type = typename ProductType<
1533 Number,
1535
1541 Number,
1543 };
1544
1549 struct ShapeFunctionData
1550 {
1559 bool is_nonzero_shape_function_component
1560 [value_type::n_independent_components];
1561
1571 unsigned int row_index[value_type::n_independent_components];
1572
1582
1587 };
1588
1593
1603 SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1604 const unsigned int first_tensor_component);
1605
1611
1615 // NOLINTNEXTLINE OSX does not compile with noexcept
1617
1624
1630
1649 value(const unsigned int shape_function, const unsigned int q_point) const;
1650
1665 divergence(const unsigned int shape_function,
1666 const unsigned int q_point) const;
1667
1685 template <class InputVector>
1686 void
1687 get_function_values(
1688 const InputVector &fe_function,
1689 std::vector<solution_value_type<typename InputVector::value_type>>
1690 &values) const;
1691
1726 template <class InputVector>
1727 void
1728 get_function_values_from_local_dof_values(
1729 const InputVector &dof_values,
1730 std::vector<solution_value_type<typename InputVector::value_type>>
1731 &values) const;
1732
1754 template <class InputVector>
1755 void
1756 get_function_divergences(
1757 const InputVector &fe_function,
1758 std::vector<solution_divergence_type<typename InputVector::value_type>>
1759 &divergences) const;
1760
1767 template <class InputVector>
1768 void
1769 get_function_divergences_from_local_dof_values(
1770 const InputVector &dof_values,
1771 std::vector<solution_divergence_type<typename InputVector::value_type>>
1772 &divergences) const;
1773
1774 private:
1778 const SmartPointer<const FEValuesBase<dim, spacedim>> fe_values;
1779
1784 const unsigned int first_tensor_component;
1785
1789 std::vector<ShapeFunctionData> shape_function_data;
1790 };
1791
1792
1793 template <int rank, int dim, int spacedim = dim>
1794 class Tensor;
1795
1814 template <int dim, int spacedim>
1815 class Tensor<2, dim, spacedim>
1816 {
1817 public:
1823
1828
1834
1841 template <typename Number>
1843
1850 template <typename Number>
1853
1860 template <typename Number>
1863
1864
1871 template <typename Number>
1872 struct DEAL_II_DEPRECATED OutputType
1873 {
1878 using value_type = typename ProductType<
1879 Number,
1881
1887 Number,
1889
1894 using gradient_type = typename ProductType<
1895 Number,
1897 };
1898
1903 struct ShapeFunctionData
1904 {
1913 bool is_nonzero_shape_function_component
1914 [value_type::n_independent_components];
1915
1925 unsigned int row_index[value_type::n_independent_components];
1926
1936
1941 };
1942
1946 Tensor();
1947
1953
1957 // NOLINTNEXTLINE OSX does not compile with noexcept
1959
1963 ~Tensor() = default;
1964
1974 Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1975 const unsigned int first_tensor_component);
1976
1977
1982 Tensor &
1984
1988 Tensor &
1989 operator=(Tensor<2, dim, spacedim> &&) = default; // NOLINT
1990
2008 value(const unsigned int shape_function, const unsigned int q_point) const;
2009
2024 divergence(const unsigned int shape_function,
2025 const unsigned int q_point) const;
2026
2041 gradient(const unsigned int shape_function,
2042 const unsigned int q_point) const;
2043
2061 template <class InputVector>
2062 void
2063 get_function_values(
2064 const InputVector &fe_function,
2066 &values) const;
2067
2102 template <class InputVector>
2103 void
2104 get_function_values_from_local_dof_values(
2105 const InputVector &dof_values,
2107 &values) const;
2108
2130 template <class InputVector>
2131 void
2132 get_function_divergences(
2133 const InputVector &fe_function,
2135 &divergences) const;
2136
2143 template <class InputVector>
2144 void
2145 get_function_divergences_from_local_dof_values(
2146 const InputVector &dof_values,
2148 &divergences) const;
2149
2166 template <class InputVector>
2167 void
2168 get_function_gradients(
2169 const InputVector &fe_function,
2171 &gradients) const;
2172
2179 template <class InputVector>
2180 void
2181 get_function_gradients_from_local_dof_values(
2182 const InputVector &dof_values,
2184 &gradients) const;
2185
2186 private:
2191
2196 const unsigned int first_tensor_component;
2197
2201 std::vector<ShapeFunctionData> shape_function_data;
2202 };
2203
2204} // namespace FEValuesViews
2205
2206
2207namespace internal
2208{
2210 {
2215 template <int dim, int spacedim, typename Extractor>
2217 {};
2218
2226 template <int dim, int spacedim>
2227 struct ViewType<dim, spacedim, FEValuesExtractors::Scalar>
2228 {
2229 using type = typename ::FEValuesViews::Scalar<dim, spacedim>;
2230 };
2231
2239 template <int dim, int spacedim>
2240 struct ViewType<dim, spacedim, FEValuesExtractors::Vector>
2241 {
2242 using type = typename ::FEValuesViews::Vector<dim, spacedim>;
2243 };
2244
2252 template <int dim, int spacedim, int rank>
2253 struct ViewType<dim, spacedim, FEValuesExtractors::Tensor<rank>>
2254 {
2255 using type = typename ::FEValuesViews::Tensor<rank, dim, spacedim>;
2256 };
2257
2265 template <int dim, int spacedim, int rank>
2266 struct ViewType<dim, spacedim, FEValuesExtractors::SymmetricTensor<rank>>
2267 {
2268 using type =
2269 typename ::FEValuesViews::SymmetricTensor<rank, dim, spacedim>;
2270 };
2271
2279 template <int dim, int spacedim>
2280 struct Cache
2281 {
2286 std::vector<::FEValuesViews::Scalar<dim, spacedim>> scalars;
2287 std::vector<::FEValuesViews::Vector<dim, spacedim>> vectors;
2288 std::vector<::FEValuesViews::SymmetricTensor<2, dim, spacedim>>
2290 std::vector<::FEValuesViews::Tensor<2, dim, spacedim>>
2292
2296 Cache(const FEValuesBase<dim, spacedim> &fe_values);
2297 };
2298 } // namespace FEValuesViews
2299} // namespace internal
2300
2301namespace FEValuesViews
2302{
2307 template <int dim, int spacedim, typename Extractor>
2308 using View = typename ::internal::FEValuesViews::
2309 ViewType<dim, spacedim, Extractor>::type;
2310} // namespace FEValuesViews
2311
2312
2412template <int dim, int spacedim>
2414{
2415public:
2419 static constexpr unsigned int dimension = dim;
2420
2424 static constexpr unsigned int space_dimension = spacedim;
2425
2433 const unsigned int n_quadrature_points;
2434
2444 const unsigned int max_n_quadrature_points;
2445
2451 const unsigned int dofs_per_cell;
2452
2453
2461 FEValuesBase(const unsigned int n_q_points,
2462 const unsigned int dofs_per_cell,
2463 const UpdateFlags update_flags,
2464 const Mapping<dim, spacedim> & mapping,
2466
2471 FEValuesBase &
2472 operator=(const FEValuesBase &) = delete;
2473
2478 FEValuesBase(const FEValuesBase &) = delete;
2479
2483 virtual ~FEValuesBase() override;
2484
2485
2489
2511 const double &
2512 shape_value(const unsigned int i, const unsigned int q_point) const;
2513
2534 double
2535 shape_value_component(const unsigned int i,
2536 const unsigned int q_point,
2537 const unsigned int component) const;
2538
2564 const Tensor<1, spacedim> &
2565 shape_grad(const unsigned int i, const unsigned int q_point) const;
2566
2584 shape_grad_component(const unsigned int i,
2585 const unsigned int q_point,
2586 const unsigned int component) const;
2587
2607 const Tensor<2, spacedim> &
2608 shape_hessian(const unsigned int i, const unsigned int q_point) const;
2609
2627 shape_hessian_component(const unsigned int i,
2628 const unsigned int q_point,
2629 const unsigned int component) const;
2630
2650 const Tensor<3, spacedim> &
2651 shape_3rd_derivative(const unsigned int i, const unsigned int q_point) const;
2652
2671 const unsigned int q_point,
2672 const unsigned int component) const;
2673
2676
2728 template <class InputVector>
2729 void
2730 get_function_values(
2731 const InputVector & fe_function,
2732 std::vector<typename InputVector::value_type> &values) const;
2733
2747 template <class InputVector>
2748 void
2749 get_function_values(
2750 const InputVector & fe_function,
2751 std::vector<Vector<typename InputVector::value_type>> &values) const;
2752
2809 template <class InputVector>
2810 void
2811 get_function_values(
2812 const InputVector & fe_function,
2814 std::vector<typename InputVector::value_type> & values) const;
2815
2824 template <class InputVector>
2825 void
2826 get_function_values(
2827 const InputVector & fe_function,
2829 std::vector<Vector<typename InputVector::value_type>> &values) const;
2830
2831
2853 template <class InputVector>
2854 void
2855 get_function_values(
2856 const InputVector & fe_function,
2858 ArrayView<std::vector<typename InputVector::value_type>> values,
2859 const bool quadrature_points_fastest) const;
2860
2863
2907 template <class InputVector>
2908 void
2909 get_function_gradients(
2910 const InputVector &fe_function,
2912 &gradients) const;
2913
2930 template <class InputVector>
2931 void
2932 get_function_gradients(
2933 const InputVector &fe_function,
2934 std::vector<
2936 &gradients) const;
2937
2946 template <class InputVector>
2947 void
2948 get_function_gradients(
2949 const InputVector & fe_function,
2952 &gradients) const;
2953
2962 template <class InputVector>
2963 void
2964 get_function_gradients(
2965 const InputVector & fe_function,
2967 ArrayView<
2969 gradients,
2970 const bool quadrature_points_fastest = false) const;
2971
2976
3015 template <class InputVector>
3016 void
3017 get_function_hessians(
3018 const InputVector &fe_function,
3020 &hessians) const;
3021
3039 template <class InputVector>
3040 void
3041 get_function_hessians(
3042 const InputVector &fe_function,
3043 std::vector<
3045 & hessians,
3046 const bool quadrature_points_fastest = false) const;
3047
3056 template <class InputVector>
3057 void
3058 get_function_hessians(
3059 const InputVector & fe_function,
3062 &hessians) const;
3063
3072 template <class InputVector>
3073 void
3074 get_function_hessians(
3075 const InputVector & fe_function,
3077 ArrayView<
3079 hessians,
3080 const bool quadrature_points_fastest = false) const;
3081
3122 template <class InputVector>
3123 void
3124 get_function_laplacians(
3125 const InputVector & fe_function,
3126 std::vector<typename InputVector::value_type> &laplacians) const;
3127
3147 template <class InputVector>
3148 void
3149 get_function_laplacians(
3150 const InputVector & fe_function,
3151 std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
3152
3161 template <class InputVector>
3162 void
3163 get_function_laplacians(
3164 const InputVector & fe_function,
3166 std::vector<typename InputVector::value_type> & laplacians) const;
3167
3176 template <class InputVector>
3177 void
3178 get_function_laplacians(
3179 const InputVector & fe_function,
3181 std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
3182
3191 template <class InputVector>
3192 void
3193 get_function_laplacians(
3194 const InputVector & fe_function,
3196 std::vector<std::vector<typename InputVector::value_type>> &laplacians,
3197 const bool quadrature_points_fastest = false) const;
3198
3201
3241 template <class InputVector>
3242 void
3243 get_function_third_derivatives(
3244 const InputVector &fe_function,
3246 &third_derivatives) const;
3247
3266 template <class InputVector>
3267 void
3268 get_function_third_derivatives(
3269 const InputVector &fe_function,
3270 std::vector<
3272 & third_derivatives,
3273 const bool quadrature_points_fastest = false) const;
3274
3283 template <class InputVector>
3284 void
3285 get_function_third_derivatives(
3286 const InputVector & fe_function,
3289 &third_derivatives) const;
3290
3299 template <class InputVector>
3300 void
3301 get_function_third_derivatives(
3302 const InputVector & fe_function,
3304 ArrayView<
3306 third_derivatives,
3307 const bool quadrature_points_fastest = false) const;
3311
3338
3372 dof_indices_starting_at(const unsigned int start_dof_index) const;
3373
3405 dof_indices_ending_at(const unsigned int end_dof_index) const;
3406
3410
3434
3441 const Point<spacedim> &
3442 quadrature_point(const unsigned int q_point) const;
3443
3449 const std::vector<Point<spacedim>> &
3451
3468 double
3469 JxW(const unsigned int q_point) const;
3470
3474 const std::vector<double> &
3476
3484 jacobian(const unsigned int q_point) const;
3485
3492 const std::vector<DerivativeForm<1, dim, spacedim>> &
3494
3503 jacobian_grad(const unsigned int q_point) const;
3504
3511 const std::vector<DerivativeForm<2, dim, spacedim>> &
3513
3522 const Tensor<3, spacedim> &
3523 jacobian_pushed_forward_grad(const unsigned int q_point) const;
3524
3531 const std::vector<Tensor<3, spacedim>> &
3533
3542 jacobian_2nd_derivative(const unsigned int q_point) const;
3543
3550 const std::vector<DerivativeForm<3, dim, spacedim>> &
3552
3562 const Tensor<4, spacedim> &
3563 jacobian_pushed_forward_2nd_derivative(const unsigned int q_point) const;
3564
3571 const std::vector<Tensor<4, spacedim>> &
3573
3583 jacobian_3rd_derivative(const unsigned int q_point) const;
3584
3591 const std::vector<DerivativeForm<4, dim, spacedim>> &
3593
3603 const Tensor<5, spacedim> &
3604 jacobian_pushed_forward_3rd_derivative(const unsigned int q_point) const;
3605
3612 const std::vector<Tensor<5, spacedim>> &
3614
3622 inverse_jacobian(const unsigned int q_point) const;
3623
3630 const std::vector<DerivativeForm<1, spacedim, dim>> &
3632
3652 const Tensor<1, spacedim> &
3653 normal_vector(const unsigned int q_point) const;
3654
3662 const std::vector<Tensor<1, spacedim>> &
3663 get_normal_vectors() const;
3664
3668
3680
3691
3703
3704
3715
3719
3726
3731 get_fe() const;
3732
3738
3743 get_cell() const;
3744
3751 get_cell_similarity() const;
3752
3757 std::size_t
3758 memory_consumption() const;
3769 ExcAccessToUninitializedField,
3770 std::string,
3771 << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3772 << "object for which this kind of information has not been computed. What "
3773 << "information these objects compute is determined by the update_* flags you "
3774 << "pass to the constructor. Here, the operation you are attempting requires "
3775 << "the <" << arg1
3776 << "> flag to be set, but it was apparently not specified "
3777 << "upon construction.");
3778
3784 DeclExceptionMsg(ExcNotReinited,
3785 "FEValues object is not reinit'ed to any cell");
3786
3794 ExcFEDontMatch,
3795 "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3796 "to the DoFHandler that provided the cell iterator do not match.");
3802 DeclException1(ExcShapeFunctionNotPrimitive,
3803 int,
3804 << "The shape function with index " << arg1
3805 << " is not primitive, i.e. it is vector-valued and "
3806 << "has more than one non-zero vector component. This "
3807 << "function cannot be called for these shape functions. "
3808 << "Maybe you want to use the same function with the "
3809 << "_component suffix?");
3810
3818 ExcFENotPrimitive,
3819 "The given FiniteElement is not a primitive element but the requested operation "
3820 "only works for those. See FiniteElement::is_primitive() for more information.");
3821
3822protected:
3830 {
3831 public:
3833 ExcNeedsDoFHandler,
3834 "You have previously called the FEValues::reinit() function with a "
3835 "cell iterator of type Triangulation<dim,spacedim>::cell_iterator. However, "
3836 "when you do this, you cannot call some functions in the FEValues "
3837 "class, such as the get_function_values/gradients/hessians/third_derivatives "
3838 "functions. If you need these functions, then you need to call "
3839 "FEValues::reinit() with an iterator type that allows to extract "
3840 "degrees of freedom, such as DoFHandler<dim,spacedim>::cell_iterator.");
3841
3846
3850 template <bool lda>
3853
3858 const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3859
3863 bool
3864 is_initialized() const;
3865
3872 operator typename Triangulation<dim, spacedim>::cell_iterator() const;
3873
3879 n_dofs_for_dof_handler() const;
3880
3885 template <typename VectorType>
3886 void
3887 get_interpolated_dof_values(
3888 const VectorType & in,
3890
3895 void
3896 get_interpolated_dof_values(const IndexSet & in,
3897 Vector<IndexSet::value_type> &out) const;
3898
3899 private:
3904 };
3905
3912
3920 boost::signals2::connection tria_listener_refinement;
3921
3929 boost::signals2::connection tria_listener_mesh_transform;
3930
3936 void
3937 invalidate_present_cell();
3938
3948 void
3949 maybe_invalidate_previous_present_cell(
3950 const typename Triangulation<dim, spacedim>::cell_iterator &cell);
3951
3957
3963 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3965
3972
3980
3986 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3988
3994 spacedim>
3996
3997
4002
4012 compute_update_flags(const UpdateFlags update_flags) const;
4013
4020
4026 void
4027 check_cell_similarity(
4028 const typename Triangulation<dim, spacedim>::cell_iterator &cell);
4029
4030private:
4035
4036 // Make the view classes friends of this class, since they access internal
4037 // data.
4038 template <int, int>
4040 template <int, int>
4042 template <int, int, int>
4044 template <int, int, int>
4046};
4047
4048
4049
4059template <int dim, int spacedim = dim>
4060class FEValues : public FEValuesBase<dim, spacedim>
4061{
4062public:
4067 static constexpr unsigned int integral_dimension = dim;
4068
4075 const Quadrature<dim> & quadrature,
4076 const UpdateFlags update_flags);
4077
4086 const hp::QCollection<dim> & quadrature,
4087 const UpdateFlags update_flags);
4088
4095 const Quadrature<dim> & quadrature,
4096 const UpdateFlags update_flags);
4097
4105 const hp::QCollection<dim> & quadrature,
4106 const UpdateFlags update_flags);
4107
4114 template <bool level_dof_access>
4115 void
4118
4132 void
4134
4139 const Quadrature<dim> &
4141
4146 std::size_t
4148
4165
4166private:
4171
4175 void
4176 initialize(const UpdateFlags update_flags);
4177
4184 void
4186};
4187
4188
4198template <int dim, int spacedim = dim>
4199class FEFaceValuesBase : public FEValuesBase<dim, spacedim>
4200{
4201public:
4206 static constexpr unsigned int integral_dimension = dim - 1;
4207
4219 FEFaceValuesBase(const unsigned int dofs_per_cell,
4220 const UpdateFlags update_flags,
4221 const Mapping<dim, spacedim> & mapping,
4223 const Quadrature<dim - 1> & quadrature);
4224
4231 FEFaceValuesBase(const unsigned int dofs_per_cell,
4232 const UpdateFlags update_flags,
4233 const Mapping<dim, spacedim> & mapping,
4235 const hp::QCollection<dim - 1> & quadrature);
4236
4244 const Tensor<1, spacedim> &
4245 boundary_form(const unsigned int q_point) const;
4246
4253 const std::vector<Tensor<1, spacedim>> &
4255
4260 unsigned int
4262
4267 unsigned int
4269
4274 const Quadrature<dim - 1> &
4276
4281 std::size_t
4283
4284protected:
4289 unsigned int present_face_no;
4290
4296
4301};
4302
4303
4304
4318template <int dim, int spacedim = dim>
4319class FEFaceValues : public FEFaceValuesBase<dim, spacedim>
4320{
4321public:
4326 static constexpr unsigned int dimension = dim;
4327
4328 static constexpr unsigned int space_dimension = spacedim;
4329
4334 static constexpr unsigned int integral_dimension = dim - 1;
4335
4342 const Quadrature<dim - 1> & quadrature,
4343 const UpdateFlags update_flags);
4344
4353 const hp::QCollection<dim - 1> & quadrature,
4354 const UpdateFlags update_flags);
4355
4362 const Quadrature<dim - 1> & quadrature,
4363 const UpdateFlags update_flags);
4364
4372 const hp::QCollection<dim - 1> & quadrature,
4373 const UpdateFlags update_flags);
4374
4379 template <bool level_dof_access>
4380 void
4383 const unsigned int face_no);
4384
4391 template <bool level_dof_access>
4392 void
4395 const typename Triangulation<dim, spacedim>::face_iterator & face);
4396
4410 void
4412 const unsigned int face_no);
4413
4414 /*
4415 * Reinitialize the gradients, Jacobi determinants, etc for the given face
4416 * on a given cell of type "iterator into a Triangulation object", and the
4417 * given finite element. Since iterators into a triangulation alone only
4418 * convey information about the geometry of a cell, but not about degrees of
4419 * freedom possibly associated with this cell, you will not be able to call
4420 * some functions of this class if they need information about degrees of
4421 * freedom. These functions are, above all, the
4422 * <tt>get_function_value/gradients/hessians/third_derivatives</tt>
4423 * functions. If you want to call these functions, you have to call the @p
4424 * reinit variants that take iterators into DoFHandler or other DoF handler
4425 * type objects.
4426 *
4427 * @note @p face must be one of @p cell's face iterators.
4428 */
4429 void
4431 const typename Triangulation<dim, spacedim>::face_iterator &face);
4432
4449
4450private:
4454 void
4455 initialize(const UpdateFlags update_flags);
4456
4463 void
4464 do_reinit(const unsigned int face_no);
4465};
4466
4467
4484template <int dim, int spacedim = dim>
4485class FESubfaceValues : public FEFaceValuesBase<dim, spacedim>
4486{
4487public:
4491 static constexpr unsigned int dimension = dim;
4492
4496 static constexpr unsigned int space_dimension = spacedim;
4497
4502 static constexpr unsigned int integral_dimension = dim - 1;
4503
4510 const Quadrature<dim - 1> & face_quadrature,
4511 const UpdateFlags update_flags);
4512
4521 const hp::QCollection<dim - 1> & face_quadrature,
4522 const UpdateFlags update_flags);
4523
4530 const Quadrature<dim - 1> & face_quadrature,
4531 const UpdateFlags update_flags);
4532
4540 const hp::QCollection<dim - 1> & face_quadrature,
4541 const UpdateFlags update_flags);
4542
4549 template <bool level_dof_access>
4550 void
4553 const unsigned int face_no,
4554 const unsigned int subface_no);
4555
4560 template <bool level_dof_access>
4561 void
4564 const typename Triangulation<dim, spacedim>::face_iterator & face,
4565 const typename Triangulation<dim, spacedim>::face_iterator &subface);
4566
4580 void
4582 const unsigned int face_no,
4583 const unsigned int subface_no);
4584
4604 void
4607 const typename Triangulation<dim, spacedim>::face_iterator &subface);
4608
4625
4631 DeclException0(ExcReinitCalledWithBoundaryFace);
4632
4638 DeclException0(ExcFaceHasNoSubfaces);
4639
4640private:
4644 void
4645 initialize(const UpdateFlags update_flags);
4646
4653 void
4654 do_reinit(const unsigned int face_no, const unsigned int subface_no);
4655};
4656
4657
4658#ifndef DOXYGEN
4659
4660
4661/*------------------------ Inline functions: namespace FEValuesViews --------*/
4662
4663namespace FEValuesViews
4664{
4665 template <int dim, int spacedim>
4666 inline typename Scalar<dim, spacedim>::value_type
4667 Scalar<dim, spacedim>::value(const unsigned int shape_function,
4668 const unsigned int q_point) const
4669 {
4670 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4671 Assert(
4672 fe_values->update_flags & update_values,
4674 "update_values"))));
4675
4676 // an adaptation of the FEValuesBase::shape_value_component function
4677 // except that here we know the component as fixed and we have
4678 // pre-computed and cached a bunch of information. See the comments there.
4679 if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4680 return fe_values->finite_element_output.shape_values(
4681 shape_function_data[shape_function].row_index, q_point);
4682 else
4683 return 0;
4684 }
4685
4686
4687
4688 template <int dim, int spacedim>
4689 inline typename Scalar<dim, spacedim>::gradient_type
4690 Scalar<dim, spacedim>::gradient(const unsigned int shape_function,
4691 const unsigned int q_point) const
4692 {
4693 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4694 Assert(fe_values->update_flags & update_gradients,
4696 "update_gradients")));
4697
4698 // an adaptation of the FEValuesBase::shape_grad_component
4699 // function except that here we know the component as fixed and we have
4700 // pre-computed and cached a bunch of information. See the comments there.
4701 if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4702 return fe_values->finite_element_output
4703 .shape_gradients[shape_function_data[shape_function].row_index]
4704 [q_point];
4705 else
4706 return gradient_type();
4707 }
4708
4709
4710
4711 template <int dim, int spacedim>
4712 inline typename Scalar<dim, spacedim>::hessian_type
4713 Scalar<dim, spacedim>::hessian(const unsigned int shape_function,
4714 const unsigned int q_point) const
4715 {
4716 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4717 Assert(fe_values->update_flags & update_hessians,
4719 "update_hessians")));
4720
4721 // an adaptation of the FEValuesBase::shape_hessian_component
4722 // function except that here we know the component as fixed and we have
4723 // pre-computed and cached a bunch of information. See the comments there.
4724 if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4725 return fe_values->finite_element_output
4726 .shape_hessians[shape_function_data[shape_function].row_index][q_point];
4727 else
4728 return hessian_type();
4729 }
4730
4731
4732
4733 template <int dim, int spacedim>
4734 inline typename Scalar<dim, spacedim>::third_derivative_type
4735 Scalar<dim, spacedim>::third_derivative(const unsigned int shape_function,
4736 const unsigned int q_point) const
4737 {
4738 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4739 Assert(fe_values->update_flags & update_3rd_derivatives,
4741 "update_3rd_derivatives")));
4742
4743 // an adaptation of the FEValuesBase::shape_3rdderivative_component
4744 // function except that here we know the component as fixed and we have
4745 // pre-computed and cached a bunch of information. See the comments there.
4746 if (shape_function_data[shape_function].is_nonzero_shape_function_component)
4747 return fe_values->finite_element_output
4748 .shape_3rd_derivatives[shape_function_data[shape_function].row_index]
4749 [q_point];
4750 else
4751 return third_derivative_type();
4752 }
4753
4754
4755
4756 template <int dim, int spacedim>
4757 inline typename Vector<dim, spacedim>::value_type
4758 Vector<dim, spacedim>::value(const unsigned int shape_function,
4759 const unsigned int q_point) const
4760 {
4761 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4762 Assert(fe_values->update_flags & update_values,
4764 "update_values")));
4765
4766 // same as for the scalar case except that we have one more index
4767 const int snc =
4768 shape_function_data[shape_function].single_nonzero_component;
4769 if (snc == -2)
4770 return value_type();
4771 else if (snc != -1)
4772 {
4773 value_type return_value;
4774 return_value[shape_function_data[shape_function]
4775 .single_nonzero_component_index] =
4776 fe_values->finite_element_output.shape_values(snc, q_point);
4777 return return_value;
4778 }
4779 else
4780 {
4781 value_type return_value;
4782 for (unsigned int d = 0; d < dim; ++d)
4783 if (shape_function_data[shape_function]
4784 .is_nonzero_shape_function_component[d])
4785 return_value[d] = fe_values->finite_element_output.shape_values(
4786 shape_function_data[shape_function].row_index[d], q_point);
4787
4788 return return_value;
4789 }
4790 }
4791
4792
4793
4794 template <int dim, int spacedim>
4796 Vector<dim, spacedim>::gradient(const unsigned int shape_function,
4797 const unsigned int q_point) const
4798 {
4799 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4800 Assert(fe_values->update_flags & update_gradients,
4802 "update_gradients")));
4803
4804 // same as for the scalar case except that we have one more index
4805 const int snc =
4806 shape_function_data[shape_function].single_nonzero_component;
4807 if (snc == -2)
4808 return gradient_type();
4809 else if (snc != -1)
4810 {
4811 gradient_type return_value;
4812 return_value[shape_function_data[shape_function]
4813 .single_nonzero_component_index] =
4814 fe_values->finite_element_output.shape_gradients[snc][q_point];
4815 return return_value;
4816 }
4817 else
4818 {
4819 gradient_type return_value;
4820 for (unsigned int d = 0; d < dim; ++d)
4821 if (shape_function_data[shape_function]
4822 .is_nonzero_shape_function_component[d])
4823 return_value[d] =
4824 fe_values->finite_element_output.shape_gradients
4825 [shape_function_data[shape_function].row_index[d]][q_point];
4826
4827 return return_value;
4828 }
4829 }
4830
4831
4832
4833 template <int dim, int spacedim>
4835 Vector<dim, spacedim>::divergence(const unsigned int shape_function,
4836 const unsigned int q_point) const
4837 {
4838 // this function works like in the case above
4839 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4840 Assert(fe_values->update_flags & update_gradients,
4842 "update_gradients")));
4843
4844 // same as for the scalar case except that we have one more index
4845 const int snc =
4846 shape_function_data[shape_function].single_nonzero_component;
4847 if (snc == -2)
4848 return divergence_type();
4849 else if (snc != -1)
4850 return fe_values->finite_element_output
4851 .shape_gradients[snc][q_point][shape_function_data[shape_function]
4852 .single_nonzero_component_index];
4853 else
4854 {
4855 divergence_type return_value = 0;
4856 for (unsigned int d = 0; d < dim; ++d)
4857 if (shape_function_data[shape_function]
4858 .is_nonzero_shape_function_component[d])
4859 return_value +=
4860 fe_values->finite_element_output.shape_gradients
4861 [shape_function_data[shape_function].row_index[d]][q_point][d];
4862
4863 return return_value;
4864 }
4865 }
4866
4867
4868
4869 template <int dim, int spacedim>
4870 inline typename Vector<dim, spacedim>::curl_type
4871 Vector<dim, spacedim>::curl(const unsigned int shape_function,
4872 const unsigned int q_point) const
4873 {
4874 // this function works like in the case above
4875
4876 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
4877 Assert(fe_values->update_flags & update_gradients,
4879 "update_gradients")));
4880 // same as for the scalar case except that we have one more index
4881 const int snc =
4882 shape_function_data[shape_function].single_nonzero_component;
4883
4884 if (snc == -2)
4885 return curl_type();
4886
4887 else
4888 switch (dim)
4889 {
4890 case 1:
4891 {
4892 Assert(false,
4893 ExcMessage(
4894 "Computing the curl in 1d is not a useful operation"));
4895 return curl_type();
4896 }
4897
4898 case 2:
4899 {
4900 if (snc != -1)
4901 {
4902 curl_type return_value;
4903
4904 // the single nonzero component can only be zero or one in 2d
4905 if (shape_function_data[shape_function]
4906 .single_nonzero_component_index == 0)
4907 return_value[0] =
4908 -1.0 * fe_values->finite_element_output
4909 .shape_gradients[snc][q_point][1];
4910 else
4911 return_value[0] = fe_values->finite_element_output
4912 .shape_gradients[snc][q_point][0];
4913
4914 return return_value;
4915 }
4916
4917 else
4918 {
4919 curl_type return_value;
4920
4921 return_value[0] = 0.0;
4922
4923 if (shape_function_data[shape_function]
4924 .is_nonzero_shape_function_component[0])
4925 return_value[0] -=
4926 fe_values->finite_element_output
4927 .shape_gradients[shape_function_data[shape_function]
4928 .row_index[0]][q_point][1];
4929
4930 if (shape_function_data[shape_function]
4931 .is_nonzero_shape_function_component[1])
4932 return_value[0] +=
4933 fe_values->finite_element_output
4934 .shape_gradients[shape_function_data[shape_function]
4935 .row_index[1]][q_point][0];
4936
4937 return return_value;
4938 }
4939 }
4940
4941 case 3:
4942 {
4943 if (snc != -1)
4944 {
4945 curl_type return_value;
4946
4947 switch (shape_function_data[shape_function]
4948 .single_nonzero_component_index)
4949 {
4950 case 0:
4951 {
4952 return_value[0] = 0;
4953 return_value[1] = fe_values->finite_element_output
4954 .shape_gradients[snc][q_point][2];
4955 return_value[2] =
4956 -1.0 * fe_values->finite_element_output
4957 .shape_gradients[snc][q_point][1];
4958 return return_value;
4959 }
4960
4961 case 1:
4962 {
4963 return_value[0] =
4964 -1.0 * fe_values->finite_element_output
4965 .shape_gradients[snc][q_point][2];
4966 return_value[1] = 0;
4967 return_value[2] = fe_values->finite_element_output
4968 .shape_gradients[snc][q_point][0];
4969 return return_value;
4970 }
4971
4972 default:
4973 {
4974 return_value[0] = fe_values->finite_element_output
4975 .shape_gradients[snc][q_point][1];
4976 return_value[1] =
4977 -1.0 * fe_values->finite_element_output
4978 .shape_gradients[snc][q_point][0];
4979 return_value[2] = 0;
4980 return return_value;
4981 }
4982 }
4983 }
4984
4985 else
4986 {
4987 curl_type return_value;
4988
4989 for (unsigned int i = 0; i < dim; ++i)
4990 return_value[i] = 0.0;
4991
4992 if (shape_function_data[shape_function]
4993 .is_nonzero_shape_function_component[0])
4994 {
4995 return_value[1] +=
4996 fe_values->finite_element_output
4997 .shape_gradients[shape_function_data[shape_function]
4998 .row_index[0]][q_point][2];
4999 return_value[2] -=
5000 fe_values->finite_element_output
5001 .shape_gradients[shape_function_data[shape_function]
5002 .row_index[0]][q_point][1];
5003 }
5004
5005 if (shape_function_data[shape_function]
5006 .is_nonzero_shape_function_component[1])
5007 {
5008 return_value[0] -=
5009 fe_values->finite_element_output
5010 .shape_gradients[shape_function_data[shape_function]
5011 .row_index[1]][q_point][2];
5012 return_value[2] +=
5013 fe_values->finite_element_output
5014 .shape_gradients[shape_function_data[shape_function]
5015 .row_index[1]][q_point][0];
5016 }
5017
5018 if (shape_function_data[shape_function]
5019 .is_nonzero_shape_function_component[2])
5020 {
5021 return_value[0] +=
5022 fe_values->finite_element_output
5023 .shape_gradients[shape_function_data[shape_function]
5024 .row_index[2]][q_point][1];
5025 return_value[1] -=
5026 fe_values->finite_element_output
5027 .shape_gradients[shape_function_data[shape_function]
5028 .row_index[2]][q_point][0];
5029 }
5030
5031 return return_value;
5032 }
5033 }
5034 }
5035 // should not end up here
5036 Assert(false, ExcInternalError());
5037 return curl_type();
5038 }
5039
5040
5041
5042 template <int dim, int spacedim>
5044 Vector<dim, spacedim>::hessian(const unsigned int shape_function,
5045 const unsigned int q_point) const
5046 {
5047 // this function works like in the case above
5048 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5049 Assert(fe_values->update_flags & update_hessians,
5051 "update_hessians")));
5052
5053 // same as for the scalar case except that we have one more index
5054 const int snc =
5055 shape_function_data[shape_function].single_nonzero_component;
5056 if (snc == -2)
5057 return hessian_type();
5058 else if (snc != -1)
5059 {
5060 hessian_type return_value;
5061 return_value[shape_function_data[shape_function]
5062 .single_nonzero_component_index] =
5063 fe_values->finite_element_output.shape_hessians[snc][q_point];
5064 return return_value;
5065 }
5066 else
5067 {
5068 hessian_type return_value;
5069 for (unsigned int d = 0; d < dim; ++d)
5070 if (shape_function_data[shape_function]
5071 .is_nonzero_shape_function_component[d])
5072 return_value[d] =
5073 fe_values->finite_element_output.shape_hessians
5074 [shape_function_data[shape_function].row_index[d]][q_point];
5075
5076 return return_value;
5077 }
5078 }
5079
5080
5081
5082 template <int dim, int spacedim>
5084 Vector<dim, spacedim>::third_derivative(const unsigned int shape_function,
5085 const unsigned int q_point) const
5086 {
5087 // this function works like in the case above
5088 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5089 Assert(fe_values->update_flags & update_3rd_derivatives,
5091 "update_3rd_derivatives")));
5092
5093 // same as for the scalar case except that we have one more index
5094 const int snc =
5095 shape_function_data[shape_function].single_nonzero_component;
5096 if (snc == -2)
5097 return third_derivative_type();
5098 else if (snc != -1)
5099 {
5100 third_derivative_type return_value;
5101 return_value[shape_function_data[shape_function]
5102 .single_nonzero_component_index] =
5103 fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
5104 return return_value;
5105 }
5106 else
5107 {
5108 third_derivative_type return_value;
5109 for (unsigned int d = 0; d < dim; ++d)
5110 if (shape_function_data[shape_function]
5111 .is_nonzero_shape_function_component[d])
5112 return_value[d] =
5113 fe_values->finite_element_output.shape_3rd_derivatives
5114 [shape_function_data[shape_function].row_index[d]][q_point];
5115
5116 return return_value;
5117 }
5118 }
5119
5120
5121
5122 namespace internal
5123 {
5128 inline ::SymmetricTensor<2, 1>
5129 symmetrize_single_row(const unsigned int n, const ::Tensor<1, 1> &t)
5130 {
5131 AssertIndexRange(n, 1);
5132 (void)n;
5133
5134 return {{t[0]}};
5135 }
5136
5137
5138
5139 inline ::SymmetricTensor<2, 2>
5140 symmetrize_single_row(const unsigned int n, const ::Tensor<1, 2> &t)
5141 {
5142 switch (n)
5143 {
5144 case 0:
5145 {
5146 return {{t[0], 0, t[1] / 2}};
5147 }
5148 case 1:
5149 {
5150 return {{0, t[1], t[0] / 2}};
5151 }
5152 default:
5153 {
5154 AssertIndexRange(n, 2);
5155 return {};
5156 }
5157 }
5158 }
5159
5160
5161
5162 inline ::SymmetricTensor<2, 3>
5163 symmetrize_single_row(const unsigned int n, const ::Tensor<1, 3> &t)
5164 {
5165 switch (n)
5166 {
5167 case 0:
5168 {
5169 return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
5170 }
5171 case 1:
5172 {
5173 return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
5174 }
5175 case 2:
5176 {
5177 return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
5178 }
5179 default:
5180 {
5181 AssertIndexRange(n, 3);
5182 return {};
5183 }
5184 }
5185 }
5186 } // namespace internal
5187
5188
5189
5190 template <int dim, int spacedim>
5192 Vector<dim, spacedim>::symmetric_gradient(const unsigned int shape_function,
5193 const unsigned int q_point) const
5194 {
5195 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5196 Assert(fe_values->update_flags & update_gradients,
5198 "update_gradients")));
5199
5200 // same as for the scalar case except that we have one more index
5201 const int snc =
5202 shape_function_data[shape_function].single_nonzero_component;
5203 if (snc == -2)
5204 return symmetric_gradient_type();
5205 else if (snc != -1)
5206 return internal::symmetrize_single_row(
5207 shape_function_data[shape_function].single_nonzero_component_index,
5208 fe_values->finite_element_output.shape_gradients[snc][q_point]);
5209 else
5210 {
5211 gradient_type return_value;
5212 for (unsigned int d = 0; d < dim; ++d)
5213 if (shape_function_data[shape_function]
5214 .is_nonzero_shape_function_component[d])
5215 return_value[d] =
5216 fe_values->finite_element_output.shape_gradients
5217 [shape_function_data[shape_function].row_index[d]][q_point];
5218
5219 return symmetrize(return_value);
5220 }
5221 }
5222
5223
5224
5225 template <int dim, int spacedim>
5227 SymmetricTensor<2, dim, spacedim>::value(const unsigned int shape_function,
5228 const unsigned int q_point) const
5229 {
5230 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5231 Assert(fe_values->update_flags & update_values,
5233 "update_values")));
5234
5235 // similar to the vector case where we have more then one index and we need
5236 // to convert between unrolled and component indexing for tensors
5237 const int snc =
5238 shape_function_data[shape_function].single_nonzero_component;
5239
5240 if (snc == -2)
5241 {
5242 // shape function is zero for the selected components
5243 return value_type();
5244 }
5245 else if (snc != -1)
5246 {
5247 value_type return_value;
5248 const unsigned int comp =
5249 shape_function_data[shape_function].single_nonzero_component_index;
5250 return_value[value_type::unrolled_to_component_indices(comp)] =
5251 fe_values->finite_element_output.shape_values(snc, q_point);
5252 return return_value;
5253 }
5254 else
5255 {
5256 value_type return_value;
5257 for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
5258 if (shape_function_data[shape_function]
5259 .is_nonzero_shape_function_component[d])
5260 return_value[value_type::unrolled_to_component_indices(d)] =
5261 fe_values->finite_element_output.shape_values(
5262 shape_function_data[shape_function].row_index[d], q_point);
5263 return return_value;
5264 }
5265 }
5266
5267
5268
5269 template <int dim, int spacedim>
5272 const unsigned int shape_function,
5273 const unsigned int q_point) const
5274 {
5275 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5276 Assert(fe_values->update_flags & update_gradients,
5278 "update_gradients")));
5279
5280 const int snc =
5281 shape_function_data[shape_function].single_nonzero_component;
5282
5283 if (snc == -2)
5284 {
5285 // shape function is zero for the selected components
5286 return divergence_type();
5287 }
5288 else if (snc != -1)
5289 {
5290 // we have a single non-zero component when the symmetric tensor is
5291 // represented in unrolled form. this implies we potentially have
5292 // two non-zero components when represented in component form! we
5293 // will only have one non-zero entry if the non-zero component lies on
5294 // the diagonal of the tensor.
5295 //
5296 // the divergence of a second-order tensor is a first order tensor.
5297 //
5298 // assume the second-order tensor is A with components A_{ij}. then
5299 // A_{ij} = A_{ji} and there is only one (if diagonal) or two non-zero
5300 // entries in the tensorial representation. define the
5301 // divergence as:
5302 // b_i \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_j}.
5303 // (which is incidentally also
5304 // b_j \dealcoloneq \dfrac{\partial phi_{ij}}{\partial x_i}).
5305 // In both cases, a sum is implied.
5306 //
5307 // Now, we know the nonzero component in unrolled form: it is indicated
5308 // by 'snc'. we can figure out which tensor components belong to this:
5309 const unsigned int comp =
5310 shape_function_data[shape_function].single_nonzero_component_index;
5311 const unsigned int ii =
5312 value_type::unrolled_to_component_indices(comp)[0];
5313 const unsigned int jj =
5314 value_type::unrolled_to_component_indices(comp)[1];
5315
5316 // given the form of the divergence above, if ii=jj there is only a
5317 // single nonzero component of the full tensor and the gradient
5318 // equals
5319 // b_ii \dealcoloneq \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
5320 // all other entries of 'b' are zero
5321 //
5322 // on the other hand, if ii!=jj, then there are two nonzero entries in
5323 // the full tensor and
5324 // b_ii \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
5325 // b_jj \dealcoloneq \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
5326 // again, all other entries of 'b' are zero
5327 const ::Tensor<1, spacedim> &phi_grad =
5328 fe_values->finite_element_output.shape_gradients[snc][q_point];
5329
5330 divergence_type return_value;
5331 return_value[ii] = phi_grad[jj];
5332
5333 if (ii != jj)
5334 return_value[jj] = phi_grad[ii];
5335
5336 return return_value;
5337 }
5338 else
5339 {
5340 Assert(false, ExcNotImplemented());
5341 divergence_type return_value;
5342 return return_value;
5343 }
5344 }
5345
5346
5347
5348 template <int dim, int spacedim>
5350 Tensor<2, dim, spacedim>::value(const unsigned int shape_function,
5351 const unsigned int q_point) const
5352 {
5353 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5354 Assert(fe_values->update_flags & update_values,
5356 "update_values")));
5357
5358 // similar to the vector case where we have more then one index and we need
5359 // to convert between unrolled and component indexing for tensors
5360 const int snc =
5361 shape_function_data[shape_function].single_nonzero_component;
5362
5363 if (snc == -2)
5364 {
5365 // shape function is zero for the selected components
5366 return value_type();
5367 }
5368 else if (snc != -1)
5369 {
5370 value_type return_value;
5371 const unsigned int comp =
5372 shape_function_data[shape_function].single_nonzero_component_index;
5373 const TableIndices<2> indices =
5375 return_value[indices] =
5376 fe_values->finite_element_output.shape_values(snc, q_point);
5377 return return_value;
5378 }
5379 else
5380 {
5381 value_type return_value;
5382 for (unsigned int d = 0; d < dim * dim; ++d)
5383 if (shape_function_data[shape_function]
5384 .is_nonzero_shape_function_component[d])
5385 {
5386 const TableIndices<2> indices =
5388 return_value[indices] =
5389 fe_values->finite_element_output.shape_values(
5390 shape_function_data[shape_function].row_index[d], q_point);
5391 }
5392 return return_value;
5393 }
5394 }
5395
5396
5397
5398 template <int dim, int spacedim>
5400 Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
5401 const unsigned int q_point) const
5402 {
5403 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5404 Assert(fe_values->update_flags & update_gradients,
5406 "update_gradients")));
5407
5408 const int snc =
5409 shape_function_data[shape_function].single_nonzero_component;
5410
5411 if (snc == -2)
5412 {
5413 // shape function is zero for the selected components
5414 return divergence_type();
5415 }
5416 else if (snc != -1)
5417 {
5418 // we have a single non-zero component when the tensor is
5419 // represented in unrolled form.
5420 //
5421 // the divergence of a second-order tensor is a first order tensor.
5422 //
5423 // assume the second-order tensor is A with components A_{ij},
5424 // then divergence is d_i := \frac{\partial A_{ij}}{\partial x_j}
5425 //
5426 // Now, we know the nonzero component in unrolled form: it is indicated
5427 // by 'snc'. we can figure out which tensor components belong to this:
5428 const unsigned int comp =
5429 shape_function_data[shape_function].single_nonzero_component_index;
5430 const TableIndices<2> indices =
5432 const unsigned int ii = indices[0];
5433 const unsigned int jj = indices[1];
5434
5435 const ::Tensor<1, spacedim> &phi_grad =
5436 fe_values->finite_element_output.shape_gradients[snc][q_point];
5437
5438 divergence_type return_value;
5439 // note that we contract \nabla from the right
5440 return_value[ii] = phi_grad[jj];
5441
5442 return return_value;
5443 }
5444 else
5445 {
5446 Assert(false, ExcNotImplemented());
5447 divergence_type return_value;
5448 return return_value;
5449 }
5450 }
5451
5452
5453
5454 template <int dim, int spacedim>
5456 Tensor<2, dim, spacedim>::gradient(const unsigned int shape_function,
5457 const unsigned int q_point) const
5458 {
5459 AssertIndexRange(shape_function, fe_values->fe->n_dofs_per_cell());
5460 Assert(fe_values->update_flags & update_gradients,
5462 "update_gradients")));
5463
5464 const int snc =
5465 shape_function_data[shape_function].single_nonzero_component;
5466
5467 if (snc == -2)
5468 {
5469 // shape function is zero for the selected components
5470 return gradient_type();
5471 }
5472 else if (snc != -1)
5473 {
5474 // we have a single non-zero component when the tensor is
5475 // represented in unrolled form.
5476 //
5477 // the gradient of a second-order tensor is a third order tensor.
5478 //
5479 // assume the second-order tensor is A with components A_{ij},
5480 // then gradient is B_{ijk} := \frac{\partial A_{ij}}{\partial x_k}
5481 //
5482 // Now, we know the nonzero component in unrolled form: it is indicated
5483 // by 'snc'. we can figure out which tensor components belong to this:
5484 const unsigned int comp =
5485 shape_function_data[shape_function].single_nonzero_component_index;
5486 const TableIndices<2> indices =
5488 const unsigned int ii = indices[0];
5489 const unsigned int jj = indices[1];
5490
5491 const ::Tensor<1, spacedim> &phi_grad =
5492 fe_values->finite_element_output.shape_gradients[snc][q_point];
5493
5494 gradient_type return_value;
5495 return_value[ii][jj] = phi_grad;
5496
5497 return return_value;
5498 }
5499 else
5500 {
5501 Assert(false, ExcNotImplemented());
5502 gradient_type return_value;
5503 return return_value;
5504 }
5505 }
5506
5507} // namespace FEValuesViews
5508
5509
5510
5511/*---------------------- Inline functions: FEValuesBase ---------------------*/
5512
5513
5514
5515template <int dim, int spacedim>
5516template <bool lda>
5520 : initialized(true)
5521 , cell(cell)
5522 , dof_handler(&cell->get_dof_handler())
5523 , level_dof_access(lda)
5524{}
5525
5526
5527
5528template <int dim, int spacedim>
5531 const FEValuesExtractors::Scalar &scalar) const
5532{
5533 AssertIndexRange(scalar.component, fe_values_views_cache.scalars.size());
5534
5535 return fe_values_views_cache.scalars[scalar.component];
5536}
5537
5538
5539
5540template <int dim, int spacedim>
5543 const FEValuesExtractors::Vector &vector) const
5544{
5546 fe_values_views_cache.vectors.size());
5547
5548 return fe_values_views_cache.vectors[vector.first_vector_component];
5549}
5550
5551
5552
5553template <int dim, int spacedim>
5556 const FEValuesExtractors::SymmetricTensor<2> &tensor) const
5557{
5558 Assert(
5559 tensor.first_tensor_component <
5560 fe_values_views_cache.symmetric_second_order_tensors.size(),
5562 0,
5563 fe_values_views_cache.symmetric_second_order_tensors.size()));
5564
5565 return fe_values_views_cache
5566 .symmetric_second_order_tensors[tensor.first_tensor_component];
5567}
5568
5569
5570
5571template <int dim, int spacedim>
5574 const FEValuesExtractors::Tensor<2> &tensor) const
5575{
5577 fe_values_views_cache.second_order_tensors.size());
5578
5579 return fe_values_views_cache
5580 .second_order_tensors[tensor.first_tensor_component];
5581}
5582
5583
5584
5585template <int dim, int spacedim>
5586inline const double &
5587FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
5588 const unsigned int q_point) const
5589{
5590 AssertIndexRange(i, fe->n_dofs_per_cell());
5591 Assert(this->update_flags & update_values,
5592 ExcAccessToUninitializedField("update_values"));
5593 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5594 Assert(present_cell.is_initialized(), ExcNotReinited());
5595 // if the entire FE is primitive,
5596 // then we can take a short-cut:
5597 if (fe->is_primitive())
5598 return this->finite_element_output.shape_values(i, q_point);
5599 else
5600 {
5601 // otherwise, use the mapping
5602 // between shape function
5603 // numbers and rows. note that
5604 // by the assertions above, we
5605 // know that this particular
5606 // shape function is primitive,
5607 // so we can call
5608 // system_to_component_index
5609 const unsigned int row =
5610 this->finite_element_output
5611 .shape_function_to_row_table[i * fe->n_components() +
5612 fe->system_to_component_index(i).first];
5613 return this->finite_element_output.shape_values(row, q_point);
5614 }
5615}
5616
5617
5618
5619template <int dim, int spacedim>
5620inline double
5622 const unsigned int i,
5623 const unsigned int q_point,
5624 const unsigned int component) const
5625{
5626 AssertIndexRange(i, fe->n_dofs_per_cell());
5627 Assert(this->update_flags & update_values,
5628 ExcAccessToUninitializedField("update_values"));
5629 AssertIndexRange(component, fe->n_components());
5630 Assert(present_cell.is_initialized(), ExcNotReinited());
5631
5632 // check whether the shape function
5633 // is non-zero at all within
5634 // this component:
5635 if (fe->get_nonzero_components(i)[component] == false)
5636 return 0;
5637
5638 // look up the right row in the
5639 // table and take the data from
5640 // there
5641 const unsigned int row =
5642 this->finite_element_output
5643 .shape_function_to_row_table[i * fe->n_components() + component];
5644 return this->finite_element_output.shape_values(row, q_point);
5645}
5646
5647
5648
5649template <int dim, int spacedim>
5650inline const Tensor<1, spacedim> &
5651FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
5652 const unsigned int q_point) const
5653{
5654 AssertIndexRange(i, fe->n_dofs_per_cell());
5655 Assert(this->update_flags & update_gradients,
5656 ExcAccessToUninitializedField("update_gradients"));
5657 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5658 Assert(present_cell.is_initialized(), ExcNotReinited());
5659 // if the entire FE is primitive,
5660 // then we can take a short-cut:
5661 if (fe->is_primitive())
5662 return this->finite_element_output.shape_gradients[i][q_point];
5663 else
5664 {
5665 // otherwise, use the mapping
5666 // between shape function
5667 // numbers and rows. note that
5668 // by the assertions above, we
5669 // know that this particular
5670 // shape function is primitive,
5671 // so we can call
5672 // system_to_component_index
5673 const unsigned int row =
5674 this->finite_element_output
5675 .shape_function_to_row_table[i * fe->n_components() +
5676 fe->system_to_component_index(i).first];
5677 return this->finite_element_output.shape_gradients[row][q_point];
5678 }
5679}
5680
5681
5682
5683template <int dim, int spacedim>
5686 const unsigned int i,
5687 const unsigned int q_point,
5688 const unsigned int component) const
5689{
5690 AssertIndexRange(i, fe->n_dofs_per_cell());
5691 Assert(this->update_flags & update_gradients,
5692 ExcAccessToUninitializedField("update_gradients"));
5693 AssertIndexRange(component, fe->n_components());
5694 Assert(present_cell.is_initialized(), ExcNotReinited());
5695 // check whether the shape function
5696 // is non-zero at all within
5697 // this component:
5698 if (fe->get_nonzero_components(i)[component] == false)
5699 return Tensor<1, spacedim>();
5700
5701 // look up the right row in the
5702 // table and take the data from
5703 // there
5704 const unsigned int row =
5705 this->finite_element_output
5706 .shape_function_to_row_table[i * fe->n_components() + component];
5707 return this->finite_element_output.shape_gradients[row][q_point];
5708}
5709
5710
5711
5712template <int dim, int spacedim>
5713inline const Tensor<2, spacedim> &
5715 const unsigned int q_point) const
5716{
5717 AssertIndexRange(i, fe->n_dofs_per_cell());
5718 Assert(this->update_flags & update_hessians,
5719 ExcAccessToUninitializedField("update_hessians"));
5720 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5721 Assert(present_cell.is_initialized(), ExcNotReinited());
5722 // if the entire FE is primitive,
5723 // then we can take a short-cut:
5724 if (fe->is_primitive())
5725 return this->finite_element_output.shape_hessians[i][q_point];
5726 else
5727 {
5728 // otherwise, use the mapping
5729 // between shape function
5730 // numbers and rows. note that
5731 // by the assertions above, we
5732 // know that this particular
5733 // shape function is primitive,
5734 // so we can call
5735 // system_to_component_index
5736 const unsigned int row =
5737 this->finite_element_output
5738 .shape_function_to_row_table[i * fe->n_components() +
5739 fe->system_to_component_index(i).first];
5740 return this->finite_element_output.shape_hessians[row][q_point];
5741 }
5742}
5743
5744
5745
5746template <int dim, int spacedim>
5749 const unsigned int i,
5750 const unsigned int q_point,
5751 const unsigned int component) const
5752{
5753 AssertIndexRange(i, fe->n_dofs_per_cell());
5754 Assert(this->update_flags & update_hessians,
5755 ExcAccessToUninitializedField("update_hessians"));
5756 AssertIndexRange(component, fe->n_components());
5757 Assert(present_cell.is_initialized(), ExcNotReinited());
5758 // check whether the shape function
5759 // is non-zero at all within
5760 // this component:
5761 if (fe->get_nonzero_components(i)[component] == false)
5762 return Tensor<2, spacedim>();
5763
5764 // look up the right row in the
5765 // table and take the data from
5766 // there
5767 const unsigned int row =
5768 this->finite_element_output
5769 .shape_function_to_row_table[i * fe->n_components() + component];
5770 return this->finite_element_output.shape_hessians[row][q_point];
5771}
5772
5773
5774
5775template <int dim, int spacedim>
5776inline const Tensor<3, spacedim> &
5778 const unsigned int i,
5779 const unsigned int q_point) const
5780{
5781 AssertIndexRange(i, fe->n_dofs_per_cell());
5782 Assert(this->update_flags & update_3rd_derivatives,
5783 ExcAccessToUninitializedField("update_3rd_derivatives"));
5784 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5785 Assert(present_cell.is_initialized(), ExcNotReinited());
5786 // if the entire FE is primitive,
5787 // then we can take a short-cut:
5788 if (fe->is_primitive())
5789 return this->finite_element_output.shape_3rd_derivatives[i][q_point];
5790 else
5791 {
5792 // otherwise, use the mapping
5793 // between shape function
5794 // numbers and rows. note that
5795 // by the assertions above, we
5796 // know that this particular
5797 // shape function is primitive,
5798 // so we can call
5799 // system_to_component_index
5800 const unsigned int row =
5801 this->finite_element_output
5802 .shape_function_to_row_table[i * fe->n_components() +
5803 fe->system_to_component_index(i).first];
5804 return this->finite_element_output.shape_3rd_derivatives[row][q_point];
5805 }
5806}
5807
5808
5809
5810template <int dim, int spacedim>
5813 const unsigned int i,
5814 const unsigned int q_point,
5815 const unsigned int component) const
5816{
5817 AssertIndexRange(i, fe->n_dofs_per_cell());
5818 Assert(this->update_flags & update_3rd_derivatives,
5819 ExcAccessToUninitializedField("update_3rd_derivatives"));
5820 AssertIndexRange(component, fe->n_components());
5821 Assert(present_cell.is_initialized(), ExcNotReinited());
5822 // check whether the shape function
5823 // is non-zero at all within
5824 // this component:
5825 if (fe->get_nonzero_components(i)[component] == false)
5826 return Tensor<3, spacedim>();
5827
5828 // look up the right row in the
5829 // table and take the data from
5830 // there
5831 const unsigned int row =
5832 this->finite_element_output
5833 .shape_function_to_row_table[i * fe->n_components() + component];
5834 return this->finite_element_output.shape_3rd_derivatives[row][q_point];
5835}
5836
5837
5838
5839template <int dim, int spacedim>
5840inline const FiniteElement<dim, spacedim> &
5842{
5843 return *fe;
5844}
5845
5846
5847
5848template <int dim, int spacedim>
5849inline const Mapping<dim, spacedim> &
5851{
5852 return *mapping;
5853}
5854
5855
5856
5857template <int dim, int spacedim>
5858inline UpdateFlags
5860{
5861 return this->update_flags;
5862}
5863
5864
5865
5866template <int dim, int spacedim>
5867inline const std::vector<Point<spacedim>> &
5869{
5870 Assert(this->update_flags & update_quadrature_points,
5871 ExcAccessToUninitializedField("update_quadrature_points"));
5872 Assert(present_cell.is_initialized(), ExcNotReinited());
5873 return this->mapping_output.quadrature_points;
5874}
5875
5876
5877
5878template <int dim, int spacedim>
5879inline const std::vector<double> &
5881{
5882 Assert(this->update_flags & update_JxW_values,
5883 ExcAccessToUninitializedField("update_JxW_values"));
5884 Assert(present_cell.is_initialized(), ExcNotReinited());
5885 return this->mapping_output.JxW_values;
5886}
5887
5888
5889
5890template <int dim, int spacedim>
5891inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5893{
5894 Assert(this->update_flags & update_jacobians,
5895 ExcAccessToUninitializedField("update_jacobians"));
5896 Assert(present_cell.is_initialized(), ExcNotReinited());
5897 return this->mapping_output.jacobians;
5898}
5899
5900
5901
5902template <int dim, int spacedim>
5903inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5905{
5906 Assert(this->update_flags & update_jacobian_grads,
5907 ExcAccessToUninitializedField("update_jacobians_grads"));
5908 Assert(present_cell.is_initialized(), ExcNotReinited());
5909 return this->mapping_output.jacobian_grads;
5910}
5911
5912
5913
5914template <int dim, int spacedim>
5915inline const Tensor<3, spacedim> &
5917 const unsigned int q_point) const
5918{
5919 Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5920 ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5921 Assert(present_cell.is_initialized(), ExcNotReinited());
5922 return this->mapping_output.jacobian_pushed_forward_grads[q_point];
5923}
5924
5925
5926
5927template <int dim, int spacedim>
5928inline const std::vector<Tensor<3, spacedim>> &
5930{
5931 Assert(this->update_flags & update_jacobian_pushed_forward_grads,
5932 ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
5933 Assert(present_cell.is_initialized(), ExcNotReinited());
5934 return this->mapping_output.jacobian_pushed_forward_grads;
5935}
5936
5937
5938
5939template <int dim, int spacedim>
5942 const unsigned int q_point) const
5943{
5944 Assert(this->update_flags & update_jacobian_2nd_derivatives,
5945 ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5946 Assert(present_cell.is_initialized(), ExcNotReinited());
5947 return this->mapping_output.jacobian_2nd_derivatives[q_point];
5948}
5949
5950
5951
5952template <int dim, int spacedim>
5953inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5955{
5956 Assert(this->update_flags & update_jacobian_2nd_derivatives,
5957 ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
5958 Assert(present_cell.is_initialized(), ExcNotReinited());
5959 return this->mapping_output.jacobian_2nd_derivatives;
5960}
5961
5962
5963
5964template <int dim, int spacedim>
5965inline const Tensor<4, spacedim> &
5967 const unsigned int q_point) const
5968{
5970 ExcAccessToUninitializedField(
5971 "update_jacobian_pushed_forward_2nd_derivatives"));
5972 Assert(present_cell.is_initialized(), ExcNotReinited());
5973 return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[q_point];
5974}
5975
5976
5977
5978template <int dim, int spacedim>
5979inline const std::vector<Tensor<4, spacedim>> &
5981{
5983 ExcAccessToUninitializedField(
5984 "update_jacobian_pushed_forward_2nd_derivatives"));
5985 Assert(present_cell.is_initialized(), ExcNotReinited());
5986 return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5987}
5988
5989
5990
5991template <int dim, int spacedim>
5994 const unsigned int q_point) const
5995{
5996 Assert(this->update_flags & update_jacobian_3rd_derivatives,
5997 ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
5998 Assert(present_cell.is_initialized(), ExcNotReinited());
5999 return this->mapping_output.jacobian_3rd_derivatives[q_point];
6000}
6001
6002
6003
6004template <int dim, int spacedim>
6005inline const std::vector<DerivativeForm<4, dim, spacedim>> &
6007{
6008 Assert(this->update_flags & update_jacobian_3rd_derivatives,
6009 ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
6010 Assert(present_cell.is_initialized(), ExcNotReinited());
6011 return this->mapping_output.jacobian_3rd_derivatives;
6012}
6013
6014
6015
6016template <int dim, int spacedim>
6017inline const Tensor<5, spacedim> &
6019 const unsigned int q_point) const
6020{
6022 ExcAccessToUninitializedField(
6023 "update_jacobian_pushed_forward_3rd_derivatives"));
6024 Assert(present_cell.is_initialized(), ExcNotReinited());
6025 return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[q_point];
6026}
6027
6028
6029
6030template <int dim, int spacedim>
6031inline const std::vector<Tensor<5, spacedim>> &
6033{
6035 ExcAccessToUninitializedField(
6036 "update_jacobian_pushed_forward_3rd_derivatives"));
6037 Assert(present_cell.is_initialized(), ExcNotReinited());
6038 return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
6039}
6040
6041
6042
6043template <int dim, int spacedim>
6044inline const std::vector<DerivativeForm<1, spacedim, dim>> &
6046{
6047 Assert(this->update_flags & update_inverse_jacobians,
6048 ExcAccessToUninitializedField("update_inverse_jacobians"));
6049 Assert(present_cell.is_initialized(), ExcNotReinited());
6050 return this->mapping_output.inverse_jacobians;
6051}
6052
6053
6054
6055template <int dim, int spacedim>
6058{
6059 return {0U, dofs_per_cell};
6060}
6061
6062
6063
6064template <int dim, int spacedim>
6067 const unsigned int start_dof_index) const
6068{
6069 Assert(start_dof_index <= dofs_per_cell,
6070 ExcIndexRange(start_dof_index, 0, dofs_per_cell + 1));
6071 return {start_dof_index, dofs_per_cell};
6072}
6073
6074
6075
6076template <int dim, int spacedim>
6079 const unsigned int end_dof_index) const
6080{
6081 Assert(end_dof_index < dofs_per_cell,
6082 ExcIndexRange(end_dof_index, 0, dofs_per_cell));
6083 return {0U, end_dof_index + 1};
6084}
6085
6086
6087
6088template <int dim, int spacedim>
6091{
6092 return {0U, n_quadrature_points};
6093}
6094
6095
6096
6097template <int dim, int spacedim>
6098inline const Point<spacedim> &
6099FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int q_point) const
6100{
6101 Assert(this->update_flags & update_quadrature_points,
6102 ExcAccessToUninitializedField("update_quadrature_points"));
6103 AssertIndexRange(q_point, this->mapping_output.quadrature_points.size());
6104 Assert(present_cell.is_initialized(), ExcNotReinited());
6105
6106 return this->mapping_output.quadrature_points[q_point];
6107}
6108
6109
6110
6111template <int dim, int spacedim>
6112inline double
6113FEValuesBase<dim, spacedim>::JxW(const unsigned int q_point) const
6114{
6115 Assert(this->update_flags & update_JxW_values,
6116 ExcAccessToUninitializedField("update_JxW_values"));
6117 AssertIndexRange(q_point, this->mapping_output.JxW_values.size());
6118 Assert(present_cell.is_initialized(), ExcNotReinited());
6119
6120 return this->mapping_output.JxW_values[q_point];
6121}
6122
6123
6124
6125template <int dim, int spacedim>
6127FEValuesBase<dim, spacedim>::jacobian(const unsigned int q_point) const
6128{
6129 Assert(this->update_flags & update_jacobians,
6130 ExcAccessToUninitializedField("update_jacobians"));
6131 AssertIndexRange(q_point, this->mapping_output.jacobians.size());
6132 Assert(present_cell.is_initialized(), ExcNotReinited());
6133
6134 return this->mapping_output.jacobians[q_point];
6135}
6136
6137
6138
6139template <int dim, int spacedim>
6141FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int q_point) const
6142{
6143 Assert(this->update_flags & update_jacobian_grads,
6144 ExcAccessToUninitializedField("update_jacobians_grads"));
6145 AssertIndexRange(q_point, this->mapping_output.jacobian_grads.size());
6146 Assert(present_cell.is_initialized(), ExcNotReinited());
6147
6148 return this->mapping_output.jacobian_grads[q_point];
6149}
6150
6151
6152
6153template <int dim, int spacedim>
6155FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int q_point) const
6156{
6157 Assert(this->update_flags & update_inverse_jacobians,
6158 ExcAccessToUninitializedField("update_inverse_jacobians"));
6159 AssertIndexRange(q_point, this->mapping_output.inverse_jacobians.size());
6160 Assert(present_cell.is_initialized(), ExcNotReinited());
6161
6162 return this->mapping_output.inverse_jacobians[q_point];
6163}
6164
6165
6166
6167template <int dim, int spacedim>
6168inline const Tensor<1, spacedim> &
6169FEValuesBase<dim, spacedim>::normal_vector(const unsigned int q_point) const
6170{
6171 Assert(this->update_flags & update_normal_vectors,
6173 "update_normal_vectors")));
6174 AssertIndexRange(q_point, this->mapping_output.normal_vectors.size());
6175 Assert(present_cell.is_initialized(), ExcNotReinited());
6176
6177 return this->mapping_output.normal_vectors[q_point];
6178}
6179
6180
6181
6182/*--------------------- Inline functions: FEValues --------------------------*/
6183
6184
6185template <int dim, int spacedim>
6186inline const Quadrature<dim> &
6188{
6189 return quadrature;
6190}
6191
6192
6193
6194template <int dim, int spacedim>
6195inline const FEValues<dim, spacedim> &
6197{
6198 return *this;
6199}
6200
6201
6202/*---------------------- Inline functions: FEFaceValuesBase -----------------*/
6203
6204
6205template <int dim, int spacedim>
6206inline unsigned int
6208{
6209 return present_face_no;
6210}
6211
6212
6213template <int dim, int spacedim>
6214inline unsigned int
6216{
6217 return present_face_index;
6218}
6219
6220
6221/*----------------------- Inline functions: FE*FaceValues -------------------*/
6222
6223template <int dim, int spacedim>
6224inline const Quadrature<dim - 1> &
6226{
6227 return quadrature[quadrature.size() == 1 ? 0 : present_face_no];
6228}
6229
6230
6231
6232template <int dim, int spacedim>
6233inline const FEFaceValues<dim, spacedim> &
6235{
6236 return *this;
6237}
6238
6239
6240
6241template <int dim, int spacedim>
6242inline const FESubfaceValues<dim, spacedim> &
6244{
6245 return *this;
6246}
6247
6248
6249
6250template <int dim, int spacedim>
6251inline const Tensor<1, spacedim> &
6252FEFaceValuesBase<dim, spacedim>::boundary_form(const unsigned int q_point) const
6253{
6254 AssertIndexRange(q_point, this->mapping_output.boundary_forms.size());
6255 Assert(this->update_flags & update_boundary_forms,
6257 "update_boundary_forms")));
6258
6259 return this->mapping_output.boundary_forms[q_point];
6260}
6261
6262#endif // DOXYGEN
6263
6265
6266#endif
const Tensor< 1, spacedim > & boundary_form(const unsigned int q_point) const
const Quadrature< dim - 1 > & get_quadrature() const
unsigned int get_face_index() const
unsigned int present_face_no
Definition fe_values.h:4289
FEFaceValuesBase(const unsigned int dofs_per_cell, const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const hp::QCollection< dim - 1 > &quadrature)
unsigned int present_face_index
Definition fe_values.h:4295
FEFaceValuesBase(const unsigned int dofs_per_cell, const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim - 1 > &quadrature)
std::size_t memory_consumption() const
const std::vector< Tensor< 1, spacedim > > & get_boundary_forms() const
unsigned int get_face_number() const
const hp::QCollection< dim - 1 > quadrature
Definition fe_values.h:4300
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access > > &cell, const typename Triangulation< dim, spacedim >::face_iterator &face)
const FEFaceValues< dim, spacedim > & get_present_fe_values() const
FEFaceValues(const FiniteElement< dim, spacedim > &fe, const Quadrature< dim - 1 > &quadrature, const UpdateFlags update_flags)
void initialize(const UpdateFlags update_flags)
void reinit(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no)
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access > > &cell, const unsigned int face_no)
void do_reinit(const unsigned int face_no)
FEFaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim - 1 > &quadrature, const UpdateFlags update_flags)
FEFaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const hp::QCollection< dim - 1 > &quadrature, const UpdateFlags update_flags)
FEFaceValues(const FiniteElement< dim, spacedim > &fe, const hp::QCollection< dim - 1 > &quadrature, const UpdateFlags update_flags)
void reinit(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const typename Triangulation< dim, spacedim >::face_iterator &face)
FESubfaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim - 1 > &face_quadrature, const UpdateFlags update_flags)
void reinit(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no)
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access > > &cell, const unsigned int face_no, const unsigned int subface_no)
void initialize(const UpdateFlags update_flags)
const FESubfaceValues< dim, spacedim > & get_present_fe_values() const
FESubfaceValues(const FiniteElement< dim, spacedim > &fe, const Quadrature< dim - 1 > &face_quadrature, const UpdateFlags update_flags)
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access > > &cell, const typename Triangulation< dim, spacedim >::face_iterator &face, const typename Triangulation< dim, spacedim >::face_iterator &subface)
void do_reinit(const unsigned int face_no, const unsigned int subface_no)
FESubfaceValues(const FiniteElement< dim, spacedim > &fe, const hp::QCollection< dim - 1 > &face_quadrature, const UpdateFlags update_flags)
FESubfaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const hp::QCollection< dim - 1 > &face_quadrature, const UpdateFlags update_flags)
void reinit(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const typename Triangulation< dim, spacedim >::face_iterator &face, const typename Triangulation< dim, spacedim >::face_iterator &subface)
const DoFHandler< dim, spacedim > * dof_handler
Definition fe_values.h:3902
CellIteratorContainer(const TriaIterator< DoFCellAccessor< dim, spacedim, lda > > &cell)
Triangulation< dim, spacedim >::cell_iterator cell
Definition fe_values.h:3901
CellSimilarity::Similarity cell_similarity
Definition fe_values.h:4019
const Tensor< 5, spacedim > & jacobian_pushed_forward_3rd_derivative(const unsigned int q_point) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_ending_at(const unsigned int end_dof_index) const
const FEValuesViews::Vector< dim, spacedim > & operator[](const FEValuesExtractors::Vector &vector) const
const std::vector< double > & get_JxW_values() const
CellIteratorContainer present_cell
Definition fe_values.h:3911
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int q_point) const
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition fe_values.h:4034
FEValuesBase(const FEValuesBase &)=delete
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
boost::signals2::connection tria_listener_mesh_transform
Definition fe_values.h:3929
const std::vector< Point< spacedim > > & get_quadrature_points() const
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
const DerivativeForm< 2, dim, spacedim > & jacobian_grad(const unsigned int q_point) const
internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > mapping_output
Definition fe_values.h:3971
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition fe_values.h:3956
const std::vector< DerivativeForm< 3, dim, spacedim > > & get_jacobian_2nd_derivatives() const
UpdateFlags get_update_flags() const
const FEValuesViews::Tensor< 2, dim, spacedim > & operator[](const FEValuesExtractors::Tensor< 2 > &tensor) const
const unsigned int dofs_per_cell
Definition fe_values.h:2451
Tensor< 3, spacedim > shape_3rd_derivative_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
const Tensor< 1, spacedim > & normal_vector(const unsigned int q_point) const
UpdateFlags update_flags
Definition fe_values.h:4001
Tensor< 2, spacedim > shape_hessian_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
Definition fe_values.h:3979
const Point< spacedim > & quadrature_point(const unsigned int q_point) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_starting_at(const unsigned int start_dof_index) const
double shape_value_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
const Mapping< dim, spacedim > & get_mapping() const
const unsigned int n_quadrature_points
Definition fe_values.h:2433
const Tensor< 3, spacedim > & jacobian_pushed_forward_grad(const unsigned int q_point) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > quadrature_point_indices() const
const Tensor< 2, spacedim > & shape_hessian(const unsigned int i, const unsigned int q_point) const
const std::vector< Tensor< 5, spacedim > > & get_jacobian_pushed_forward_3rd_derivatives() const
const Tensor< 4, spacedim > & jacobian_pushed_forward_2nd_derivative(const unsigned int q_point) const
const Tensor< 3, spacedim > & shape_3rd_derivative(const unsigned int i, const unsigned int q_point) const
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition fe_values.h:3964
const std::vector< DerivativeForm< 2, dim, spacedim > > & get_jacobian_grads() const
boost::signals2::connection tria_listener_refinement
Definition fe_values.h:3920
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int q_point) const
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int q_point) const
const std::vector< Tensor< 4, spacedim > > & get_jacobian_pushed_forward_2nd_derivatives() const
const DerivativeForm< 3, dim, spacedim > & jacobian_2nd_derivative(const unsigned int q_point) const
Tensor< 1, spacedim > shape_grad_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition fe_values.h:3995
const FEValuesViews::SymmetricTensor< 2, dim, spacedim > & operator[](const FEValuesExtractors::SymmetricTensor< 2 > &tensor) const
const Tensor< 1, spacedim > & shape_grad(const unsigned int i, const unsigned int q_point) const
const std::vector< DerivativeForm< 4, dim, spacedim > > & get_jacobian_3rd_derivatives() const
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > fe_data
Definition fe_values.h:3987
const FiniteElement< dim, spacedim > & get_fe() const
double JxW(const unsigned int q_point) const
const double & shape_value(const unsigned int i, const unsigned int q_point) const
FEValuesBase & operator=(const FEValuesBase &)=delete
const std::vector< DerivativeForm< 1, spacedim, dim > > & get_inverse_jacobians() const
const unsigned int max_n_quadrature_points
Definition fe_values.h:2444
void get_function_third_derivatives_from_local_dof_values(const InputVector &dof_values, std::vector< solution_third_derivative_type< typename InputVector::value_type > > &third_derivatives) const
Scalar & operator=(const Scalar< dim, spacedim > &)=delete
typename ProductType< Number, hessian_type >::type solution_hessian_type
Definition fe_values.h:215
value_type value(const unsigned int shape_function, const unsigned int q_point) const
const unsigned int component
Definition fe_values.h:635
void get_function_hessians_from_local_dof_values(const InputVector &dof_values, std::vector< solution_hessian_type< typename InputVector::value_type > > &hessians) const
void get_function_values_from_local_dof_values(const InputVector &dof_values, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
std::vector< ShapeFunctionData > shape_function_data
Definition fe_values.h:640
void get_function_laplacians(const InputVector &fe_function, std::vector< solution_laplacian_type< typename InputVector::value_type > > &laplacians) const
typename ProductType< Number, value_type >::type solution_laplacian_type
Definition fe_values.h:205
void get_function_third_derivatives(const InputVector &fe_function, std::vector< solution_third_derivative_type< typename InputVector::value_type > > &third_derivatives) const
Scalar & operator=(Scalar< dim, spacedim > &&) noexcept=default
third_derivative_type third_derivative(const unsigned int shape_function, const unsigned int q_point) const
void get_function_hessians(const InputVector &fe_function, std::vector< solution_hessian_type< typename InputVector::value_type > > &hessians) const
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
Definition fe_values.h:225
typename ProductType< Number, value_type >::type solution_value_type
Definition fe_values.h:185
typename ProductType< Number, gradient_type >::type solution_gradient_type
Definition fe_values.h:195
void get_function_laplacians_from_local_dof_values(const InputVector &dof_values, std::vector< solution_laplacian_type< typename InputVector::value_type > > &laplacians) const
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition fe_values.h:629
void get_function_gradients(const InputVector &fe_function, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
void get_function_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
Scalar(Scalar< dim, spacedim > &&)=default
Scalar(const Scalar< dim, spacedim > &)=delete
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
void get_function_values(const InputVector &fe_function, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
SymmetricTensor(const SymmetricTensor< 2, dim, spacedim > &)=delete
SymmetricTensor(SymmetricTensor< 2, dim, spacedim > &&)=default
SymmetricTensor & operator=(SymmetricTensor< 2, dim, spacedim > &&) noexcept=default
typename ProductType< Number, value_type >::type solution_value_type
Definition fe_values.h:1506
SymmetricTensor & operator=(const SymmetricTensor< 2, dim, spacedim > &)=delete
typename ProductType< Number, divergence_type >::type solution_divergence_type
Definition fe_values.h:1516
typename ProductType< Number, value_type >::type solution_value_type
Definition fe_values.h:1842
typename ProductType< Number, divergence_type >::type solution_divergence_type
Definition fe_values.h:1852
typename ProductType< Number, gradient_type >::type solution_gradient_type
Definition fe_values.h:1862
divergence_type divergence(const unsigned int shape_function, const unsigned int q_point) const
Tensor(const Tensor< 2, dim, spacedim > &)=delete
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition fe_values.h:2190
Tensor & operator=(const Tensor< 2, dim, spacedim > &)=delete
Tensor & operator=(Tensor< 2, dim, spacedim > &&)=default
value_type value(const unsigned int shape_function, const unsigned int q_point) const
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
Tensor(Tensor< 2, dim, spacedim > &&)=default
std::vector< ShapeFunctionData > shape_function_data
Definition fe_values.h:2201
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
Definition fe_values.h:812
typename ProductType< Number, divergence_type >::type solution_divergence_type
Definition fe_values.h:773
Vector(const Vector< dim, spacedim > &)=delete
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
typename ProductType< Number, hessian_type >::type solution_hessian_type
Definition fe_values.h:802
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition fe_values.h:1436
Vector & operator=(Vector< dim, spacedim > &&)=default
typename ProductType< Number, symmetric_gradient_type >::type solution_symmetric_gradient_type
Definition fe_values.h:763
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
divergence_type divergence(const unsigned int shape_function, const unsigned int q_point) const
Vector & operator=(const Vector< dim, spacedim > &)=delete
third_derivative_type third_derivative(const unsigned int shape_function, const unsigned int q_point) const
typename ProductType< Number, gradient_type >::type solution_gradient_type
Definition fe_values.h:753
typename ProductType< Number, value_type >::type solution_value_type
Definition fe_values.h:743
symmetric_gradient_type symmetric_gradient(const unsigned int shape_function, const unsigned int q_point) const
Vector(Vector< dim, spacedim > &&)=default
typename ::internal::CurlType< spacedim >::type curl_type
Definition fe_values.h:720
const unsigned int first_vector_component
Definition fe_values.h:1442
typename ProductType< Number, curl_type >::type solution_curl_type
Definition fe_values.h:792
std::vector< ShapeFunctionData > shape_function_data
Definition fe_values.h:1447
value_type value(const unsigned int shape_function, const unsigned int q_point) const
curl_type curl(const unsigned int shape_function, const unsigned int q_point) const
typename ProductType< Number, value_type >::type solution_laplacian_type
Definition fe_values.h:783
FEValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, const UpdateFlags update_flags)
void do_reinit()
const Quadrature< dim > & get_quadrature() const
const Quadrature< dim > quadrature
Definition fe_values.h:4170
FEValues(const FiniteElement< dim, spacedim > &fe, const hp::QCollection< dim > &quadrature, const UpdateFlags update_flags)
const FEValues< dim, spacedim > & get_present_fe_values() const
FEValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const hp::QCollection< dim > &quadrature, const UpdateFlags update_flags)
FEValues(const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, const UpdateFlags update_flags)
void initialize(const UpdateFlags update_flags)
void reinit(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access > > &cell)
std::size_t memory_consumption() const
Abstract base class for mapping classes.
Definition mapping.h:317
Definition point.h:112
typename Tensor< rank_ - 1, dim, Number >::tensor_type value_type
Definition tensor.h:549
static constexpr TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
Number value_type
Definition vector.h:128
#define DEAL_II_DEPRECATED
Definition config.h:172
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define DeclException0(Exception0)
Definition exceptions.h:465
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:488
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:510
static ::ExceptionBase & ExcMessage(std::string arg1)
TriaIterator< CellAccessor< dim, spacedim > > cell_iterator
Definition tria.h:1370
UpdateFlags
@ update_jacobian_pushed_forward_2nd_derivatives
@ update_jacobian_pushed_forward_grads
@ update_hessians
Second derivatives of shape functions.
@ update_jacobian_3rd_derivatives
@ update_values
Shape function values.
@ update_jacobian_grads
Gradient of volume element.
@ update_normal_vectors
Normal vectors.
@ update_3rd_derivatives
Third derivatives of shape functions.
@ update_JxW_values
Transformed quadrature weights.
@ update_jacobians
Volume element.
@ update_inverse_jacobians
Volume element.
@ update_gradients
Shape function gradients.
@ update_quadrature_points
Transformed quadrature points.
@ update_jacobian_pushed_forward_3rd_derivatives
@ update_boundary_forms
Outer normal vector, not normalized.
@ update_jacobian_2nd_derivatives
typename ::internal::FEValuesViews::ViewType< dim, spacedim, Extractor >::type View
Definition fe_values.h:2309
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
boost::integer_range< IncrementableType > iota_view
Definition iota_view.h:46
STL namespace.
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type laplacian_type
Definition fe_values.h:258
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type value_type
Definition fe_values.h:242
typename ProductType< Number, typename Scalar< dim, spacedim >::gradient_type >::type gradient_type
Definition fe_values.h:250
typename ProductType< Number, typename Scalar< dim, spacedim >::hessian_type >::type hessian_type
Definition fe_values.h:266
typename ProductType< Number, typename Scalar< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition fe_values.h:274
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition fe_values.h:1542
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::value_type >::type value_type
Definition fe_values.h:1534
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::gradient_type >::type gradient_type
Definition fe_values.h:1896
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::value_type >::type value_type
Definition fe_values.h:1880
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition fe_values.h:1888
typename ProductType< Number, typename Vector< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition fe_values.h:885
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type laplacian_type
Definition fe_values.h:861
typename ProductType< Number, typename Vector< dim, spacedim >::hessian_type >::type hessian_type
Definition fe_values.h:877
typename ProductType< Number, typename Vector< dim, spacedim >::symmetric_gradient_type >::type symmetric_gradient_type
Definition fe_values.h:845
typename ProductType< Number, typename Vector< dim, spacedim >::gradient_type >::type gradient_type
Definition fe_values.h:837
typename ProductType< Number, typename Vector< dim, spacedim >::curl_type >::type curl_type
Definition fe_values.h:869
typename ProductType< Number, typename Vector< dim, spacedim >::divergence_type >::type divergence_type
Definition fe_values.h:853
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type value_type
Definition fe_values.h:829
typename internal::ProductTypeImpl< typename std::decay< T >::type, typename std::decay< U >::type >::type type
std::vector<::FEValuesViews::Scalar< dim, spacedim > > scalars
Definition fe_values.h:2286
std::vector<::FEValuesViews::Vector< dim, spacedim > > vectors
Definition fe_values.h:2287
std::vector<::FEValuesViews::SymmetricTensor< 2, dim, spacedim > > symmetric_second_order_tensors
Definition fe_values.h:2289
std::vector<::FEValuesViews::Tensor< 2, dim, spacedim > > second_order_tensors
Definition fe_values.h:2291
typename ::FEValuesViews::SymmetricTensor< rank, dim, spacedim > type
Definition fe_values.h:2269
DEAL_II_HOST constexpr SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)