Loading [MathJax]/extensions/TeX/AMSsymbols.js
 deal.II version GIT relicensing-3075-gc235bd4825 2025-04-15 08:10:00+00:00
\(\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\}}\)
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
fe_interface_values.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2019 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_fe_interface_values_h
16#define dealii_fe_interface_values_h
17
18#include <deal.II/base/config.h>
19
21
23#include <deal.II/fe/mapping.h>
24
28
30
31#ifndef DOXYGEN
32template <int dim, int spacedim>
34#endif
35
41{
45 template <int dim, int spacedim = dim>
46 class Base
47 {
48 public:
53
54 protected:
59
64 template <class InputVector, class OutputVector>
65 void
66 get_local_dof_values(const InputVector &dof_values,
67 OutputVector &local_dof_values) const;
68 };
69
70
71
75 template <int dim, int spacedim = dim>
76 class Scalar : public Base<dim, spacedim>
77 {
78 public:
82 using value_type = double;
83
90
97
104
111 template <typename Number>
113
120 template <typename Number>
123
130 template <typename Number>
133
140 template <typename Number>
143
148 const unsigned int component);
149
175 value(const bool here_or_there,
176 const unsigned int interface_dof_index,
177 const unsigned int q_point) const;
178
197 jump_in_values(const unsigned int interface_dof_index,
198 const unsigned int q_point) const;
199
211 jump_in_gradients(const unsigned int interface_dof_index,
212 const unsigned int q_point) const;
213
225 jump_in_hessians(const unsigned int interface_dof_index,
226 const unsigned int q_point) const;
227
239 jump_in_third_derivatives(const unsigned int interface_dof_index,
240 const unsigned int q_point) const;
241
260 average_of_values(const unsigned int interface_dof_index,
261 const unsigned int q_point) const;
262
274 average_of_gradients(const unsigned int interface_dof_index,
275 const unsigned int q_point) const;
276
289 average_of_hessians(const unsigned int interface_dof_index,
290 const unsigned int q_point) const;
291
318 template <class InputVector>
319 void
321 const bool here_or_there,
322 const InputVector &fe_function,
324 &values) const;
325
348 template <class InputVector>
349 void
351 const bool here_or_there,
352 const InputVector &local_dof_values,
354 &values) const;
355
376 template <class InputVector>
377 void
379 const InputVector &fe_function,
381 &values) const;
382
390 template <class InputVector>
391 void
393 const InputVector &local_dof_values,
395 &values) const;
396
410 template <class InputVector>
411 void
413 const InputVector &fe_function,
415 &gradients) const;
416
424 template <class InputVector>
425 void
427 const InputVector &local_dof_values,
429 &gradients) const;
430
444 template <class InputVector>
445 void
447 const InputVector &fe_function,
449 &hessians) const;
450
458 template <class InputVector>
459 void
461 const InputVector &local_dof_values,
463 &hessians) const;
464
479 template <class InputVector>
480 void
482 const InputVector &fe_function,
483 std::vector<
485 &third_derivatives) const;
486
494 template <class InputVector>
495 void
497 const InputVector &local_dof_values,
498 std::vector<
500 &third_derivatives) const;
501
522 template <class InputVector>
523 void
525 const InputVector &fe_function,
527 &values) const;
528
536 template <class InputVector>
537 void
539 const InputVector &local_dof_values,
541 &values) const;
542
556 template <class InputVector>
557 void
559 const InputVector &fe_function,
561 &gradients) const;
562
570 template <class InputVector>
571 void
573 const InputVector &local_dof_values,
575 &gradients) const;
576
590 template <class InputVector>
591 void
593 const InputVector &fe_function,
595 &hessians) const;
596
604 template <class InputVector>
605 void
607 const InputVector &local_dof_values,
609 &hessians) const;
610
613 private:
618 };
619
620
621
625 template <int dim, int spacedim = dim>
626 class Vector : public Base<dim, spacedim>
627 {
628 public:
634
641
649
657
664 template <typename Number>
666
673 template <typename Number>
676
683 template <typename Number>
686
693 template <typename Number>
696
701 const unsigned int first_vector_component);
702
728 value(const bool here_or_there,
729 const unsigned int interface_dof_index,
730 const unsigned int q_point) const;
731
749 jump_in_values(const unsigned int interface_dof_index,
750 const unsigned int q_point) const;
751
762 jump_in_gradients(const unsigned int interface_dof_index,
763 const unsigned int q_point) const;
764
771 jump_gradient(const unsigned int interface_dof_index,
772 const unsigned int q_point) const;
773
785 jump_in_hessians(const unsigned int interface_dof_index,
786 const unsigned int q_point) const;
787
794 jump_hessian(const unsigned int interface_dof_index,
795 const unsigned int q_point) const;
796
808 jump_in_third_derivatives(const unsigned int interface_dof_index,
809 const unsigned int q_point) const;
810
828 average_of_values(const unsigned int interface_dof_index,
829 const unsigned int q_point) const;
830
841 average_of_gradients(const unsigned int interface_dof_index,
842 const unsigned int q_point) const;
843
856 average_of_hessians(const unsigned int interface_dof_index,
857 const unsigned int q_point) const;
858
865 average_hessian(const unsigned int interface_dof_index,
866 const unsigned int q_point) const;
867
894 template <class InputVector>
895 void
897 const bool here_or_there,
898 const InputVector &fe_function,
900 &values) const;
901
924 template <class InputVector>
925 void
927 const bool here_or_there,
928 const InputVector &local_dof_values,
930 &values) const;
931
952 template <class InputVector>
953 void
955 const InputVector &fe_function,
957 &values) const;
958
966 template <class InputVector>
967 void
969 const InputVector &local_dof_values,
971 &values) const;
972
986 template <class InputVector>
987 void
989 const InputVector &fe_function,
991 &gradients) const;
992
1000 template <class InputVector>
1001 void
1003 const InputVector &local_dof_values,
1005 &gradients) const;
1006
1020 template <class InputVector>
1021 void
1023 const InputVector &fe_function,
1025 &hessians) const;
1026
1034 template <class InputVector>
1035 void
1037 const InputVector &local_dof_values,
1039 &hessians) const;
1040
1055 template <class InputVector>
1056 void
1058 const InputVector &fe_function,
1059 std::vector<
1061 &third_derivatives) const;
1062
1070 template <class InputVector>
1071 void
1073 const InputVector &local_dof_values,
1074 std::vector<
1076 &third_derivatives) const;
1077
1098 template <class InputVector>
1099 void
1101 const InputVector &fe_function,
1103 &values) const;
1104
1112 template <class InputVector>
1113 void
1115 const InputVector &local_dof_values,
1117 &values) const;
1118
1132 template <class InputVector>
1133 void
1135 const InputVector &fe_function,
1137 &gradients) const;
1138
1146 template <class InputVector>
1147 void
1149 const InputVector &local_dof_values,
1151 &gradients) const;
1152
1166 template <class InputVector>
1167 void
1169 const InputVector &fe_function,
1171 &hessians) const;
1172
1180 template <class InputVector>
1181 void
1183 const InputVector &local_dof_values,
1185 &hessians) const;
1186
1189 private:
1194 };
1195} // namespace FEInterfaceViews
1196
1197
1198namespace internal
1199{
1201 {
1206 template <int dim, int spacedim, typename Extractor>
1208 {};
1209
1217 template <int dim, int spacedim>
1218 struct ViewType<dim, spacedim, FEValuesExtractors::Scalar>
1219 {
1220 using type = typename ::FEInterfaceViews::Scalar<dim, spacedim>;
1221 };
1222
1230 template <int dim, int spacedim>
1231 struct ViewType<dim, spacedim, FEValuesExtractors::Vector>
1232 {
1233 using type = typename ::FEInterfaceViews::Vector<dim, spacedim>;
1234 };
1235 } // namespace FEInterfaceViews
1236} // namespace internal
1237
1238namespace FEInterfaceViews
1239{
1244 template <int dim, int spacedim, typename Extractor>
1245 using View = typename ::internal::FEInterfaceViews::
1246 ViewType<dim, spacedim, Extractor>::type;
1247} // namespace FEInterfaceViews
1248
1249
1250
1275template <int dim, int spacedim = dim>
1277{
1278public:
1282 const unsigned int n_quadrature_points;
1283
1291 const Quadrature<dim - 1> &quadrature,
1292 const UpdateFlags update_flags);
1293
1301 const hp::QCollection<dim - 1> &quadrature,
1302 const UpdateFlags update_flags);
1303
1311 const Quadrature<dim - 1> &quadrature,
1312 const UpdateFlags update_flags);
1313
1319 const hp::MappingCollection<dim, spacedim> &mapping_collection,
1320 const hp::FECollection<dim, spacedim> &fe_collection,
1321 const hp::QCollection<dim - 1> &quadrature_collection,
1322 const UpdateFlags update_flags);
1323
1328 const hp::QCollection<dim - 1> &quadrature_collection,
1329 const UpdateFlags update_flags);
1330
1447 template <typename CellIteratorType, typename CellNeighborIteratorType>
1448 void
1449 reinit(const CellIteratorType &cell,
1450 const unsigned int face_no,
1451 const unsigned int sub_face_no,
1452 const CellNeighborIteratorType &cell_neighbor,
1453 const unsigned int face_no_neighbor,
1454 const unsigned int sub_face_no_neighbor,
1455 const unsigned int q_index = numbers::invalid_unsigned_int,
1456 const unsigned int mapping_index = numbers::invalid_unsigned_int,
1457 const unsigned int fe_index = numbers::invalid_unsigned_int,
1458 const unsigned int fe_index_neighbor = numbers::invalid_unsigned_int);
1459
1486 template <typename CellIteratorType>
1487 void
1488 reinit(const CellIteratorType &cell,
1489 const unsigned int face_no,
1490 const unsigned int q_index = numbers::invalid_unsigned_int,
1491 const unsigned int mapping_index = numbers::invalid_unsigned_int,
1492 const unsigned int fe_index = numbers::invalid_unsigned_int);
1493
1502 get_fe_face_values(const unsigned int cell_index) const;
1503
1509
1514 get_fe() const;
1515
1519 const Quadrature<dim - 1> &
1521
1527
1533
1537 const hp::QCollection<dim - 1> &
1539
1544 bool
1546
1552
1560 get_cell(const unsigned int cell_index) const;
1561
1569 unsigned int
1570 get_face_number(const unsigned int cell_index) const;
1571
1583 bool
1585
1597 double
1598 JxW(const unsigned int quadrature_point) const;
1599
1605 const std::vector<double> &
1607
1616 DEAL_II_DEPRECATED_EARLY_WITH_COMMENT("Use the function normal_vector().")
1617 Tensor<1, spacedim>
1618 normal(const unsigned int q_point_index) const;
1619
1628 Tensor<1, spacedim>
1629 normal_vector(const unsigned int q_point_index) const;
1630
1640 const std::vector<Tensor<1, spacedim>> &
1642
1650 std_cxx20::ranges::iota_view<unsigned int, unsigned int>
1652
1659 const Point<spacedim> &
1660 quadrature_point(const unsigned int q_point) const;
1661
1667 const std::vector<Point<spacedim>> &
1669
1680 unsigned
1682
1710 std_cxx20::ranges::iota_view<unsigned int, unsigned int>
1712
1723 std::vector<types::global_dof_index>
1725
1738 std::array<unsigned int, 2>
1739 interface_dof_to_dof_indices(const unsigned int interface_dof_index) const;
1740
1771 double
1772 shape_value(const bool here_or_there,
1773 const unsigned int interface_dof_index,
1774 const unsigned int q_point,
1775 const unsigned int component = 0) const;
1776
1804 double
1805 jump_in_shape_values(const unsigned int interface_dof_index,
1806 const unsigned int q_point,
1807 const unsigned int component = 0) const;
1808
1822 Tensor<1, spacedim>
1823 jump_in_shape_gradients(const unsigned int interface_dof_index,
1824 const unsigned int q_point,
1825 const unsigned int component = 0) const;
1826
1841 Tensor<2, spacedim>
1842 jump_in_shape_hessians(const unsigned int interface_dof_index,
1843 const unsigned int q_point,
1844 const unsigned int component = 0) const;
1845
1859 Tensor<3, spacedim>
1860 jump_in_shape_3rd_derivatives(const unsigned int interface_dof_index,
1861 const unsigned int q_point,
1862 const unsigned int component = 0) const;
1863
1886 double
1887 average_of_shape_values(const unsigned int interface_dof_index,
1888 const unsigned int q_point,
1889 const unsigned int component = 0) const;
1890
1904 Tensor<1, spacedim>
1905 average_of_shape_gradients(const unsigned int interface_dof_index,
1906 const unsigned int q_point,
1907 const unsigned int component = 0) const;
1908
1923 Tensor<2, spacedim>
1924 average_of_shape_hessians(const unsigned int interface_dof_index,
1925 const unsigned int q_point,
1926 const unsigned int component = 0) const;
1927
1947 template <class InputVector>
1948 void
1950 const InputVector &fe_function,
1951 std::vector<typename InputVector::value_type> &values) const;
1952
1961 template <class InputVector>
1962 void
1964 const InputVector &fe_function,
1965 std::vector<Tensor<1, spacedim, typename InputVector::value_type>>
1966 &gradients) const;
1967
1975 template <class InputVector>
1976 void
1978 const InputVector &fe_function,
1979 std::vector<Tensor<2, spacedim, typename InputVector::value_type>>
1980 &hessians) const;
1981
1990 template <class InputVector>
1991 void
1993 const InputVector &fe_function,
1994 std::vector<Tensor<3, spacedim, typename InputVector::value_type>>
1995 &third_derivatives) const;
1996
2012 template <class InputVector>
2013 void
2015 const InputVector &fe_function,
2016 std::vector<typename InputVector::value_type> &values) const;
2017
2025 template <class InputVector>
2026 void
2028 const InputVector &fe_function,
2029 std::vector<Tensor<1, spacedim, typename InputVector::value_type>>
2030 &gradients) const;
2031
2039 template <class InputVector>
2040 void
2042 const InputVector &fe_function,
2043 std::vector<Tensor<2, spacedim, typename InputVector::value_type>>
2044 &hessians) const;
2045
2063 FEInterfaceViews::Scalar<dim, spacedim>
2064 operator[](const FEValuesExtractors::Scalar &scalar) const;
2065
2072 FEInterfaceViews::Vector<dim, spacedim>
2073 operator[](const FEValuesExtractors::Vector &vector) const;
2074
2079private:
2083 std::vector<types::global_dof_index> interface_dof_indices;
2084
2090 std::vector<std::array<unsigned int, 2>> dofmap;
2091
2097
2104
// non-hp data
2109
2113 std::unique_ptr<FEFaceValues<dim, spacedim>> internal_fe_face_values;
2114
2119
2124
2128 std::unique_ptr<FESubfaceValues<dim, spacedim>>
2130
// non-hp data
2132
// hp data
2137
2142 std::unique_ptr<hp::FEFaceValues<dim, spacedim>> internal_hp_fe_face_values;
2143
2148 std::unique_ptr<hp::FESubfaceValues<dim, spacedim>>
2150
2155 std::unique_ptr<hp::FEFaceValues<dim, spacedim>>
2157
2162 std::unique_ptr<hp::FESubfaceValues<dim, spacedim>>
2164
2172 "The current function doesn't make sense when used with a "
2173 "FEInterfaceValues object with hp-capabilities.");
2174
2182 "The current function doesn't make sense when used with a "
2183 "FEInterfaceValues object without hp-capabilities.");
2184
// hp data
2186
2187 /*
2188 * Make the view classes friends of this class, since they
2189 * access internal data.
2190 */
2191 template <int, int>
2193 template <int, int>
2195};
2196
2197
2198
2199#ifndef DOXYGEN
2200
2201/*---------------------- Inline functions ---------------------*/
2202
2203template <int dim, int spacedim>
2205 const Mapping<dim, spacedim> &mapping,
2207 const Quadrature<dim - 1> &quadrature,
2208 const UpdateFlags update_flags)
2209 : n_quadrature_points(quadrature.size())
2210 , fe_face_values(nullptr)
2211 , fe_face_values_neighbor(nullptr)
2213 std::make_unique<FEFaceValues<dim, spacedim>>(mapping,
2214 fe,
2215 quadrature,
2216 update_flags))
2218 std::make_unique<FESubfaceValues<dim, spacedim>>(mapping,
2219 fe,
2220 quadrature,
2221 update_flags))
2223 std::make_unique<FEFaceValues<dim, spacedim>>(mapping,
2224 fe,
2225 quadrature,
2226 update_flags))
2228 std::make_unique<FESubfaceValues<dim, spacedim>>(mapping,
2229 fe,
2230 quadrature,
2231 update_flags))
2232{}
2233
2234
2235
2236template <int dim, int spacedim>
2239 const Quadrature<dim - 1> &quadrature,
2240 const UpdateFlags update_flags)
2242 fe.reference_cell().template get_default_linear_mapping<dim, spacedim>(),
2243 fe,
2244 quadrature,
2245 update_flags)
2246{}
2247
2248
2249
2250template <int dim, int spacedim>
2252 const Mapping<dim, spacedim> &mapping,
2254 const hp::QCollection<dim - 1> &quadrature,
2255 const UpdateFlags update_flags)
2256 : n_quadrature_points(quadrature.max_n_quadrature_points())
2257 , fe_face_values(nullptr)
2258 , fe_face_values_neighbor(nullptr)
2259 , internal_fe_face_values(
2260 std::make_unique<FEFaceValues<dim, spacedim>>(mapping,
2261 fe,
2262 quadrature,
2263 update_flags))
2264 , internal_fe_subface_values(
2265 std::make_unique<FESubfaceValues<dim, spacedim>>(mapping,
2266 fe,
2267 quadrature,
2268 update_flags))
2269 , internal_fe_face_values_neighbor(
2270 std::make_unique<FEFaceValues<dim, spacedim>>(mapping,
2271 fe,
2272 quadrature[0],
2273 update_flags))
2274 , internal_fe_subface_values_neighbor(
2275 std::make_unique<FESubfaceValues<dim, spacedim>>(mapping,
2276 fe,
2277 quadrature[0],
2278 update_flags))
2279{}
2280
2281
2282
2283template <int dim, int spacedim>
2285 const hp::MappingCollection<dim, spacedim> &mapping_collection,
2286 const hp::FECollection<dim, spacedim> &fe_collection,
2287 const hp::QCollection<dim - 1> &quadrature_collection,
2288 const UpdateFlags update_flags)
2289 : n_quadrature_points(quadrature_collection.max_n_quadrature_points())
2290 , fe_face_values(nullptr)
2291 , fe_face_values_neighbor(nullptr)
2292 , internal_hp_fe_face_values(
2293 std::make_unique<hp::FEFaceValues<dim, spacedim>>(mapping_collection,
2294 fe_collection,
2295 quadrature_collection,
2296 update_flags))
2297 , internal_hp_fe_subface_values(
2298 std::make_unique<hp::FESubfaceValues<dim, spacedim>>(
2299 mapping_collection,
2300 fe_collection,
2301 quadrature_collection,
2302 update_flags))
2303 , internal_hp_fe_face_values_neighbor(
2304 std::make_unique<hp::FEFaceValues<dim, spacedim>>(mapping_collection,
2305 fe_collection,
2306 quadrature_collection,
2307 update_flags))
2308 , internal_hp_fe_subface_values_neighbor(
2309 std::make_unique<hp::FESubfaceValues<dim, spacedim>>(
2310 mapping_collection,
2311 fe_collection,
2312 quadrature_collection,
2313 update_flags))
2314{
2315 AssertDimension(dim, spacedim);
2316}
2317
2318
2319
2320template <int dim, int spacedim>
2322 const hp::FECollection<dim, spacedim> &fe_collection,
2323 const hp::QCollection<dim - 1> &quadrature_collection,
2324 const UpdateFlags update_flags)
2325 : FEInterfaceValues(fe_collection.get_reference_cell_default_linear_mapping(),
2326 fe_collection,
2327 quadrature_collection,
2328 update_flags)
2329{}
2330
2331
2332
2333template <int dim, int spacedim>
2334template <typename CellIteratorType, typename CellNeighborIteratorType>
2335void
2337 const CellIteratorType &cell,
2338 const unsigned int face_no,
2339 const unsigned int sub_face_no,
2340 const CellNeighborIteratorType &cell_neighbor,
2341 const unsigned int face_no_neighbor,
2342 const unsigned int sub_face_no_neighbor,
2343 const unsigned int q_index,
2344 const unsigned int mapping_index,
2345 const unsigned int fe_index_in,
2346 const unsigned int fe_index_neighbor_in)
2347{
2348 Assert(internal_fe_face_values || internal_hp_fe_face_values,
2350
2351 constexpr bool is_dof_cell_accessor =
2352 std::is_same_v<DoFCellAccessor<dim, spacedim, true>,
2353 typename CellIteratorType::AccessorType> ||
2354 std::is_same_v<DoFCellAccessor<dim, spacedim, false>,
2355 typename CellIteratorType::AccessorType>;
2356
2357 constexpr bool is_dof_cell_accessor_neighbor =
2358 std::is_same_v<DoFCellAccessor<dim, spacedim, true>,
2359 typename CellNeighborIteratorType::AccessorType> ||
2360 std::is_same_v<DoFCellAccessor<dim, spacedim, false>,
2361 typename CellNeighborIteratorType::AccessorType>;
2362
2363 unsigned int active_fe_index = 0;
2364 unsigned int active_fe_index_neighbor = 0;
2365
2366 if (internal_fe_face_values)
2367 {
2368 if (sub_face_no == numbers::invalid_unsigned_int)
2369 {
2370 internal_fe_face_values->reinit(cell, face_no);
2371 fe_face_values = internal_fe_face_values.get();
2372 }
2373 else
2374 {
2375 internal_fe_subface_values->reinit(cell, face_no, sub_face_no);
2376 fe_face_values = internal_fe_subface_values.get();
2377 }
2378 if (sub_face_no_neighbor == numbers::invalid_unsigned_int)
2379 {
2380 internal_fe_face_values_neighbor->reinit(cell_neighbor,
2381 face_no_neighbor);
2382 fe_face_values_neighbor = internal_fe_face_values_neighbor.get();
2383 }
2384 else
2385 {
2386 internal_fe_subface_values_neighbor->reinit(cell_neighbor,
2387 face_no_neighbor,
2388 sub_face_no_neighbor);
2389 fe_face_values_neighbor = internal_fe_subface_values_neighbor.get();
2390 }
2391
2392 AssertDimension(fe_face_values->n_quadrature_points,
2393 fe_face_values_neighbor->n_quadrature_points);
2394
2395 const_cast<unsigned int &>(this->n_quadrature_points) =
2396 fe_face_values->n_quadrature_points;
2397 }
2398 else if (internal_hp_fe_face_values)
2399 {
2400 active_fe_index = fe_index_in;
2401 active_fe_index_neighbor =
2402 (fe_index_neighbor_in != numbers::invalid_unsigned_int) ?
2403 fe_index_neighbor_in :
2404 fe_index_in;
2405
2406 if (active_fe_index == numbers::invalid_unsigned_int)
2407 {
2408 if constexpr (is_dof_cell_accessor)
2409 active_fe_index = cell->active_fe_index();
2410 else
2411 active_fe_index = 0;
2412 }
2413
2414 if (active_fe_index_neighbor == numbers::invalid_unsigned_int)
2415 {
2416 if constexpr (is_dof_cell_accessor_neighbor)
2417 active_fe_index_neighbor = cell_neighbor->active_fe_index();
2418 else
2419 active_fe_index_neighbor = 0;
2420 }
2421
2422 unsigned int used_q_index = q_index;
2423 unsigned int used_mapping_index = mapping_index;
2424
2425 // First check. If there is only one element in a collection, and if none
2426 // had been specified explicitly, then that's clearly the one to take:
2427 if (used_q_index == numbers::invalid_unsigned_int)
2428 if (internal_hp_fe_face_values->get_quadrature_collection().size() == 1)
2429 used_q_index = 0;
2430
2431 if (used_mapping_index == numbers::invalid_unsigned_int)
2432 if (internal_hp_fe_face_values->get_mapping_collection().size() == 1)
2433 used_mapping_index = 0;
2434
2435 // Second check: See if the two quadrature objects are the same, because
2436 // in that case it does not matter which one we use. Unfortunately, we
2437 // currently have no way of testing that two mapping objects are the
2438 // same :-(
2439 if (used_q_index == numbers::invalid_unsigned_int)
2440 if (internal_hp_fe_face_values
2441 ->get_quadrature_collection()[active_fe_index] ==
2442 internal_hp_fe_face_values
2443 ->get_quadrature_collection()[active_fe_index_neighbor])
2444 used_q_index = active_fe_index;
2445
2446 // Third check, if the above did not already suffice. We see if we
2447 // can get somewhere via the dominated's finite element index.
2448 if ((used_q_index == numbers::invalid_unsigned_int) ||
2449 (used_mapping_index == numbers::invalid_unsigned_int))
2450 {
2451 const unsigned int dominated_fe_index =
2452 ((used_q_index == numbers::invalid_unsigned_int) ||
2453 (used_mapping_index == numbers::invalid_unsigned_int) ?
2454 internal_hp_fe_face_values->get_fe_collection()
2455 .find_dominated_fe(
2456 {active_fe_index, active_fe_index_neighbor}) :
2458
2459 if (used_q_index == numbers::invalid_unsigned_int)
2460 {
2461 Assert(dominated_fe_index != numbers::invalid_fe_index,
2462 ExcMessage(
2463 "You called this function with 'q_index' left at its "
2464 "default value, but this can only work if one of "
2465 "the two finite elements adjacent to this face "
2466 "dominates the other. See the documentation "
2467 "of this function for more information of how "
2468 "to deal with this situation."));
2469 used_q_index = dominated_fe_index;
2470 }
2471
2472 if (used_mapping_index == numbers::invalid_unsigned_int)
2473 {
2474 Assert(dominated_fe_index != numbers::invalid_fe_index,
2475 ExcMessage(
2476 "You called this function with 'mapping_index' left "
2477 "at its default value, but this can only work if one "
2478 "of the two finite elements adjacent to this face "
2479 "dominates the other. See the documentation "
2480 "of this function for more information of how "
2481 "to deal with this situation."));
2482 used_mapping_index = dominated_fe_index;
2483 }
2484 }
2485
2486 // Same as if above, but when hp is enabled.
2487 if (sub_face_no == numbers::invalid_unsigned_int)
2488 {
2489 internal_hp_fe_face_values->reinit(
2490 cell, face_no, used_q_index, used_mapping_index, active_fe_index);
2491 fe_face_values = &const_cast<FEFaceValues<dim, spacedim> &>(
2492 internal_hp_fe_face_values->get_present_fe_values());
2493 }
2494 else
2495 {
2496 internal_hp_fe_subface_values->reinit(cell,
2497 face_no,
2498 sub_face_no,
2499 used_q_index,
2500 used_mapping_index,
2501 active_fe_index);
2502
2503 fe_face_values = &const_cast<FESubfaceValues<dim, spacedim> &>(
2504 internal_hp_fe_subface_values->get_present_fe_values());
2505 }
2506 if (sub_face_no_neighbor == numbers::invalid_unsigned_int)
2507 {
2508 internal_hp_fe_face_values_neighbor->reinit(cell_neighbor,
2509 face_no_neighbor,
2510 used_q_index,
2511 used_mapping_index,
2512 active_fe_index_neighbor);
2513
2514 fe_face_values_neighbor = &const_cast<FEFaceValues<dim, spacedim> &>(
2515 internal_hp_fe_face_values_neighbor->get_present_fe_values());
2516 }
2517 else
2518 {
2519 internal_hp_fe_subface_values_neighbor->reinit(
2520 cell_neighbor,
2521 face_no_neighbor,
2522 sub_face_no_neighbor,
2523 used_q_index,
2524 used_mapping_index,
2525 active_fe_index_neighbor);
2526
2527 fe_face_values_neighbor =
2528 &const_cast<FESubfaceValues<dim, spacedim> &>(
2529 internal_hp_fe_subface_values_neighbor->get_present_fe_values());
2530 }
2531
2532 AssertDimension(fe_face_values->n_quadrature_points,
2533 fe_face_values_neighbor->n_quadrature_points);
2534
2535 const_cast<unsigned int &>(this->n_quadrature_points) =
2536 fe_face_values->n_quadrature_points;
2537 }
2538
2539 // Set up dof mapping and remove duplicates (for continuous elements).
2540 if constexpr (is_dof_cell_accessor_neighbor && is_dof_cell_accessor)
2541 {
2542 // Get dof indices first:
2543 std::vector<types::global_dof_index> v(
2544 fe_face_values->get_fe().n_dofs_per_cell());
2545 cell->get_active_or_mg_dof_indices(v);
2546 std::vector<types::global_dof_index> v2(
2547 fe_face_values_neighbor->get_fe().n_dofs_per_cell());
2548 cell_neighbor->get_active_or_mg_dof_indices(v2);
2549
2550 // Fill a map from the global dof index to the left and right
2551 // local index.
2552 std::map<types::global_dof_index, std::pair<unsigned int, unsigned int>>
2553 tempmap;
2554 std::pair<unsigned int, unsigned int> invalid_entry(
2556
2557 for (unsigned int i = 0; i < v.size(); ++i)
2558 {
2559 // If not already existing, add an invalid entry:
2560 auto result = tempmap.insert(std::make_pair(v[i], invalid_entry));
2561 result.first->second.first = i;
2562 }
2563
2564 for (unsigned int i = 0; i < v2.size(); ++i)
2565 {
2566 // If not already existing, add an invalid entry:
2567 auto result = tempmap.insert(std::make_pair(v2[i], invalid_entry));
2568 result.first->second.second = i;
2569 }
2570
2571 // Transfer from the map to the sorted std::vectors.
2572 dofmap.resize(tempmap.size());
2573 interface_dof_indices.resize(tempmap.size());
2574 unsigned int idx = 0;
2575 for (auto &x : tempmap)
2576 {
2577 interface_dof_indices[idx] = x.first;
2578 dofmap[idx] = {{x.second.first, x.second.second}};
2579 ++idx;
2580 }
2581 }
2582 else
2583 {
2584 const unsigned int n_dofs_per_cell_1 = fe_face_values->dofs_per_cell;
2585 const unsigned int n_dofs_per_cell_2 =
2586 fe_face_values_neighbor->dofs_per_cell;
2587
2588 interface_dof_indices.resize(n_dofs_per_cell_1 + n_dofs_per_cell_2);
2589 dofmap.resize(n_dofs_per_cell_1 + n_dofs_per_cell_2);
2590
2591 for (unsigned int i = 0; i < n_dofs_per_cell_1; ++i)
2592 {
2593 interface_dof_indices[i] = numbers::invalid_dof_index;
2594 dofmap[i] = {{i, numbers::invalid_unsigned_int}};
2595 }
2596
2597 for (unsigned int i = 0; i < n_dofs_per_cell_2; ++i)
2598 {
2599 interface_dof_indices[i + n_dofs_per_cell_1] =
2601 dofmap[i + n_dofs_per_cell_1] = {{numbers::invalid_unsigned_int, i}};
2602 }
2603 }
2604}
2605
2606
2607
2608template <int dim, int spacedim>
2609template <typename CellIteratorType>
2610void
2611FEInterfaceValues<dim, spacedim>::reinit(const CellIteratorType &cell,
2612 const unsigned int face_no,
2613 const unsigned int q_index,
2614 const unsigned int mapping_index,
2615 const unsigned int fe_index)
2616{
2617 Assert(internal_fe_face_values || internal_hp_fe_face_values,
2619
2620 if (internal_fe_face_values)
2621 {
2622 Assert((q_index == 0 || q_index == numbers::invalid_unsigned_int),
2624 Assert((mapping_index == 0 ||
2625 mapping_index == numbers::invalid_unsigned_int),
2627 Assert((fe_index == 0 || fe_index == numbers::invalid_unsigned_int),
2629
2630 internal_fe_face_values->reinit(cell, face_no);
2631 fe_face_values = internal_fe_face_values.get();
2632 fe_face_values_neighbor = nullptr;
2633 }
2634 else if (internal_hp_fe_face_values)
2635 {
2636 internal_hp_fe_face_values->reinit(
2637 cell, face_no, q_index, mapping_index, fe_index);
2638 fe_face_values = &const_cast<FEFaceValues<dim> &>(
2639 internal_hp_fe_face_values->get_present_fe_values());
2640 fe_face_values_neighbor = nullptr;
2641 }
2642
2643 interface_dof_indices.resize(fe_face_values->get_fe().n_dofs_per_cell());
2644
2645 if constexpr (std::is_same_v<typename CellIteratorType::AccessorType,
2647 std::is_same_v<typename CellIteratorType::AccessorType,
2649 {
2650 cell->get_active_or_mg_dof_indices(interface_dof_indices);
2651 }
2652 else
2653 {
2654 for (auto &i : interface_dof_indices)
2656 }
2657
2658 dofmap.resize(interface_dof_indices.size());
2659 for (unsigned int i = 0; i < interface_dof_indices.size(); ++i)
2660 dofmap[i] = {{i, numbers::invalid_unsigned_int}};
2661}
2662
2663
2664
2665template <int dim, int spacedim>
2666inline double
2667FEInterfaceValues<dim, spacedim>::JxW(const unsigned int q) const
2668{
2669 Assert(fe_face_values != nullptr,
2670 ExcMessage("This call requires a call to reinit() first."));
2671 return fe_face_values->JxW(q);
2672}
2673
2674
2675
2676template <int dim, int spacedim>
2677const std::vector<double> &
2679{
2680 Assert(fe_face_values != nullptr,
2681 ExcMessage("This call requires a call to reinit() first."));
2682 return fe_face_values->get_JxW_values();
2683}
2684
2685
2686
2687template <int dim, int spacedim>
2688const std::vector<Tensor<1, spacedim>> &
2690{
2691 Assert(fe_face_values != nullptr,
2692 ExcMessage("This call requires a call to reinit() first."));
2693 return fe_face_values->get_normal_vectors();
2694}
2695
2696
2697
2698template <int dim, int spacedim>
2701{
2702 Assert(!has_hp_capabilities(), ExcOnlyAvailableWithoutHP());
2703 return internal_fe_face_values->get_mapping();
2704}
2705
2706
2707
2708template <int dim, int spacedim>
2711{
2712 Assert(!has_hp_capabilities(), ExcOnlyAvailableWithoutHP());
2713 return internal_fe_face_values->get_fe();
2714}
2715
2716
2717
2718template <int dim, int spacedim>
2719const Quadrature<dim - 1> &
2721{
2722 Assert(!has_hp_capabilities(), ExcOnlyAvailableWithoutHP());
2723 return internal_fe_face_values->get_quadrature();
2724}
2725
2726
2727
2728template <int dim, int spacedim>
2731{
2732 Assert(has_hp_capabilities(), ExcOnlyAvailableWithHP());
2733 return internal_hp_fe_face_values->get_mapping_collection();
2734}
2735
2736
2737
2738template <int dim, int spacedim>
2741{
2742 Assert(has_hp_capabilities(), ExcOnlyAvailableWithHP());
2743 return internal_hp_fe_face_values->get_fe_collection();
2744}
2745
2746
2747
2748template <int dim, int spacedim>
2749const hp::QCollection<dim - 1> &
2751{
2752 Assert(has_hp_capabilities(), ExcOnlyAvailableWithHP());
2753 return internal_hp_fe_face_values->get_quadrature_collection();
2754}
2755
2756
2757
2758template <int dim, int spacedim>
2759bool
2761{
2762 if (internal_hp_fe_face_values || internal_hp_fe_subface_values ||
2763 internal_hp_fe_face_values_neighbor ||
2764 internal_hp_fe_subface_values_neighbor)
2765 {
2766 Assert(!internal_fe_face_values, ExcInternalError());
2767 Assert(!internal_fe_subface_values, ExcInternalError());
2768 Assert(!internal_fe_face_values_neighbor, ExcInternalError());
2769 Assert(!internal_fe_subface_values_neighbor, ExcInternalError());
2770
2771 return true;
2772 }
2773
2774 Assert(internal_fe_face_values || internal_fe_subface_values ||
2775 internal_fe_face_values_neighbor ||
2776 internal_fe_subface_values_neighbor,
2778 Assert(!internal_hp_fe_face_values, ExcInternalError());
2779 Assert(!internal_hp_fe_subface_values, ExcInternalError());
2780 Assert(!internal_hp_fe_face_values_neighbor, ExcInternalError());
2781 Assert(!internal_hp_fe_subface_values_neighbor, ExcInternalError());
2782
2783 return false;
2784}
2785
2786
2787
2788template <int dim, int spacedim>
2791{
2793 0U, n_quadrature_points);
2794}
2795
2796
2797
2798template <int dim, int spacedim>
2799const Point<spacedim> &
2801 const unsigned int q_point) const
2802{
2803 Assert(fe_face_values != nullptr,
2804 ExcMessage("This call requires a call to reinit() first."));
2805 return fe_face_values->quadrature_point(q_point);
2806}
2807
2808
2809
2810template <int dim, int spacedim>
2811const std::vector<Point<spacedim>> &
2813{
2814 Assert(fe_face_values != nullptr,
2815 ExcMessage("This call requires a call to reinit() first."));
2816 return fe_face_values->get_quadrature_points();
2817}
2818
2819
2820
2821template <int dim, int spacedim>
2824{
2825 if (has_hp_capabilities())
2826 return internal_hp_fe_face_values->get_update_flags();
2827 else
2828 return internal_fe_face_values->get_update_flags();
2829}
2830
2831
2832
2833template <int dim, int spacedim>
2836{
2837 return get_fe_face_values(cell_index).get_cell();
2838}
2839
2840
2841
2842template <int dim, int spacedim>
2843inline unsigned int
2845 const unsigned int cell_index) const
2846{
2847 return get_fe_face_values(cell_index).get_face_number();
2848}
2849
2850
2851
2852template <int dim, int spacedim>
2853unsigned
2855{
2856 Assert(
2857 interface_dof_indices.size() > 0,
2858 ExcMessage(
2859 "n_current_interface_dofs() is only available after a call to reinit()."));
2860 return interface_dof_indices.size();
2861}
2862
2863
2864
2865template <int dim, int spacedim>
2868{
2869 return {0U, n_current_interface_dofs()};
2870}
2871
2872
2873
2874template <int dim, int spacedim>
2875bool
2877{
2878 return fe_face_values_neighbor == nullptr;
2879}
2880
2881
2882
2883template <int dim, int spacedim>
2884std::vector<types::global_dof_index>
2886{
2887 return interface_dof_indices;
2888}
2889
2890
2891
2892template <int dim, int spacedim>
2893std::array<unsigned int, 2>
2895 const unsigned int interface_dof_index) const
2896{
2897 AssertIndexRange(interface_dof_index, n_current_interface_dofs());
2898 return dofmap[interface_dof_index];
2899}
2900
2901
2902
2903template <int dim, int spacedim>
2906 const unsigned int cell_index) const
2907{
2909 Assert(
2910 cell_index == 0 || !at_boundary(),
2911 ExcMessage(
2912 "You are on a boundary, so you can only ask for the first FEFaceValues object."));
2913
2914 return (cell_index == 0) ? *fe_face_values : *fe_face_values_neighbor;
2915}
2916
2917
2918
2919template <int dim, int spacedim>
2921FEInterfaceValues<dim, spacedim>::normal(const unsigned int q_point_index) const
2922{
2923 return normal_vector(q_point_index);
2924}
2925
2926
2927
2928template <int dim, int spacedim>
2931 const unsigned int q_point_index) const
2932{
2933 return fe_face_values->normal_vector(q_point_index);
2934}
2935
2936
2937
2938template <int dim, int spacedim>
2939double
2941 const bool here_or_there,
2942 const unsigned int interface_dof_index,
2943 const unsigned int q_point,
2944 const unsigned int component) const
2945{
2946 const auto dof_pair = dofmap[interface_dof_index];
2947
2948 if (here_or_there && dof_pair[0] != numbers::invalid_unsigned_int)
2949 return get_fe_face_values(0).shape_value_component(dof_pair[0],
2950 q_point,
2951 component);
2952 if (!here_or_there && dof_pair[1] != numbers::invalid_unsigned_int)
2953 return get_fe_face_values(1).shape_value_component(dof_pair[1],
2954 q_point,
2955 component);
2956
2957 return 0.0;
2958}
2959
2960
2961
2962template <int dim, int spacedim>
2963double
2965 const unsigned int interface_dof_index,
2966 const unsigned int q_point,
2967 const unsigned int component) const
2968{
2969 const auto dof_pair = dofmap[interface_dof_index];
2970
2971 double value = 0.0;
2972
2973 if (dof_pair[0] != numbers::invalid_unsigned_int)
2974 value += get_fe_face_values(0).shape_value_component(dof_pair[0],
2975 q_point,
2976 component);
2977 if (dof_pair[1] != numbers::invalid_unsigned_int)
2978 value -= get_fe_face_values(1).shape_value_component(dof_pair[1],
2979 q_point,
2980 component);
2981 return value;
2982}
2983
2984
2985
2986template <int dim, int spacedim>
2987double
2989 const unsigned int interface_dof_index,
2990 const unsigned int q_point,
2991 const unsigned int component) const
2992{
2993 const auto dof_pair = dofmap[interface_dof_index];
2994
2995 if (at_boundary())
2996 return get_fe_face_values(0).shape_value_component(dof_pair[0],
2997 q_point,
2998 component);
2999
3000 double value = 0.0;
3001
3002 if (dof_pair[0] != numbers::invalid_unsigned_int)
3003 value += 0.5 * get_fe_face_values(0).shape_value_component(dof_pair[0],
3004 q_point,
3005 component);
3006 if (dof_pair[1] != numbers::invalid_unsigned_int)
3007 value += 0.5 * get_fe_face_values(1).shape_value_component(dof_pair[1],
3008 q_point,
3009 component);
3010
3011 return value;
3012}
3013
3014
3015
3016template <int dim, int spacedim>
3019 const unsigned int interface_dof_index,
3020 const unsigned int q_point,
3021 const unsigned int component) const
3022{
3023 const auto dof_pair = dofmap[interface_dof_index];
3024
3025 if (at_boundary())
3026 return get_fe_face_values(0).shape_grad_component(dof_pair[0],
3027 q_point,
3028 component);
3029
3031
3032 if (dof_pair[0] != numbers::invalid_unsigned_int)
3033 value += 0.5 * get_fe_face_values(0).shape_grad_component(dof_pair[0],
3034 q_point,
3035 component);
3036 if (dof_pair[1] != numbers::invalid_unsigned_int)
3037 value += 0.5 * get_fe_face_values(1).shape_grad_component(dof_pair[1],
3038 q_point,
3039 component);
3040
3041 return value;
3042}
3043
3044
3045
3046template <int dim, int spacedim>
3049 const unsigned int interface_dof_index,
3050 const unsigned int q_point,
3051 const unsigned int component) const
3052{
3053 const auto dof_pair = dofmap[interface_dof_index];
3054
3055 if (at_boundary())
3056 return get_fe_face_values(0).shape_hessian_component(dof_pair[0],
3057 q_point,
3058 component);
3059
3061
3062 if (dof_pair[0] != numbers::invalid_unsigned_int)
3063 value += 0.5 * get_fe_face_values(0).shape_hessian_component(dof_pair[0],
3064 q_point,
3065 component);
3066 if (dof_pair[1] != numbers::invalid_unsigned_int)
3067 value += 0.5 * get_fe_face_values(1).shape_hessian_component(dof_pair[1],
3068 q_point,
3069 component);
3070
3071 return value;
3072}
3073
3074
3075
3076template <int dim, int spacedim>
3079 const unsigned int interface_dof_index,
3080 const unsigned int q_point,
3081 const unsigned int component) const
3082{
3083 const auto dof_pair = dofmap[interface_dof_index];
3084
3085 if (at_boundary())
3086 return get_fe_face_values(0).shape_grad_component(dof_pair[0],
3087 q_point,
3088 component);
3089
3091
3092 if (dof_pair[0] != numbers::invalid_unsigned_int)
3093 value += get_fe_face_values(0).shape_grad_component(dof_pair[0],
3094 q_point,
3095 component);
3096 if (dof_pair[1] != numbers::invalid_unsigned_int)
3097 value -= get_fe_face_values(1).shape_grad_component(dof_pair[1],
3098 q_point,
3099 component);
3100
3101 return value;
3102}
3103
3104
3105
3106template <int dim, int spacedim>
3109 const unsigned int interface_dof_index,
3110 const unsigned int q_point,
3111 const unsigned int component) const
3112{
3113 const auto dof_pair = dofmap[interface_dof_index];
3114
3115 if (at_boundary())
3116 return get_fe_face_values(0).shape_hessian_component(dof_pair[0],
3117 q_point,
3118 component);
3119
3121
3122 if (dof_pair[0] != numbers::invalid_unsigned_int)
3123 value += get_fe_face_values(0).shape_hessian_component(dof_pair[0],
3124 q_point,
3125 component);
3126 if (dof_pair[1] != numbers::invalid_unsigned_int)
3127 value -= get_fe_face_values(1).shape_hessian_component(dof_pair[1],
3128 q_point,
3129 component);
3130
3131 return value;
3132}
3133
3134
3135
3136template <int dim, int spacedim>
3139 const unsigned int interface_dof_index,
3140 const unsigned int q_point,
3141 const unsigned int component) const
3142{
3143 const auto dof_pair = dofmap[interface_dof_index];
3144
3145 if (at_boundary())
3146 return get_fe_face_values(0).shape_3rd_derivative_component(dof_pair[0],
3147 q_point,
3148 component);
3149
3151
3152 if (dof_pair[0] != numbers::invalid_unsigned_int)
3153 value += get_fe_face_values(0).shape_3rd_derivative_component(dof_pair[0],
3154 q_point,
3155 component);
3156 if (dof_pair[1] != numbers::invalid_unsigned_int)
3157 value -= get_fe_face_values(1).shape_3rd_derivative_component(dof_pair[1],
3158 q_point,
3159 component);
3160
3161 return value;
3162}
3163
3164
3165
3166template <int dim, int spacedim>
3167template <class InputVector>
3168void
3170 const InputVector &fe_function,
3171 std::vector<typename InputVector::value_type> &values) const
3172{
3173 AssertDimension(values.size(), n_quadrature_points);
3174
3176 this->operator[](scalar).get_jump_in_function_values(fe_function, values);
3177}
3178
3179
3180
3181template <int dim, int spacedim>
3182template <class InputVector>
3183void
3185 const InputVector &fe_function,
3187 const
3188{
3189 AssertDimension(gradients.size(), n_quadrature_points);
3190
3192 this->operator[](scalar).get_jump_in_function_gradients(fe_function,
3193 gradients);
3194}
3195
3196
3197
3198template <int dim, int spacedim>
3199template <class InputVector>
3200void
3202 const InputVector &fe_function,
3204 const
3205{
3206 AssertDimension(hessians.size(), n_quadrature_points);
3207
3209 this->operator[](scalar).get_jump_in_function_hessians(fe_function, hessians);
3210}
3211
3212
3213
3214template <int dim, int spacedim>
3215template <class InputVector>
3216void
3218 const InputVector &fe_function,
3220 &third_derivatives) const
3221{
3222 AssertDimension(third_derivatives.size(), n_quadrature_points);
3223
3225 this->operator[](scalar).get_jump_in_function_third_derivatives(
3226 fe_function, third_derivatives);
3227}
3228
3229
3230
3231template <int dim, int spacedim>
3232template <class InputVector>
3233void
3235 const InputVector &fe_function,
3236 std::vector<typename InputVector::value_type> &values) const
3237{
3238 AssertDimension(values.size(), n_quadrature_points);
3239
3241 this->operator[](scalar).get_average_of_function_values(fe_function, values);
3242}
3243
3244
3245
3246template <int dim, int spacedim>
3247template <class InputVector>
3248void
3250 const InputVector &fe_function,
3252 const
3253{
3254 AssertDimension(gradients.size(), n_quadrature_points);
3255
3257 this->operator[](scalar).get_average_of_function_gradients(fe_function,
3258 gradients);
3259}
3260
3261
3262
3263template <int dim, int spacedim>
3264template <class InputVector>
3265void
3267 const InputVector &fe_function,
3269 const
3270{
3271 AssertDimension(hessians.size(), n_quadrature_points);
3272
3274 this->operator[](scalar).get_average_of_function_hessians(fe_function,
3275 hessians);
3276}
3277
3278
3279
3280/*------------ Inline functions: FEInterfaceValues------------*/
3281template <int dim, int spacedim>
3284 const FEValuesExtractors::Scalar &scalar) const
3285{
3286 const unsigned int n_components =
3287 (this->has_hp_capabilities() ? this->get_fe_collection().n_components() :
3288 this->get_fe().n_components());
3289 (void)n_components;
3290 AssertIndexRange(scalar.component, n_components);
3291 return FEInterfaceViews::Scalar<dim, spacedim>(*this, scalar.component);
3292}
3293
3294
3295
3296template <int dim, int spacedim>
3299 const FEValuesExtractors::Vector &vector) const
3300{
3301 const unsigned int n_components =
3302 (this->has_hp_capabilities() ? this->get_fe_collection().n_components() :
3303 this->get_fe().n_components());
3304 const unsigned int n_vectors =
3307 0);
3308 (void)n_components;
3309 (void)n_vectors;
3310 AssertIndexRange(vector.first_vector_component, n_vectors);
3312 vector.first_vector_component);
3313}
3314
3315
3316
3317namespace FEInterfaceViews
3318{
3319 template <int dim, int spacedim>
3321 const FEInterfaceValues<dim, spacedim> &fe_interface)
3322 : fe_interface(&fe_interface)
3323 {}
3324
3325
3326
3327 template <int dim, int spacedim>
3328 Scalar<dim, spacedim>::Scalar(
3329 const FEInterfaceValues<dim, spacedim> &fe_interface,
3330 const unsigned int component)
3331 : Base<dim, spacedim>(fe_interface)
3332 , extractor(component)
3333 {}
3334
3335
3336
3337 template <int dim, int spacedim>
3338 template <class InputVector, class OutputVector>
3339 void
3340 Base<dim, spacedim>::get_local_dof_values(
3341 const InputVector &dof_values,
3342 OutputVector &local_dof_values) const
3343 {
3344 const auto &interface_dof_indices =
3345 this->fe_interface->get_interface_dof_indices();
3346
3347 AssertDimension(interface_dof_indices.size(), local_dof_values.size());
3348
3349 for (const unsigned int i : this->fe_interface->dof_indices())
3350 local_dof_values[i] = dof_values(interface_dof_indices[i]);
3351 }
3352
3353
3354
3355 template <int dim, int spacedim>
3356 typename Scalar<dim, spacedim>::value_type
3357 Scalar<dim, spacedim>::value(const bool here_or_there,
3358 const unsigned int interface_dof_index,
3359 const unsigned int q_point) const
3360 {
3361 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3362
3363 if (here_or_there && dof_pair[0] != numbers::invalid_unsigned_int)
3364 return (*(this->fe_interface->fe_face_values))[extractor].value(
3365 dof_pair[0], q_point);
3366
3367 if (!here_or_there && dof_pair[1] != numbers::invalid_unsigned_int)
3368 return (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3369 dof_pair[1], q_point);
3370
3371 return 0.0;
3372 }
3373
3374
3375
3376 template <int dim, int spacedim>
3377 typename Scalar<dim, spacedim>::value_type
3378 Scalar<dim, spacedim>::jump_in_values(const unsigned int interface_dof_index,
3379 const unsigned int q_point) const
3380 {
3381 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3382
3383 value_type value = 0.0;
3384
3385 if (dof_pair[0] != numbers::invalid_unsigned_int)
3386 value +=
3387 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3388 q_point);
3389
3390 if (dof_pair[1] != numbers::invalid_unsigned_int)
3391 value -=
3392 (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3393 dof_pair[1], q_point);
3394
3395 return value;
3396 }
3397
3398
3399
3400 template <int dim, int spacedim>
3401 typename Scalar<dim, spacedim>::value_type
3402 Scalar<dim, spacedim>::average_of_values(
3403 const unsigned int interface_dof_index,
3404 const unsigned int q_point) const
3405 {
3406 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3407
3408 if (this->fe_interface->at_boundary())
3409 return (*(this->fe_interface->fe_face_values))[extractor].value(
3410 dof_pair[0], q_point);
3411
3412 value_type value = 0.0;
3413
3414 if (dof_pair[0] != numbers::invalid_unsigned_int)
3415 value +=
3416 0.5 *
3417 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3418 q_point);
3419
3420 if (dof_pair[1] != numbers::invalid_unsigned_int)
3421 value +=
3422 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3423 dof_pair[1], q_point);
3424
3425 return value;
3426 }
3427
3428
3429
3430 template <int dim, int spacedim>
3431 typename Scalar<dim, spacedim>::gradient_type
3432 Scalar<dim, spacedim>::average_of_gradients(
3433 const unsigned int interface_dof_index,
3434 const unsigned int q_point) const
3435 {
3436 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3437
3438 if (this->fe_interface->at_boundary())
3439 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
3440 dof_pair[0], q_point);
3441
3442 gradient_type value;
3443
3444 if (dof_pair[0] != numbers::invalid_unsigned_int)
3445 value +=
3446 0.5 *
3447 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
3448 q_point);
3449
3450 if (dof_pair[1] != numbers::invalid_unsigned_int)
3451 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3452 .gradient(dof_pair[1], q_point);
3453
3454 return value;
3455 }
3456
3457
3458
3459 template <int dim, int spacedim>
3460 typename Scalar<dim, spacedim>::gradient_type
3461 Scalar<dim, spacedim>::jump_in_gradients(
3462 const unsigned int interface_dof_index,
3463 const unsigned int q_point) const
3464 {
3465 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3466
3467 if (this->fe_interface->at_boundary())
3468 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
3469 dof_pair[0], q_point);
3470
3471 gradient_type value;
3472
3473 if (dof_pair[0] != numbers::invalid_unsigned_int)
3474 value +=
3475 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
3476 q_point);
3477
3478 if (dof_pair[1] != numbers::invalid_unsigned_int)
3479 value -=
3480 (*(this->fe_interface->fe_face_values_neighbor))[extractor].gradient(
3481 dof_pair[1], q_point);
3482
3483 return value;
3484 }
3485
3486
3487
3488 template <int dim, int spacedim>
3489 typename Scalar<dim, spacedim>::hessian_type
3490 Scalar<dim, spacedim>::average_of_hessians(
3491 const unsigned int interface_dof_index,
3492 const unsigned int q_point) const
3493 {
3494 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3495
3496 if (this->fe_interface->at_boundary())
3497 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
3498 dof_pair[0], q_point);
3499
3500 hessian_type value;
3501
3502 if (dof_pair[0] != numbers::invalid_unsigned_int)
3503 value +=
3504 0.5 *
3505 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
3506 q_point);
3507
3508 if (dof_pair[1] != numbers::invalid_unsigned_int)
3509 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3510 .hessian(dof_pair[1], q_point);
3511
3512 return value;
3513 }
3514
3515
3516
3517 template <int dim, int spacedim>
3518 typename Scalar<dim, spacedim>::third_derivative_type
3519 Scalar<dim, spacedim>::jump_in_third_derivatives(
3520 const unsigned int interface_dof_index,
3521 const unsigned int q_point) const
3522 {
3523 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3524
3525 if (this->fe_interface->at_boundary())
3526 return (*(this->fe_interface->fe_face_values))[extractor]
3527 .third_derivative(dof_pair[0], q_point);
3528
3529 third_derivative_type value;
3530
3531 if (dof_pair[0] != numbers::invalid_unsigned_int)
3532 value +=
3533 (*(this->fe_interface->fe_face_values))[extractor].third_derivative(
3534 dof_pair[0], q_point);
3535
3536 if (dof_pair[1] != numbers::invalid_unsigned_int)
3537 value -= (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3538 .third_derivative(dof_pair[1], q_point);
3539
3540 return value;
3541 }
3542
3543
3544
3545 template <int dim, int spacedim>
3546 typename Scalar<dim, spacedim>::hessian_type
3547 Scalar<dim, spacedim>::jump_in_hessians(
3548 const unsigned int interface_dof_index,
3549 const unsigned int q_point) const
3550 {
3551 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3552
3553 if (this->fe_interface->at_boundary())
3554 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
3555 dof_pair[0], q_point);
3556
3557 hessian_type value;
3558
3559 if (dof_pair[0] != numbers::invalid_unsigned_int)
3560 value +=
3561 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
3562 q_point);
3563
3564 if (dof_pair[1] != numbers::invalid_unsigned_int)
3565 value -=
3566 (*(this->fe_interface->fe_face_values_neighbor))[extractor].hessian(
3567 dof_pair[1], q_point);
3568
3569 return value;
3570 }
3571
3572
3573
3574 template <int dim, int spacedim>
3575 template <class InputVector>
3576 void
3577 Scalar<dim, spacedim>::get_function_values_from_local_dof_values(
3578 const bool here_or_there,
3579 const InputVector &local_dof_values,
3580 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3581 const
3582 {
3583 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
3584
3585 for (const auto dof_index : this->fe_interface->dof_indices())
3586 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3587 {
3588 if (dof_index == 0)
3589 values[q_index] = 0.;
3590
3591 values[q_index] += local_dof_values[dof_index] *
3592 value(here_or_there, dof_index, q_index);
3593 }
3594 }
3595
3596
3597
3598 template <int dim, int spacedim>
3599 template <class InputVector>
3600 void
3601 Scalar<dim, spacedim>::get_function_values(
3602 const bool here_or_there,
3603 const InputVector &fe_function,
3604 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3605 const
3606 {
3607 std::vector<typename InputVector::value_type> local_dof_values(
3608 this->fe_interface->n_current_interface_dofs());
3609 this->get_local_dof_values(fe_function, local_dof_values);
3610
3611 get_function_values_from_local_dof_values(here_or_there,
3612 local_dof_values,
3613 values);
3614 }
3615
3616
3617
3618 template <int dim, int spacedim>
3619 template <class InputVector>
3620 void
3621 Scalar<dim, spacedim>::get_jump_in_function_values_from_local_dof_values(
3622 const InputVector &local_dof_values,
3623 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3624 const
3625 {
3626 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
3627
3628 for (const auto dof_index : this->fe_interface->dof_indices())
3629 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3630 {
3631 if (dof_index == 0)
3632 values[q_index] = 0.;
3633
3634 values[q_index] +=
3635 local_dof_values[dof_index] * jump_in_values(dof_index, q_index);
3636 }
3637 }
3638
3639
3640
3641 template <int dim, int spacedim>
3642 template <class InputVector>
3643 void
3644 Scalar<dim, spacedim>::get_jump_in_function_values(
3645 const InputVector &fe_function,
3646 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3647 const
3648 {
3649 std::vector<typename InputVector::value_type> local_dof_values(
3650 this->fe_interface->n_current_interface_dofs());
3651 this->get_local_dof_values(fe_function, local_dof_values);
3652
3653 get_jump_in_function_values_from_local_dof_values(local_dof_values, values);
3654 }
3655
3656
3657
3658 template <int dim, int spacedim>
3659 template <class InputVector>
3660 void
3661 Scalar<dim, spacedim>::get_jump_in_function_gradients_from_local_dof_values(
3662 const InputVector &local_dof_values,
3663 std::vector<solution_gradient_type<typename InputVector::value_type>>
3664 &gradients) const
3665 {
3666 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
3667
3668 for (const auto dof_index : this->fe_interface->dof_indices())
3669 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3670 {
3671 if (dof_index == 0)
3672 gradients[q_index] = 0.;
3673
3674 gradients[q_index] +=
3675 local_dof_values[dof_index] * jump_in_gradients(dof_index, q_index);
3676 }
3677 }
3678
3679
3680
3681 template <int dim, int spacedim>
3682 template <class InputVector>
3683 void
3684 Scalar<dim, spacedim>::get_jump_in_function_gradients(
3685 const InputVector &fe_function,
3686 std::vector<solution_gradient_type<typename InputVector::value_type>>
3687 &gradients) const
3688 {
3689 std::vector<typename InputVector::value_type> local_dof_values(
3690 this->fe_interface->n_current_interface_dofs());
3691 this->get_local_dof_values(fe_function, local_dof_values);
3692
3693 get_jump_in_function_gradients_from_local_dof_values(local_dof_values,
3694 gradients);
3695 }
3696
3697
3698
3699 template <int dim, int spacedim>
3700 template <class InputVector>
3701 void
3702 Scalar<dim, spacedim>::get_average_of_function_values_from_local_dof_values(
3703 const InputVector &local_dof_values,
3704 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3705 const
3706 {
3707 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
3708
3709 for (const auto dof_index : this->fe_interface->dof_indices())
3710 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3711 {
3712 if (dof_index == 0)
3713 values[q_index] = 0.;
3714
3715 values[q_index] +=
3716 local_dof_values[dof_index] * average_of_values(dof_index, q_index);
3717 }
3718 }
3719
3720
3721
3722 template <int dim, int spacedim>
3723 template <class InputVector>
3724 void
3725 Scalar<dim, spacedim>::get_average_of_function_values(
3726 const InputVector &fe_function,
3727 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3728 const
3729 {
3730 std::vector<typename InputVector::value_type> local_dof_values(
3731 this->fe_interface->n_current_interface_dofs());
3732 this->get_local_dof_values(fe_function, local_dof_values);
3733
3734 get_average_of_function_values_from_local_dof_values(local_dof_values,
3735 values);
3736 }
3737
3738
3739
3740 template <int dim, int spacedim>
3741 template <class InputVector>
3742 void
3743 Scalar<dim, spacedim>::
3744 get_average_of_function_gradients_from_local_dof_values(
3745 const InputVector &local_dof_values,
3746 std::vector<solution_gradient_type<typename InputVector::value_type>>
3747 &gradients) const
3748 {
3749 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
3750
3751 for (const auto dof_index : this->fe_interface->dof_indices())
3752 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3753 {
3754 if (dof_index == 0)
3755 gradients[q_index] = 0.;
3756
3757 gradients[q_index] += local_dof_values[dof_index] *
3758 average_of_gradients(dof_index, q_index);
3759 }
3760 }
3761
3762
3763
3764 template <int dim, int spacedim>
3765 template <class InputVector>
3766 void
3767 Scalar<dim, spacedim>::get_average_of_function_gradients(
3768 const InputVector &fe_function,
3769 std::vector<solution_gradient_type<typename InputVector::value_type>>
3770 &gradients) const
3771 {
3772 std::vector<typename InputVector::value_type> local_dof_values(
3773 this->fe_interface->n_current_interface_dofs());
3774 this->get_local_dof_values(fe_function, local_dof_values);
3775
3776 get_average_of_function_gradients_from_local_dof_values(local_dof_values,
3777 gradients);
3778 }
3779
3780
3781
3782 template <int dim, int spacedim>
3783 template <class InputVector>
3784 void
3785 Scalar<dim, spacedim>::get_jump_in_function_hessians_from_local_dof_values(
3786 const InputVector &local_dof_values,
3787 std::vector<solution_hessian_type<typename InputVector::value_type>>
3788 &hessians) const
3789 {
3790 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
3791
3792 for (const auto dof_index : this->fe_interface->dof_indices())
3793 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3794 {
3795 if (dof_index == 0)
3796 hessians[q_index] = 0.;
3797
3798 hessians[q_index] +=
3799 local_dof_values[dof_index] * jump_in_hessians(dof_index, q_index);
3800 }
3801 }
3802
3803
3804
3805 template <int dim, int spacedim>
3806 template <class InputVector>
3807 void
3808 Scalar<dim, spacedim>::get_jump_in_function_hessians(
3809 const InputVector &fe_function,
3810 std::vector<solution_hessian_type<typename InputVector::value_type>>
3811 &hessians) const
3812 {
3813 std::vector<typename InputVector::value_type> local_dof_values(
3814 this->fe_interface->n_current_interface_dofs());
3815 this->get_local_dof_values(fe_function, local_dof_values);
3816
3817 get_jump_in_function_hessians_from_local_dof_values(local_dof_values,
3818 hessians);
3819 }
3820
3821
3822
3823 template <int dim, int spacedim>
3824 template <class InputVector>
3825 void
3826 Scalar<dim, spacedim>::get_average_of_function_hessians_from_local_dof_values(
3827 const InputVector &local_dof_values,
3828 std::vector<solution_hessian_type<typename InputVector::value_type>>
3829 &hessians) const
3830 {
3831 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
3832
3833 for (const unsigned int dof_index : this->fe_interface->dof_indices())
3834 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3835 {
3836 if (dof_index == 0)
3837 hessians[q_index] = 0.;
3838
3839 hessians[q_index] += local_dof_values[dof_index] *
3840 average_of_hessians(dof_index, q_index);
3841 }
3842 }
3843
3844
3845
3846 template <int dim, int spacedim>
3847 template <class InputVector>
3848 void
3849 Scalar<dim, spacedim>::get_average_of_function_hessians(
3850 const InputVector &fe_function,
3851 std::vector<solution_hessian_type<typename InputVector::value_type>>
3852 &hessians) const
3853 {
3854 std::vector<typename InputVector::value_type> local_dof_values(
3855 this->fe_interface->n_current_interface_dofs());
3856 this->get_local_dof_values(fe_function, local_dof_values);
3857
3858 get_average_of_function_hessians_from_local_dof_values(local_dof_values,
3859 hessians);
3860 }
3861
3862
3863
3864 template <int dim, int spacedim>
3865 template <class InputVector>
3866 void
3867 Scalar<dim, spacedim>::
3868 get_jump_in_function_third_derivatives_from_local_dof_values(
3869 const InputVector &local_dof_values,
3870 std::vector<
3871 solution_third_derivative_type<typename InputVector::value_type>>
3872 &third_derivatives) const
3873 {
3874 AssertDimension(third_derivatives.size(),
3875 this->fe_interface->n_quadrature_points);
3876
3877 for (const unsigned int dof_index : this->fe_interface->dof_indices())
3878 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3879 {
3880 if (dof_index == 0)
3881 third_derivatives[q_index] = 0.;
3882
3883 third_derivatives[q_index] +=
3884 local_dof_values[dof_index] *
3885 jump_in_third_derivatives(dof_index, q_index);
3886 }
3887 }
3888
3889
3890
3891 template <int dim, int spacedim>
3892 template <class InputVector>
3893 void
3894 Scalar<dim, spacedim>::get_jump_in_function_third_derivatives(
3895 const InputVector &fe_function,
3896 std::vector<
3897 solution_third_derivative_type<typename InputVector::value_type>>
3898 &third_derivatives) const
3899 {
3900 std::vector<typename InputVector::value_type> local_dof_values(
3901 this->fe_interface->n_current_interface_dofs());
3902 this->get_local_dof_values(fe_function, local_dof_values);
3903
3904 get_jump_in_function_third_derivatives_from_local_dof_values(
3905 local_dof_values, third_derivatives);
3906 }
3907
3908
3909
3910 template <int dim, int spacedim>
3912 const FEInterfaceValues<dim, spacedim> &fe_interface,
3913 const unsigned int first_vector_component)
3914 : Base<dim, spacedim>(fe_interface)
3915 , extractor(first_vector_component)
3916 {}
3917
3918
3919
3920 template <int dim, int spacedim>
3922 Vector<dim, spacedim>::value(const bool here_or_there,
3923 const unsigned int interface_dof_index,
3924 const unsigned int q_point) const
3925 {
3926 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3927
3928 if (here_or_there && dof_pair[0] != numbers::invalid_unsigned_int)
3929 return (*(this->fe_interface->fe_face_values))[extractor].value(
3930 dof_pair[0], q_point);
3931
3932 if (!here_or_there && dof_pair[1] != numbers::invalid_unsigned_int)
3933 return (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3934 dof_pair[1], q_point);
3935
3936 return value_type();
3937 }
3938
3939
3940
3941 template <int dim, int spacedim>
3943 Vector<dim, spacedim>::jump_in_values(const unsigned int interface_dof_index,
3944 const unsigned int q_point) const
3945 {
3946 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3947
3948 value_type value;
3949
3950 if (dof_pair[0] != numbers::invalid_unsigned_int)
3951 value +=
3952 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3953 q_point);
3954
3955 if (dof_pair[1] != numbers::invalid_unsigned_int)
3956 value -=
3957 (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3958 dof_pair[1], q_point);
3959
3960 return value;
3961 }
3962
3963
3964
3965 template <int dim, int spacedim>
3968 const unsigned int interface_dof_index,
3969 const unsigned int q_point) const
3970 {
3971 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3972
3973 if (this->fe_interface->at_boundary())
3974 return (*(this->fe_interface->fe_face_values))[extractor].value(
3975 dof_pair[0], q_point);
3976
3977 value_type value;
3978
3979 if (dof_pair[0] != numbers::invalid_unsigned_int)
3980 value +=
3981 0.5 *
3982 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3983 q_point);
3984
3985 if (dof_pair[1] != numbers::invalid_unsigned_int)
3986 value +=
3987 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3988 dof_pair[1], q_point);
3989
3990 return value;
3991 }
3992
3993
3994
3995 template <int dim, int spacedim>
3998 const unsigned int interface_dof_index,
3999 const unsigned int q_point) const
4000 {
4001 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4002
4003 if (this->fe_interface->at_boundary())
4004 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
4005 dof_pair[0], q_point);
4006
4007 gradient_type value;
4008
4009 if (dof_pair[0] != numbers::invalid_unsigned_int)
4010 value +=
4011 0.5 *
4012 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
4013 q_point);
4014
4015 if (dof_pair[1] != numbers::invalid_unsigned_int)
4016 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
4017 .gradient(dof_pair[1], q_point);
4018
4019 return value;
4020 }
4021
4022
4023
4024 template <int dim, int spacedim>
4027 const unsigned int interface_dof_index,
4028 const unsigned int q_point) const
4029 {
4030 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4031
4032 if (this->fe_interface->at_boundary())
4033 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
4034 dof_pair[0], q_point);
4035
4036 gradient_type value;
4037
4038 if (dof_pair[0] != numbers::invalid_unsigned_int)
4039 value +=
4040 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
4041 q_point);
4042
4043 if (dof_pair[1] != numbers::invalid_unsigned_int)
4044 value -=
4045 (*(this->fe_interface->fe_face_values_neighbor))[extractor].gradient(
4046 dof_pair[1], q_point);
4047
4048 return value;
4049 }
4050
4051
4052
4053 template <int dim, int spacedim>
4055 Vector<dim, spacedim>::jump_gradient(const unsigned int interface_dof_index,
4056 const unsigned int q_point) const
4057 {
4058 return jump_in_gradients(interface_dof_index, q_point);
4059 }
4060
4061
4062
4063 template <int dim, int spacedim>
4066 const unsigned int interface_dof_index,
4067 const unsigned int q_point) const
4068 {
4069 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4070
4071 if (this->fe_interface->at_boundary())
4072 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
4073 dof_pair[0], q_point);
4074
4075 hessian_type value;
4076
4077 if (dof_pair[0] != numbers::invalid_unsigned_int)
4078 value +=
4079 0.5 *
4080 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
4081 q_point);
4082
4083 if (dof_pair[1] != numbers::invalid_unsigned_int)
4084 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
4085 .hessian(dof_pair[1], q_point);
4086
4087 return value;
4088 }
4089
4090
4091
4092 template <int dim, int spacedim>
4094 Vector<dim, spacedim>::average_hessian(const unsigned int interface_dof_index,
4095 const unsigned int q_point) const
4096 {
4097 return average_of_hessians(interface_dof_index, q_point);
4098 }
4099
4100
4101
4102 template <int dim, int spacedim>
4105 const unsigned int interface_dof_index,
4106 const unsigned int q_point) const
4107 {
4108 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4109
4110 if (this->fe_interface->at_boundary())
4111 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
4112 dof_pair[0], q_point);
4113
4114 hessian_type value;
4115
4116 if (dof_pair[0] != numbers::invalid_unsigned_int)
4117 value +=
4118 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
4119 q_point);
4120
4121 if (dof_pair[1] != numbers::invalid_unsigned_int)
4122 value -=
4123 (*(this->fe_interface->fe_face_values_neighbor))[extractor].hessian(
4124 dof_pair[1], q_point);
4125
4126 return value;
4127 }
4128
4129
4130
4131 template <int dim, int spacedim>
4133 Vector<dim, spacedim>::jump_hessian(const unsigned int interface_dof_index,
4134 const unsigned int q_point) const
4135 {
4136 return jump_in_hessians(interface_dof_index, q_point);
4137 }
4138
4139
4140
4141 template <int dim, int spacedim>
4144 const unsigned int interface_dof_index,
4145 const unsigned int q_point) const
4146 {
4147 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4148
4149 if (this->fe_interface->at_boundary())
4150 return (*(this->fe_interface->fe_face_values))[extractor]
4151 .third_derivative(dof_pair[0], q_point);
4152
4153 third_derivative_type value;
4154
4155 if (dof_pair[0] != numbers::invalid_unsigned_int)
4156 value +=
4157 (*(this->fe_interface->fe_face_values))[extractor].third_derivative(
4158 dof_pair[0], q_point);
4159
4160 if (dof_pair[1] != numbers::invalid_unsigned_int)
4161 value -= (*(this->fe_interface->fe_face_values_neighbor))[extractor]
4162 .third_derivative(dof_pair[1], q_point);
4163
4164 return value;
4165 }
4166
4167
4168
4169 template <int dim, int spacedim>
4170 template <class InputVector>
4171 void
4173 const bool here_or_there,
4174 const InputVector &local_dof_values,
4175 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4176 const
4177 {
4178 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
4179
4180 for (const auto dof_index : this->fe_interface->dof_indices())
4181 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4182 {
4183 if (dof_index == 0)
4184 values[q_index] = 0.;
4185
4186 values[q_index] += local_dof_values[dof_index] *
4187 value(here_or_there, dof_index, q_index);
4188 }
4189 }
4190
4191
4192
4193 template <int dim, int spacedim>
4194 template <class InputVector>
4195 void
4197 const bool here_or_there,
4198 const InputVector &fe_function,
4199 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4200 const
4201 {
4202 std::vector<typename InputVector::value_type> local_dof_values(
4203 this->fe_interface->n_current_interface_dofs());
4204 this->get_local_dof_values(fe_function, local_dof_values);
4205
4206 get_function_values_from_local_dof_values(here_or_there,
4207 local_dof_values,
4208 values);
4209 }
4210
4211
4212
4213 template <int dim, int spacedim>
4214 template <class InputVector>
4215 void
4217 const InputVector &local_dof_values,
4218 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4219 const
4220 {
4221 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
4222
4223 for (const auto dof_index : this->fe_interface->dof_indices())
4224 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4225 {
4226 if (dof_index == 0)
4227 values[q_index] = 0.;
4228
4229 values[q_index] +=
4230 local_dof_values[dof_index] * jump_in_values(dof_index, q_index);
4231 }
4232 }
4233
4234
4235
4236 template <int dim, int spacedim>
4237 template <class InputVector>
4238 void
4240 const InputVector &fe_function,
4241 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4242 const
4243 {
4244 std::vector<typename InputVector::value_type> local_dof_values(
4245 this->fe_interface->n_current_interface_dofs());
4246 this->get_local_dof_values(fe_function, local_dof_values);
4247
4248 get_jump_in_function_values_from_local_dof_values(local_dof_values, values);
4249 }
4250
4251
4252
4253 template <int dim, int spacedim>
4254 template <class InputVector>
4255 void
4257 const InputVector &local_dof_values,
4258 std::vector<solution_gradient_type<typename InputVector::value_type>>
4259 &gradients) const
4260 {
4261 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
4262
4263 for (const auto dof_index : this->fe_interface->dof_indices())
4264 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4265 {
4266 if (dof_index == 0)
4267 gradients[q_index] = 0.;
4268
4269 gradients[q_index] +=
4270 local_dof_values[dof_index] * jump_in_gradients(dof_index, q_index);
4271 }
4272 }
4273
4274
4275
4276 template <int dim, int spacedim>
4277 template <class InputVector>
4278 void
4280 const InputVector &fe_function,
4281 std::vector<solution_gradient_type<typename InputVector::value_type>>
4282 &gradients) const
4283 {
4284 std::vector<typename InputVector::value_type> local_dof_values(
4285 this->fe_interface->n_current_interface_dofs());
4286 this->get_local_dof_values(fe_function, local_dof_values);
4287
4288 get_jump_in_function_gradients_from_local_dof_values(local_dof_values,
4289 gradients);
4290 }
4291
4292
4293
4294 template <int dim, int spacedim>
4295 template <class InputVector>
4296 void
4298 const InputVector &local_dof_values,
4299 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4300 const
4301 {
4302 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
4303
4304 for (const auto dof_index : this->fe_interface->dof_indices())
4305 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4306 {
4307 if (dof_index == 0)
4308 values[q_index] = 0.;
4309
4310 values[q_index] +=
4311 local_dof_values[dof_index] * average_of_values(dof_index, q_index);
4312 }
4313 }
4314
4315
4316
4317 template <int dim, int spacedim>
4318 template <class InputVector>
4319 void
4321 const InputVector &fe_function,
4322 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4323 const
4324 {
4325 std::vector<typename InputVector::value_type> local_dof_values(
4326 this->fe_interface->n_current_interface_dofs());
4327 this->get_local_dof_values(fe_function, local_dof_values);
4328
4329 get_average_of_function_values_from_local_dof_values(local_dof_values,
4330 values);
4331 }
4332
4333
4334
4335 template <int dim, int spacedim>
4336 template <class InputVector>
4337 void
4340 const InputVector &local_dof_values,
4341 std::vector<solution_gradient_type<typename InputVector::value_type>>
4342 &gradients) const
4343 {
4344 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
4345
4346 for (const auto dof_index : this->fe_interface->dof_indices())
4347 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4348 {
4349 if (dof_index == 0)
4350 gradients[q_index] = 0.;
4351
4352 gradients[q_index] += local_dof_values[dof_index] *
4353 average_of_gradients(dof_index, q_index);
4354 }
4355 }
4356
4357
4358
4359 template <int dim, int spacedim>
4360 template <class InputVector>
4361 void
4363 const InputVector &fe_function,
4364 std::vector<solution_gradient_type<typename InputVector::value_type>>
4365 &gradients) const
4366 {
4367 std::vector<typename InputVector::value_type> local_dof_values(
4368 this->fe_interface->n_current_interface_dofs());
4369 this->get_local_dof_values(fe_function, local_dof_values);
4370
4371 get_average_of_function_gradients_from_local_dof_values(local_dof_values,
4372 gradients);
4373 }
4374
4375
4376
4377 template <int dim, int spacedim>
4378 template <class InputVector>
4379 void
4381 const InputVector &local_dof_values,
4382 std::vector<solution_hessian_type<typename InputVector::value_type>>
4383 &hessians) const
4384 {
4385 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
4386
4387 for (const auto dof_index : this->fe_interface->dof_indices())
4388 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4389 {
4390 if (dof_index == 0)
4391 hessians[q_index] = 0.;
4392
4393 hessians[q_index] +=
4394 local_dof_values[dof_index] * jump_in_hessians(dof_index, q_index);
4395 }
4396 }
4397
4398
4399
4400 template <int dim, int spacedim>
4401 template <class InputVector>
4402 void
4404 const InputVector &fe_function,
4405 std::vector<solution_hessian_type<typename InputVector::value_type>>
4406 &hessians) const
4407 {
4408 std::vector<typename InputVector::value_type> local_dof_values(
4409 this->fe_interface->n_current_interface_dofs());
4410 this->get_local_dof_values(fe_function, local_dof_values);
4411
4412 get_jump_in_function_hessians_from_local_dof_values(local_dof_values,
4413 hessians);
4414 }
4415
4416
4417
4418 template <int dim, int spacedim>
4419 template <class InputVector>
4420 void
4422 const InputVector &local_dof_values,
4423 std::vector<solution_hessian_type<typename InputVector::value_type>>
4424 &hessians) const
4425 {
4426 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
4427
4428 for (const unsigned int dof_index : this->fe_interface->dof_indices())
4429 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4430 {
4431 if (dof_index == 0)
4432 hessians[q_index] = 0.;
4433
4434 hessians[q_index] += local_dof_values[dof_index] *
4435 average_of_hessians(dof_index, q_index);
4436 }
4437 }
4438
4439
4440
4441 template <int dim, int spacedim>
4442 template <class InputVector>
4443 void
4445 const InputVector &fe_function,
4446 std::vector<solution_hessian_type<typename InputVector::value_type>>
4447 &hessians) const
4448 {
4449 std::vector<typename InputVector::value_type> local_dof_values(
4450 this->fe_interface->n_current_interface_dofs());
4451 this->get_local_dof_values(fe_function, local_dof_values);
4452
4453 get_average_of_function_hessians_from_local_dof_values(local_dof_values,
4454 hessians);
4455 }
4456
4457
4458
4459 template <int dim, int spacedim>
4460 template <class InputVector>
4461 void
4464 const InputVector &local_dof_values,
4465 std::vector<
4466 solution_third_derivative_type<typename InputVector::value_type>>
4467 &third_derivatives) const
4468 {
4469 AssertDimension(third_derivatives.size(),
4470 this->fe_interface->n_quadrature_points);
4471
4472 for (const unsigned int dof_index : this->fe_interface->dof_indices())
4473 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4474 {
4475 if (dof_index == 0)
4476 third_derivatives[q_index] = 0.;
4477
4478 third_derivatives[q_index] +=
4479 local_dof_values[dof_index] *
4480 jump_in_third_derivatives(dof_index, q_index);
4481 }
4482 }
4483
4484
4485
4486 template <int dim, int spacedim>
4487 template <class InputVector>
4488 void
4490 const InputVector &fe_function,
4491 std::vector<
4492 solution_third_derivative_type<typename InputVector::value_type>>
4493 &third_derivatives) const
4494 {
4495 std::vector<typename InputVector::value_type> local_dof_values(
4496 this->fe_interface->n_current_interface_dofs());
4497 this->get_local_dof_values(fe_function, local_dof_values);
4498
4499 get_jump_in_function_third_derivatives_from_local_dof_values(
4500 local_dof_values, third_derivatives);
4501 }
4502} // namespace FEInterfaceViews
4503
4504#endif // DOXYGEN
4505
4507
4508#endif
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access > > &cell, const unsigned int face_no)
friend class FEInterfaceViews::Scalar
std::unique_ptr< FEFaceValues< dim, spacedim > > internal_fe_face_values
UpdateFlags get_update_flags() const
const hp::MappingCollection< dim, spacedim > & get_mapping_collection() const
Tensor< 1, spacedim > average_of_shape_gradients(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
void get_average_of_function_gradients(const InputVector &fe_function, std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > &gradients) const
void get_jump_in_function_third_derivatives(const InputVector &fe_function, std::vector< Tensor< 3, spacedim, typename InputVector::value_type > > &third_derivatives) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices() const
FEFaceValuesBase< dim, spacedim > * fe_face_values_neighbor
bool at_boundary() const
const std::vector< double > & get_JxW_values() const
void get_average_of_function_hessians(const InputVector &fe_function, std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > &hessians) const
FEInterfaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim - 1 > &quadrature, const UpdateFlags update_flags)
void reinit(const CellIteratorType &cell, const unsigned int face_no, const unsigned int q_index=numbers::invalid_unsigned_int, const unsigned int mapping_index=numbers::invalid_unsigned_int, const unsigned int fe_index=numbers::invalid_unsigned_int)
const hp::FECollection< dim, spacedim > & get_fe_collection() const
Tensor< 3, spacedim > jump_in_shape_3rd_derivatives(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
FEInterfaceValues(const hp::MappingCollection< dim, spacedim > &mapping_collection, const hp::FECollection< dim, spacedim > &fe_collection, const hp::QCollection< dim - 1 > &quadrature_collection, const UpdateFlags update_flags)
std::unique_ptr< hp::FESubfaceValues< dim, spacedim > > internal_hp_fe_subface_values
const unsigned int n_quadrature_points
bool has_hp_capabilities() const
void reinit(const CellIteratorType &cell, const unsigned int face_no, const unsigned int sub_face_no, const CellNeighborIteratorType &cell_neighbor, const unsigned int face_no_neighbor, const unsigned int sub_face_no_neighbor, const unsigned int q_index=numbers::invalid_unsigned_int, const unsigned int mapping_index=numbers::invalid_unsigned_int, const unsigned int fe_index=numbers::invalid_unsigned_int, const unsigned int fe_index_neighbor=numbers::invalid_unsigned_int)
FEInterfaceViews::Scalar< dim, spacedim > operator[](const FEValuesExtractors::Scalar &scalar) const
unsigned n_current_interface_dofs() const
const std::vector< Tensor< 1, spacedim > > & get_normal_vectors() const
const FiniteElement< dim, spacedim > & get_fe() const
double average_of_shape_values(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
const std::vector< Point< spacedim > > & get_quadrature_points() const
std::unique_ptr< FESubfaceValues< dim, spacedim > > internal_fe_subface_values
double shape_value(const bool here_or_there, const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
const FEFaceValuesBase< dim, spacedim > & get_fe_face_values(const unsigned int cell_index) const
Tensor< 2, spacedim > average_of_shape_hessians(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
std::unique_ptr< hp::FESubfaceValues< dim, spacedim > > internal_hp_fe_subface_values_neighbor
FEInterfaceValues(const hp::FECollection< dim, spacedim > &fe_collection, const hp::QCollection< dim - 1 > &quadrature_collection, const UpdateFlags update_flags)
const hp::QCollection< dim - 1 > & get_quadrature_collection() const
FEInterfaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const hp::QCollection< dim - 1 > &quadrature, const UpdateFlags update_flags)
FEFaceValuesBase< dim, spacedim > * fe_face_values
std::vector< types::global_dof_index > interface_dof_indices
Tensor< 1, spacedim > normal_vector(const unsigned int q_point_index) const
std::unique_ptr< FESubfaceValues< dim, spacedim > > internal_fe_subface_values_neighbor
std::vector< std::array< unsigned int, 2 > > dofmap
Triangulation< dim, spacedim >::cell_iterator get_cell(const unsigned int cell_index) const
void get_jump_in_function_hessians(const InputVector &fe_function, std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > &hessians) const
double JxW(const unsigned int quadrature_point) const
Tensor< 1, spacedim > normal(const unsigned int q_point_index) const
void get_average_of_function_values(const InputVector &fe_function, std::vector< typename InputVector::value_type > &values) const
Tensor< 1, spacedim > jump_in_shape_gradients(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
Tensor< 2, spacedim > jump_in_shape_hessians(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
std::unique_ptr< hp::FEFaceValues< dim, spacedim > > internal_hp_fe_face_values_neighbor
std_cxx20::ranges::iota_view< unsigned int, unsigned int > quadrature_point_indices() const
double jump_in_shape_values(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
unsigned int get_face_number(const unsigned int cell_index) const
std::vector< types::global_dof_index > get_interface_dof_indices() const
FEInterfaceValues(const FiniteElement< dim, spacedim > &fe, const Quadrature< dim - 1 > &quadrature, const UpdateFlags update_flags)
std::array< unsigned int, 2 > interface_dof_to_dof_indices(const unsigned int interface_dof_index) const
std::unique_ptr< hp::FEFaceValues< dim, spacedim > > internal_hp_fe_face_values
std::unique_ptr< FEFaceValues< dim, spacedim > > internal_fe_face_values_neighbor
void get_jump_in_function_values(const InputVector &fe_function, std::vector< typename InputVector::value_type > &values) const
void get_jump_in_function_gradients(const InputVector &fe_function, std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > &gradients) const
const Mapping< dim, spacedim > & get_mapping() const
const Point< spacedim > & quadrature_point(const unsigned int q_point) const
const Quadrature< dim - 1 > & get_quadrature() const
void get_local_dof_values(const InputVector &dof_values, OutputVector &local_dof_values) const
Base(const FEInterfaceValues< dim, spacedim > &fe_interface)
const FEInterfaceValues< dim, spacedim > * fe_interface
void get_jump_in_function_third_derivatives(const InputVector &fe_function, std::vector< solution_third_derivative_type< typename InputVector::value_type > > &third_derivatives) const
void get_jump_in_function_gradients_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
gradient_type jump_in_gradients(const unsigned int interface_dof_index, const unsigned int q_point) const
void get_jump_in_function_hessians(const InputVector &fe_function, std::vector< solution_hessian_type< typename InputVector::value_type > > &hessians) const
hessian_type jump_in_hessians(const unsigned int interface_dof_index, const unsigned int q_point) const
const FEValuesExtractors::Scalar extractor
void get_jump_in_function_third_derivatives_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_third_derivative_type< typename InputVector::value_type > > &third_derivatives) const
typename FEValuesViews::Scalar< dim, spacedim >::gradient_type gradient_type
hessian_type average_of_hessians(const unsigned int interface_dof_index, const unsigned int q_point) const
void get_jump_in_function_gradients(const InputVector &fe_function, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
value_type value(const bool here_or_there, const unsigned int interface_dof_index, const unsigned int q_point) const
void get_average_of_function_hessians_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_hessian_type< typename InputVector::value_type > > &hessians) const
void get_average_of_function_values_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
void get_average_of_function_values(const InputVector &fe_function, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
typename ProductType< Number, value_type >::type solution_value_type
void get_average_of_function_gradients(const InputVector &fe_function, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
void get_jump_in_function_values_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
void get_jump_in_function_values(const InputVector &fe_function, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
Scalar(const FEInterfaceValues< dim, spacedim > &fe_interface, const unsigned int component)
typename FEValuesViews::Scalar< dim, spacedim >::hessian_type hessian_type
third_derivative_type jump_in_third_derivatives(const unsigned int interface_dof_index, const unsigned int q_point) const
gradient_type average_of_gradients(const unsigned int interface_dof_index, const unsigned int q_point) const
void get_average_of_function_hessians(const InputVector &fe_function, std::vector< solution_hessian_type< typename InputVector::value_type > > &hessians) const
void get_average_of_function_gradients_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
value_type average_of_values(const unsigned int interface_dof_index, const unsigned int q_point) const
void get_jump_in_function_hessians_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_hessian_type< typename InputVector::value_type > > &hessians) const
value_type jump_in_values(const unsigned int interface_dof_index, const unsigned int q_point) const
void get_function_values_from_local_dof_values(const bool here_or_there, const InputVector &local_dof_values, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
typename FEValuesViews::Scalar< dim, spacedim >::third_derivative_type third_derivative_type
void get_function_values(const bool here_or_there, const InputVector &fe_function, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
typename ProductType< Number, hessian_type >::type solution_hessian_type
typename ProductType< Number, gradient_type >::type solution_gradient_type
const FEValuesExtractors::Vector extractor
void get_average_of_function_gradients_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
void get_average_of_function_hessians(const InputVector &fe_function, std::vector< solution_hessian_type< typename InputVector::value_type > > &hessians) const
void get_jump_in_function_values(const InputVector &fe_function, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
void get_average_of_function_values_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
void get_function_values_from_local_dof_values(const bool here_or_there, const InputVector &local_dof_values, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
void get_average_of_function_hessians_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_hessian_type< typename InputVector::value_type > > &hessians) const
hessian_type average_hessian(const unsigned int interface_dof_index, const unsigned int q_point) const
Vector(const FEInterfaceValues< dim, spacedim > &fe_interface, const unsigned int first_vector_component)
void get_average_of_function_gradients(const InputVector &fe_function, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
void get_function_values(const bool here_or_there, const InputVector &fe_function, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
typename ProductType< Number, gradient_type >::type solution_gradient_type
hessian_type jump_hessian(const unsigned int interface_dof_index, const unsigned int q_point) const
typename ProductType< Number, value_type >::type solution_value_type
void get_average_of_function_values(const InputVector &fe_function, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
typename FEValuesViews::Vector< dim, spacedim >::value_type value_type
value_type value(const bool here_or_there, const unsigned int interface_dof_index, const unsigned int q_point) const
void get_jump_in_function_third_derivatives_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_third_derivative_type< typename InputVector::value_type > > &third_derivatives) const
void get_jump_in_function_gradients(const InputVector &fe_function, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
third_derivative_type jump_in_third_derivatives(const unsigned int interface_dof_index, const unsigned int q_point) const
void get_jump_in_function_gradients_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
void get_jump_in_function_hessians_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_hessian_type< typename InputVector::value_type > > &hessians) const
void get_jump_in_function_hessians(const InputVector &fe_function, std::vector< solution_hessian_type< typename InputVector::value_type > > &hessians) const
value_type jump_in_values(const unsigned int interface_dof_index, const unsigned int q_point) const
typename ProductType< Number, hessian_type >::type solution_hessian_type
void get_jump_in_function_third_derivatives(const InputVector &fe_function, std::vector< solution_third_derivative_type< typename InputVector::value_type > > &third_derivatives) const
gradient_type jump_gradient(const unsigned int interface_dof_index, const unsigned int q_point) const
value_type average_of_values(const unsigned int interface_dof_index, const unsigned int q_point) const
hessian_type jump_in_hessians(const unsigned int interface_dof_index, const unsigned int q_point) const
void get_jump_in_function_values_from_local_dof_values(const InputVector &local_dof_values, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
gradient_type average_of_gradients(const unsigned int interface_dof_index, const unsigned int q_point) const
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
typename FEValuesViews::Vector< dim, spacedim >::gradient_type gradient_type
gradient_type jump_in_gradients(const unsigned int interface_dof_index, const unsigned int q_point) const
typename FEValuesViews::Vector< dim, spacedim >::third_derivative_type third_derivative_type
hessian_type average_of_hessians(const unsigned int interface_dof_index, const unsigned int q_point) const
typename FEValuesViews::Vector< dim, spacedim >::hessian_type hessian_type
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access > > &cell, const unsigned int face_no, const unsigned int subface_no)
Abstract base class for mapping classes.
Definition mapping.h:320
Definition point.h:113
friend class Vector
Definition vector.h:1124
Number value_type
Definition vector.h:137
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:35
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
unsigned int cell_index
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcOnlyAvailableWithHP()
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcNotInitialized()
static ::ExceptionBase & ExcOnlyAvailableWithoutHP()
static ::ExceptionBase & ExcMessage(std::string arg1)
UpdateFlags
const Mapping< dim, spacedim > & get_default_linear_mapping(const Triangulation< dim, spacedim > &triangulation)
Definition mapping.cc:316
std::size_t size
Definition mpi.cc:745
typename ::internal::FEInterfaceViews::ViewType< dim, spacedim, Extractor >::type View
constexpr char U
Definition hp.h:117
constexpr types::global_dof_index invalid_dof_index
Definition types.h:269
constexpr unsigned int invalid_unsigned_int
Definition types.h:238
constexpr types::fe_index invalid_fe_index
Definition types.h:260
boost::integer_range< IncrementableType > iota_view
Definition iota_view.h:45
STL namespace.
Definition types.h:32
typename internal::ProductTypeImpl< std::decay_t< T >, std::decay_t< U > >::type type