deal.II version GIT relicensing-2206-gaa53ff9447 2024-12-02 09: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\}}\)
Loading...
Searching...
No Matches
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}) :
2457 numbers::invalid_unsigned_int);
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{
2792 return {0U, n_quadrature_points};
2793}
2794
2795
2796
2797template <int dim, int spacedim>
2798const Point<spacedim> &
2800 const unsigned int q_point) const
2801{
2802 Assert(fe_face_values != nullptr,
2803 ExcMessage("This call requires a call to reinit() first."));
2804 return fe_face_values->quadrature_point(q_point);
2805}
2806
2807
2808
2809template <int dim, int spacedim>
2810const std::vector<Point<spacedim>> &
2812{
2813 Assert(fe_face_values != nullptr,
2814 ExcMessage("This call requires a call to reinit() first."));
2815 return fe_face_values->get_quadrature_points();
2816}
2817
2818
2819
2820template <int dim, int spacedim>
2823{
2824 if (has_hp_capabilities())
2825 return internal_hp_fe_face_values->get_update_flags();
2826 else
2827 return internal_fe_face_values->get_update_flags();
2828}
2829
2830
2831
2832template <int dim, int spacedim>
2835{
2836 return get_fe_face_values(cell_index).get_cell();
2837}
2838
2839
2840
2841template <int dim, int spacedim>
2842inline unsigned int
2844 const unsigned int cell_index) const
2845{
2846 return get_fe_face_values(cell_index).get_face_number();
2847}
2848
2849
2850
2851template <int dim, int spacedim>
2852unsigned
2854{
2855 Assert(
2856 interface_dof_indices.size() > 0,
2857 ExcMessage(
2858 "n_current_interface_dofs() is only available after a call to reinit()."));
2859 return interface_dof_indices.size();
2860}
2861
2862
2863
2864template <int dim, int spacedim>
2867{
2868 return {0U, n_current_interface_dofs()};
2869}
2870
2871
2872
2873template <int dim, int spacedim>
2874bool
2876{
2877 return fe_face_values_neighbor == nullptr;
2878}
2879
2880
2881
2882template <int dim, int spacedim>
2883std::vector<types::global_dof_index>
2885{
2886 return interface_dof_indices;
2887}
2888
2889
2890
2891template <int dim, int spacedim>
2892std::array<unsigned int, 2>
2894 const unsigned int interface_dof_index) const
2895{
2896 AssertIndexRange(interface_dof_index, n_current_interface_dofs());
2897 return dofmap[interface_dof_index];
2898}
2899
2900
2901
2902template <int dim, int spacedim>
2905 const unsigned int cell_index) const
2906{
2908 Assert(
2909 cell_index == 0 || !at_boundary(),
2910 ExcMessage(
2911 "You are on a boundary, so you can only ask for the first FEFaceValues object."));
2912
2913 return (cell_index == 0) ? *fe_face_values : *fe_face_values_neighbor;
2914}
2915
2916
2917
2918template <int dim, int spacedim>
2920FEInterfaceValues<dim, spacedim>::normal(const unsigned int q_point_index) const
2921{
2922 return normal_vector(q_point_index);
2923}
2924
2925
2926
2927template <int dim, int spacedim>
2930 const unsigned int q_point_index) const
2931{
2932 return fe_face_values->normal_vector(q_point_index);
2933}
2934
2935
2936
2937template <int dim, int spacedim>
2938double
2940 const bool here_or_there,
2941 const unsigned int interface_dof_index,
2942 const unsigned int q_point,
2943 const unsigned int component) const
2944{
2945 const auto dof_pair = dofmap[interface_dof_index];
2946
2947 if (here_or_there && dof_pair[0] != numbers::invalid_unsigned_int)
2948 return get_fe_face_values(0).shape_value_component(dof_pair[0],
2949 q_point,
2950 component);
2951 if (!here_or_there && dof_pair[1] != numbers::invalid_unsigned_int)
2952 return get_fe_face_values(1).shape_value_component(dof_pair[1],
2953 q_point,
2954 component);
2955
2956 return 0.0;
2957}
2958
2959
2960
2961template <int dim, int spacedim>
2962double
2964 const unsigned int interface_dof_index,
2965 const unsigned int q_point,
2966 const unsigned int component) const
2967{
2968 const auto dof_pair = dofmap[interface_dof_index];
2969
2970 double value = 0.0;
2971
2972 if (dof_pair[0] != numbers::invalid_unsigned_int)
2973 value += get_fe_face_values(0).shape_value_component(dof_pair[0],
2974 q_point,
2975 component);
2976 if (dof_pair[1] != numbers::invalid_unsigned_int)
2977 value -= get_fe_face_values(1).shape_value_component(dof_pair[1],
2978 q_point,
2979 component);
2980 return value;
2981}
2982
2983
2984
2985template <int dim, int spacedim>
2986double
2988 const unsigned int interface_dof_index,
2989 const unsigned int q_point,
2990 const unsigned int component) const
2991{
2992 const auto dof_pair = dofmap[interface_dof_index];
2993
2994 if (at_boundary())
2995 return get_fe_face_values(0).shape_value_component(dof_pair[0],
2996 q_point,
2997 component);
2998
2999 double value = 0.0;
3000
3001 if (dof_pair[0] != numbers::invalid_unsigned_int)
3002 value += 0.5 * get_fe_face_values(0).shape_value_component(dof_pair[0],
3003 q_point,
3004 component);
3005 if (dof_pair[1] != numbers::invalid_unsigned_int)
3006 value += 0.5 * get_fe_face_values(1).shape_value_component(dof_pair[1],
3007 q_point,
3008 component);
3009
3010 return value;
3011}
3012
3013
3014
3015template <int dim, int spacedim>
3018 const unsigned int interface_dof_index,
3019 const unsigned int q_point,
3020 const unsigned int component) const
3021{
3022 const auto dof_pair = dofmap[interface_dof_index];
3023
3024 if (at_boundary())
3025 return get_fe_face_values(0).shape_grad_component(dof_pair[0],
3026 q_point,
3027 component);
3028
3030
3031 if (dof_pair[0] != numbers::invalid_unsigned_int)
3032 value += 0.5 * get_fe_face_values(0).shape_grad_component(dof_pair[0],
3033 q_point,
3034 component);
3035 if (dof_pair[1] != numbers::invalid_unsigned_int)
3036 value += 0.5 * get_fe_face_values(1).shape_grad_component(dof_pair[1],
3037 q_point,
3038 component);
3039
3040 return value;
3041}
3042
3043
3044
3045template <int dim, int spacedim>
3048 const unsigned int interface_dof_index,
3049 const unsigned int q_point,
3050 const unsigned int component) const
3051{
3052 const auto dof_pair = dofmap[interface_dof_index];
3053
3054 if (at_boundary())
3055 return get_fe_face_values(0).shape_hessian_component(dof_pair[0],
3056 q_point,
3057 component);
3058
3060
3061 if (dof_pair[0] != numbers::invalid_unsigned_int)
3062 value += 0.5 * get_fe_face_values(0).shape_hessian_component(dof_pair[0],
3063 q_point,
3064 component);
3065 if (dof_pair[1] != numbers::invalid_unsigned_int)
3066 value += 0.5 * get_fe_face_values(1).shape_hessian_component(dof_pair[1],
3067 q_point,
3068 component);
3069
3070 return value;
3071}
3072
3073
3074
3075template <int dim, int spacedim>
3078 const unsigned int interface_dof_index,
3079 const unsigned int q_point,
3080 const unsigned int component) const
3081{
3082 const auto dof_pair = dofmap[interface_dof_index];
3083
3084 if (at_boundary())
3085 return get_fe_face_values(0).shape_grad_component(dof_pair[0],
3086 q_point,
3087 component);
3088
3090
3091 if (dof_pair[0] != numbers::invalid_unsigned_int)
3092 value += get_fe_face_values(0).shape_grad_component(dof_pair[0],
3093 q_point,
3094 component);
3095 if (dof_pair[1] != numbers::invalid_unsigned_int)
3096 value -= get_fe_face_values(1).shape_grad_component(dof_pair[1],
3097 q_point,
3098 component);
3099
3100 return value;
3101}
3102
3103
3104
3105template <int dim, int spacedim>
3108 const unsigned int interface_dof_index,
3109 const unsigned int q_point,
3110 const unsigned int component) const
3111{
3112 const auto dof_pair = dofmap[interface_dof_index];
3113
3114 if (at_boundary())
3115 return get_fe_face_values(0).shape_hessian_component(dof_pair[0],
3116 q_point,
3117 component);
3118
3120
3121 if (dof_pair[0] != numbers::invalid_unsigned_int)
3122 value += get_fe_face_values(0).shape_hessian_component(dof_pair[0],
3123 q_point,
3124 component);
3125 if (dof_pair[1] != numbers::invalid_unsigned_int)
3126 value -= get_fe_face_values(1).shape_hessian_component(dof_pair[1],
3127 q_point,
3128 component);
3129
3130 return value;
3131}
3132
3133
3134
3135template <int dim, int spacedim>
3138 const unsigned int interface_dof_index,
3139 const unsigned int q_point,
3140 const unsigned int component) const
3141{
3142 const auto dof_pair = dofmap[interface_dof_index];
3143
3144 if (at_boundary())
3145 return get_fe_face_values(0).shape_3rd_derivative_component(dof_pair[0],
3146 q_point,
3147 component);
3148
3150
3151 if (dof_pair[0] != numbers::invalid_unsigned_int)
3152 value += get_fe_face_values(0).shape_3rd_derivative_component(dof_pair[0],
3153 q_point,
3154 component);
3155 if (dof_pair[1] != numbers::invalid_unsigned_int)
3156 value -= get_fe_face_values(1).shape_3rd_derivative_component(dof_pair[1],
3157 q_point,
3158 component);
3159
3160 return value;
3161}
3162
3163
3164
3165template <int dim, int spacedim>
3166template <class InputVector>
3167void
3169 const InputVector &fe_function,
3170 std::vector<typename InputVector::value_type> &values) const
3171{
3172 AssertDimension(values.size(), n_quadrature_points);
3173
3175 this->operator[](scalar).get_jump_in_function_values(fe_function, values);
3176}
3177
3178
3179
3180template <int dim, int spacedim>
3181template <class InputVector>
3182void
3184 const InputVector &fe_function,
3186 const
3187{
3188 AssertDimension(gradients.size(), n_quadrature_points);
3189
3191 this->operator[](scalar).get_jump_in_function_gradients(fe_function,
3192 gradients);
3193}
3194
3195
3196
3197template <int dim, int spacedim>
3198template <class InputVector>
3199void
3201 const InputVector &fe_function,
3203 const
3204{
3205 AssertDimension(hessians.size(), n_quadrature_points);
3206
3208 this->operator[](scalar).get_jump_in_function_hessians(fe_function, hessians);
3209}
3210
3211
3212
3213template <int dim, int spacedim>
3214template <class InputVector>
3215void
3217 const InputVector &fe_function,
3219 &third_derivatives) const
3220{
3221 AssertDimension(third_derivatives.size(), n_quadrature_points);
3222
3224 this->operator[](scalar).get_jump_in_function_third_derivatives(
3225 fe_function, third_derivatives);
3226}
3227
3228
3229
3230template <int dim, int spacedim>
3231template <class InputVector>
3232void
3234 const InputVector &fe_function,
3235 std::vector<typename InputVector::value_type> &values) const
3236{
3237 AssertDimension(values.size(), n_quadrature_points);
3238
3240 this->operator[](scalar).get_average_of_function_values(fe_function, values);
3241}
3242
3243
3244
3245template <int dim, int spacedim>
3246template <class InputVector>
3247void
3249 const InputVector &fe_function,
3251 const
3252{
3253 AssertDimension(gradients.size(), n_quadrature_points);
3254
3256 this->operator[](scalar).get_average_of_function_gradients(fe_function,
3257 gradients);
3258}
3259
3260
3261
3262template <int dim, int spacedim>
3263template <class InputVector>
3264void
3266 const InputVector &fe_function,
3268 const
3269{
3270 AssertDimension(hessians.size(), n_quadrature_points);
3271
3273 this->operator[](scalar).get_average_of_function_hessians(fe_function,
3274 hessians);
3275}
3276
3277
3278
3279/*------------ Inline functions: FEInterfaceValues------------*/
3280template <int dim, int spacedim>
3283 const FEValuesExtractors::Scalar &scalar) const
3284{
3285 const unsigned int n_components =
3286 (this->has_hp_capabilities() ? this->get_fe_collection().n_components() :
3287 this->get_fe().n_components());
3288 (void)n_components;
3289 AssertIndexRange(scalar.component, n_components);
3290 return FEInterfaceViews::Scalar<dim, spacedim>(*this, scalar.component);
3291}
3292
3293
3294
3295template <int dim, int spacedim>
3298 const FEValuesExtractors::Vector &vector) const
3299{
3300 const unsigned int n_components =
3301 (this->has_hp_capabilities() ? this->get_fe_collection().n_components() :
3302 this->get_fe().n_components());
3303 const unsigned int n_vectors =
3306 0);
3307 (void)n_components;
3308 (void)n_vectors;
3309 AssertIndexRange(vector.first_vector_component, n_vectors);
3311 vector.first_vector_component);
3312}
3313
3314
3315
3316namespace FEInterfaceViews
3317{
3318 template <int dim, int spacedim>
3320 const FEInterfaceValues<dim, spacedim> &fe_interface)
3321 : fe_interface(&fe_interface)
3322 {}
3323
3324
3325
3326 template <int dim, int spacedim>
3327 Scalar<dim, spacedim>::Scalar(
3328 const FEInterfaceValues<dim, spacedim> &fe_interface,
3329 const unsigned int component)
3330 : Base<dim, spacedim>(fe_interface)
3331 , extractor(component)
3332 {}
3333
3334
3335
3336 template <int dim, int spacedim>
3337 template <class InputVector, class OutputVector>
3338 void
3339 Base<dim, spacedim>::get_local_dof_values(
3340 const InputVector &dof_values,
3341 OutputVector &local_dof_values) const
3342 {
3343 const auto &interface_dof_indices =
3344 this->fe_interface->get_interface_dof_indices();
3345
3346 AssertDimension(interface_dof_indices.size(), local_dof_values.size());
3347
3348 for (const unsigned int i : this->fe_interface->dof_indices())
3349 local_dof_values[i] = dof_values(interface_dof_indices[i]);
3350 }
3351
3352
3353
3354 template <int dim, int spacedim>
3355 typename Scalar<dim, spacedim>::value_type
3356 Scalar<dim, spacedim>::value(const bool here_or_there,
3357 const unsigned int interface_dof_index,
3358 const unsigned int q_point) const
3359 {
3360 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3361
3362 if (here_or_there && dof_pair[0] != numbers::invalid_unsigned_int)
3363 return (*(this->fe_interface->fe_face_values))[extractor].value(
3364 dof_pair[0], q_point);
3365
3366 if (!here_or_there && dof_pair[1] != numbers::invalid_unsigned_int)
3367 return (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3368 dof_pair[1], q_point);
3369
3370 return 0.0;
3371 }
3372
3373
3374
3375 template <int dim, int spacedim>
3376 typename Scalar<dim, spacedim>::value_type
3377 Scalar<dim, spacedim>::jump_in_values(const unsigned int interface_dof_index,
3378 const unsigned int q_point) const
3379 {
3380 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3381
3382 value_type value = 0.0;
3383
3384 if (dof_pair[0] != numbers::invalid_unsigned_int)
3385 value +=
3386 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3387 q_point);
3388
3389 if (dof_pair[1] != numbers::invalid_unsigned_int)
3390 value -=
3391 (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3392 dof_pair[1], q_point);
3393
3394 return value;
3395 }
3396
3397
3398
3399 template <int dim, int spacedim>
3400 typename Scalar<dim, spacedim>::value_type
3401 Scalar<dim, spacedim>::average_of_values(
3402 const unsigned int interface_dof_index,
3403 const unsigned int q_point) const
3404 {
3405 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3406
3407 if (this->fe_interface->at_boundary())
3408 return (*(this->fe_interface->fe_face_values))[extractor].value(
3409 dof_pair[0], q_point);
3410
3411 value_type value = 0.0;
3412
3413 if (dof_pair[0] != numbers::invalid_unsigned_int)
3414 value +=
3415 0.5 *
3416 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3417 q_point);
3418
3419 if (dof_pair[1] != numbers::invalid_unsigned_int)
3420 value +=
3421 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3422 dof_pair[1], q_point);
3423
3424 return value;
3425 }
3426
3427
3428
3429 template <int dim, int spacedim>
3430 typename Scalar<dim, spacedim>::gradient_type
3431 Scalar<dim, spacedim>::average_of_gradients(
3432 const unsigned int interface_dof_index,
3433 const unsigned int q_point) const
3434 {
3435 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3436
3437 if (this->fe_interface->at_boundary())
3438 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
3439 dof_pair[0], q_point);
3440
3441 gradient_type value;
3442
3443 if (dof_pair[0] != numbers::invalid_unsigned_int)
3444 value +=
3445 0.5 *
3446 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
3447 q_point);
3448
3449 if (dof_pair[1] != numbers::invalid_unsigned_int)
3450 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3451 .gradient(dof_pair[1], q_point);
3452
3453 return value;
3454 }
3455
3456
3457
3458 template <int dim, int spacedim>
3459 typename Scalar<dim, spacedim>::gradient_type
3460 Scalar<dim, spacedim>::jump_in_gradients(
3461 const unsigned int interface_dof_index,
3462 const unsigned int q_point) const
3463 {
3464 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3465
3466 if (this->fe_interface->at_boundary())
3467 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
3468 dof_pair[0], q_point);
3469
3470 gradient_type value;
3471
3472 if (dof_pair[0] != numbers::invalid_unsigned_int)
3473 value +=
3474 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
3475 q_point);
3476
3477 if (dof_pair[1] != numbers::invalid_unsigned_int)
3478 value -=
3479 (*(this->fe_interface->fe_face_values_neighbor))[extractor].gradient(
3480 dof_pair[1], q_point);
3481
3482 return value;
3483 }
3484
3485
3486
3487 template <int dim, int spacedim>
3488 typename Scalar<dim, spacedim>::hessian_type
3489 Scalar<dim, spacedim>::average_of_hessians(
3490 const unsigned int interface_dof_index,
3491 const unsigned int q_point) const
3492 {
3493 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3494
3495 if (this->fe_interface->at_boundary())
3496 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
3497 dof_pair[0], q_point);
3498
3499 hessian_type value;
3500
3501 if (dof_pair[0] != numbers::invalid_unsigned_int)
3502 value +=
3503 0.5 *
3504 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
3505 q_point);
3506
3507 if (dof_pair[1] != numbers::invalid_unsigned_int)
3508 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3509 .hessian(dof_pair[1], q_point);
3510
3511 return value;
3512 }
3513
3514
3515
3516 template <int dim, int spacedim>
3517 typename Scalar<dim, spacedim>::third_derivative_type
3518 Scalar<dim, spacedim>::jump_in_third_derivatives(
3519 const unsigned int interface_dof_index,
3520 const unsigned int q_point) const
3521 {
3522 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3523
3524 if (this->fe_interface->at_boundary())
3525 return (*(this->fe_interface->fe_face_values))[extractor]
3526 .third_derivative(dof_pair[0], q_point);
3527
3528 third_derivative_type value;
3529
3530 if (dof_pair[0] != numbers::invalid_unsigned_int)
3531 value +=
3532 (*(this->fe_interface->fe_face_values))[extractor].third_derivative(
3533 dof_pair[0], q_point);
3534
3535 if (dof_pair[1] != numbers::invalid_unsigned_int)
3536 value -= (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3537 .third_derivative(dof_pair[1], q_point);
3538
3539 return value;
3540 }
3541
3542
3543
3544 template <int dim, int spacedim>
3545 typename Scalar<dim, spacedim>::hessian_type
3546 Scalar<dim, spacedim>::jump_in_hessians(
3547 const unsigned int interface_dof_index,
3548 const unsigned int q_point) const
3549 {
3550 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3551
3552 if (this->fe_interface->at_boundary())
3553 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
3554 dof_pair[0], q_point);
3555
3556 hessian_type value;
3557
3558 if (dof_pair[0] != numbers::invalid_unsigned_int)
3559 value +=
3560 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
3561 q_point);
3562
3563 if (dof_pair[1] != numbers::invalid_unsigned_int)
3564 value -=
3565 (*(this->fe_interface->fe_face_values_neighbor))[extractor].hessian(
3566 dof_pair[1], q_point);
3567
3568 return value;
3569 }
3570
3571
3572
3573 template <int dim, int spacedim>
3574 template <class InputVector>
3575 void
3576 Scalar<dim, spacedim>::get_function_values_from_local_dof_values(
3577 const bool here_or_there,
3578 const InputVector &local_dof_values,
3579 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3580 const
3581 {
3582 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
3583
3584 for (const auto dof_index : this->fe_interface->dof_indices())
3585 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3586 {
3587 if (dof_index == 0)
3588 values[q_index] = 0.;
3589
3590 values[q_index] += local_dof_values[dof_index] *
3591 value(here_or_there, dof_index, q_index);
3592 }
3593 }
3594
3595
3596
3597 template <int dim, int spacedim>
3598 template <class InputVector>
3599 void
3600 Scalar<dim, spacedim>::get_function_values(
3601 const bool here_or_there,
3602 const InputVector &fe_function,
3603 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3604 const
3605 {
3606 std::vector<typename InputVector::value_type> local_dof_values(
3607 this->fe_interface->n_current_interface_dofs());
3608 this->get_local_dof_values(fe_function, local_dof_values);
3609
3610 get_function_values_from_local_dof_values(here_or_there,
3611 local_dof_values,
3612 values);
3613 }
3614
3615
3616
3617 template <int dim, int spacedim>
3618 template <class InputVector>
3619 void
3620 Scalar<dim, spacedim>::get_jump_in_function_values_from_local_dof_values(
3621 const InputVector &local_dof_values,
3622 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3623 const
3624 {
3625 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
3626
3627 for (const auto dof_index : this->fe_interface->dof_indices())
3628 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3629 {
3630 if (dof_index == 0)
3631 values[q_index] = 0.;
3632
3633 values[q_index] +=
3634 local_dof_values[dof_index] * jump_in_values(dof_index, q_index);
3635 }
3636 }
3637
3638
3639
3640 template <int dim, int spacedim>
3641 template <class InputVector>
3642 void
3643 Scalar<dim, spacedim>::get_jump_in_function_values(
3644 const InputVector &fe_function,
3645 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3646 const
3647 {
3648 std::vector<typename InputVector::value_type> local_dof_values(
3649 this->fe_interface->n_current_interface_dofs());
3650 this->get_local_dof_values(fe_function, local_dof_values);
3651
3652 get_jump_in_function_values_from_local_dof_values(local_dof_values, values);
3653 }
3654
3655
3656
3657 template <int dim, int spacedim>
3658 template <class InputVector>
3659 void
3660 Scalar<dim, spacedim>::get_jump_in_function_gradients_from_local_dof_values(
3661 const InputVector &local_dof_values,
3662 std::vector<solution_gradient_type<typename InputVector::value_type>>
3663 &gradients) const
3664 {
3665 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
3666
3667 for (const auto dof_index : this->fe_interface->dof_indices())
3668 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3669 {
3670 if (dof_index == 0)
3671 gradients[q_index] = 0.;
3672
3673 gradients[q_index] +=
3674 local_dof_values[dof_index] * jump_in_gradients(dof_index, q_index);
3675 }
3676 }
3677
3678
3679
3680 template <int dim, int spacedim>
3681 template <class InputVector>
3682 void
3683 Scalar<dim, spacedim>::get_jump_in_function_gradients(
3684 const InputVector &fe_function,
3685 std::vector<solution_gradient_type<typename InputVector::value_type>>
3686 &gradients) const
3687 {
3688 std::vector<typename InputVector::value_type> local_dof_values(
3689 this->fe_interface->n_current_interface_dofs());
3690 this->get_local_dof_values(fe_function, local_dof_values);
3691
3692 get_jump_in_function_gradients_from_local_dof_values(local_dof_values,
3693 gradients);
3694 }
3695
3696
3697
3698 template <int dim, int spacedim>
3699 template <class InputVector>
3700 void
3701 Scalar<dim, spacedim>::get_average_of_function_values_from_local_dof_values(
3702 const InputVector &local_dof_values,
3703 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3704 const
3705 {
3706 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
3707
3708 for (const auto dof_index : this->fe_interface->dof_indices())
3709 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3710 {
3711 if (dof_index == 0)
3712 values[q_index] = 0.;
3713
3714 values[q_index] +=
3715 local_dof_values[dof_index] * average_of_values(dof_index, q_index);
3716 }
3717 }
3718
3719
3720
3721 template <int dim, int spacedim>
3722 template <class InputVector>
3723 void
3724 Scalar<dim, spacedim>::get_average_of_function_values(
3725 const InputVector &fe_function,
3726 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3727 const
3728 {
3729 std::vector<typename InputVector::value_type> local_dof_values(
3730 this->fe_interface->n_current_interface_dofs());
3731 this->get_local_dof_values(fe_function, local_dof_values);
3732
3733 get_average_of_function_values_from_local_dof_values(local_dof_values,
3734 values);
3735 }
3736
3737
3738
3739 template <int dim, int spacedim>
3740 template <class InputVector>
3741 void
3742 Scalar<dim, spacedim>::
3743 get_average_of_function_gradients_from_local_dof_values(
3744 const InputVector &local_dof_values,
3745 std::vector<solution_gradient_type<typename InputVector::value_type>>
3746 &gradients) const
3747 {
3748 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
3749
3750 for (const auto dof_index : this->fe_interface->dof_indices())
3751 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3752 {
3753 if (dof_index == 0)
3754 gradients[q_index] = 0.;
3755
3756 gradients[q_index] += local_dof_values[dof_index] *
3757 average_of_gradients(dof_index, q_index);
3758 }
3759 }
3760
3761
3762
3763 template <int dim, int spacedim>
3764 template <class InputVector>
3765 void
3766 Scalar<dim, spacedim>::get_average_of_function_gradients(
3767 const InputVector &fe_function,
3768 std::vector<solution_gradient_type<typename InputVector::value_type>>
3769 &gradients) const
3770 {
3771 std::vector<typename InputVector::value_type> local_dof_values(
3772 this->fe_interface->n_current_interface_dofs());
3773 this->get_local_dof_values(fe_function, local_dof_values);
3774
3775 get_average_of_function_gradients_from_local_dof_values(local_dof_values,
3776 gradients);
3777 }
3778
3779
3780
3781 template <int dim, int spacedim>
3782 template <class InputVector>
3783 void
3784 Scalar<dim, spacedim>::get_jump_in_function_hessians_from_local_dof_values(
3785 const InputVector &local_dof_values,
3786 std::vector<solution_hessian_type<typename InputVector::value_type>>
3787 &hessians) const
3788 {
3789 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
3790
3791 for (const auto dof_index : this->fe_interface->dof_indices())
3792 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3793 {
3794 if (dof_index == 0)
3795 hessians[q_index] = 0.;
3796
3797 hessians[q_index] +=
3798 local_dof_values[dof_index] * jump_in_hessians(dof_index, q_index);
3799 }
3800 }
3801
3802
3803
3804 template <int dim, int spacedim>
3805 template <class InputVector>
3806 void
3807 Scalar<dim, spacedim>::get_jump_in_function_hessians(
3808 const InputVector &fe_function,
3809 std::vector<solution_hessian_type<typename InputVector::value_type>>
3810 &hessians) const
3811 {
3812 std::vector<typename InputVector::value_type> local_dof_values(
3813 this->fe_interface->n_current_interface_dofs());
3814 this->get_local_dof_values(fe_function, local_dof_values);
3815
3816 get_jump_in_function_hessians_from_local_dof_values(local_dof_values,
3817 hessians);
3818 }
3819
3820
3821
3822 template <int dim, int spacedim>
3823 template <class InputVector>
3824 void
3825 Scalar<dim, spacedim>::get_average_of_function_hessians_from_local_dof_values(
3826 const InputVector &local_dof_values,
3827 std::vector<solution_hessian_type<typename InputVector::value_type>>
3828 &hessians) const
3829 {
3830 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
3831
3832 for (const unsigned int dof_index : this->fe_interface->dof_indices())
3833 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3834 {
3835 if (dof_index == 0)
3836 hessians[q_index] = 0.;
3837
3838 hessians[q_index] += local_dof_values[dof_index] *
3839 average_of_hessians(dof_index, q_index);
3840 }
3841 }
3842
3843
3844
3845 template <int dim, int spacedim>
3846 template <class InputVector>
3847 void
3848 Scalar<dim, spacedim>::get_average_of_function_hessians(
3849 const InputVector &fe_function,
3850 std::vector<solution_hessian_type<typename InputVector::value_type>>
3851 &hessians) const
3852 {
3853 std::vector<typename InputVector::value_type> local_dof_values(
3854 this->fe_interface->n_current_interface_dofs());
3855 this->get_local_dof_values(fe_function, local_dof_values);
3856
3857 get_average_of_function_hessians_from_local_dof_values(local_dof_values,
3858 hessians);
3859 }
3860
3861
3862
3863 template <int dim, int spacedim>
3864 template <class InputVector>
3865 void
3866 Scalar<dim, spacedim>::
3867 get_jump_in_function_third_derivatives_from_local_dof_values(
3868 const InputVector &local_dof_values,
3869 std::vector<
3870 solution_third_derivative_type<typename InputVector::value_type>>
3871 &third_derivatives) const
3872 {
3873 AssertDimension(third_derivatives.size(),
3874 this->fe_interface->n_quadrature_points);
3875
3876 for (const unsigned int dof_index : this->fe_interface->dof_indices())
3877 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3878 {
3879 if (dof_index == 0)
3880 third_derivatives[q_index] = 0.;
3881
3882 third_derivatives[q_index] +=
3883 local_dof_values[dof_index] *
3884 jump_in_third_derivatives(dof_index, q_index);
3885 }
3886 }
3887
3888
3889
3890 template <int dim, int spacedim>
3891 template <class InputVector>
3892 void
3893 Scalar<dim, spacedim>::get_jump_in_function_third_derivatives(
3894 const InputVector &fe_function,
3895 std::vector<
3896 solution_third_derivative_type<typename InputVector::value_type>>
3897 &third_derivatives) const
3898 {
3899 std::vector<typename InputVector::value_type> local_dof_values(
3900 this->fe_interface->n_current_interface_dofs());
3901 this->get_local_dof_values(fe_function, local_dof_values);
3902
3903 get_jump_in_function_third_derivatives_from_local_dof_values(
3904 local_dof_values, third_derivatives);
3905 }
3906
3907
3908
3909 template <int dim, int spacedim>
3911 const FEInterfaceValues<dim, spacedim> &fe_interface,
3912 const unsigned int first_vector_component)
3913 : Base<dim, spacedim>(fe_interface)
3914 , extractor(first_vector_component)
3915 {}
3916
3917
3918
3919 template <int dim, int spacedim>
3921 Vector<dim, spacedim>::value(const bool here_or_there,
3922 const unsigned int interface_dof_index,
3923 const unsigned int q_point) const
3924 {
3925 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3926
3927 if (here_or_there && dof_pair[0] != numbers::invalid_unsigned_int)
3928 return (*(this->fe_interface->fe_face_values))[extractor].value(
3929 dof_pair[0], q_point);
3930
3931 if (!here_or_there && dof_pair[1] != numbers::invalid_unsigned_int)
3932 return (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3933 dof_pair[1], q_point);
3934
3935 return value_type();
3936 }
3937
3938
3939
3940 template <int dim, int spacedim>
3942 Vector<dim, spacedim>::jump_in_values(const unsigned int interface_dof_index,
3943 const unsigned int q_point) const
3944 {
3945 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3946
3947 value_type value;
3948
3949 if (dof_pair[0] != numbers::invalid_unsigned_int)
3950 value +=
3951 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3952 q_point);
3953
3954 if (dof_pair[1] != numbers::invalid_unsigned_int)
3955 value -=
3956 (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3957 dof_pair[1], q_point);
3958
3959 return value;
3960 }
3961
3962
3963
3964 template <int dim, int spacedim>
3967 const unsigned int interface_dof_index,
3968 const unsigned int q_point) const
3969 {
3970 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3971
3972 if (this->fe_interface->at_boundary())
3973 return (*(this->fe_interface->fe_face_values))[extractor].value(
3974 dof_pair[0], q_point);
3975
3976 value_type value;
3977
3978 if (dof_pair[0] != numbers::invalid_unsigned_int)
3979 value +=
3980 0.5 *
3981 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3982 q_point);
3983
3984 if (dof_pair[1] != numbers::invalid_unsigned_int)
3985 value +=
3986 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3987 dof_pair[1], q_point);
3988
3989 return value;
3990 }
3991
3992
3993
3994 template <int dim, int spacedim>
3997 const unsigned int interface_dof_index,
3998 const unsigned int q_point) const
3999 {
4000 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4001
4002 if (this->fe_interface->at_boundary())
4003 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
4004 dof_pair[0], q_point);
4005
4006 gradient_type value;
4007
4008 if (dof_pair[0] != numbers::invalid_unsigned_int)
4009 value +=
4010 0.5 *
4011 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
4012 q_point);
4013
4014 if (dof_pair[1] != numbers::invalid_unsigned_int)
4015 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
4016 .gradient(dof_pair[1], q_point);
4017
4018 return value;
4019 }
4020
4021
4022
4023 template <int dim, int spacedim>
4026 const unsigned int interface_dof_index,
4027 const unsigned int q_point) const
4028 {
4029 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4030
4031 if (this->fe_interface->at_boundary())
4032 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
4033 dof_pair[0], q_point);
4034
4035 gradient_type value;
4036
4037 if (dof_pair[0] != numbers::invalid_unsigned_int)
4038 value +=
4039 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
4040 q_point);
4041
4042 if (dof_pair[1] != numbers::invalid_unsigned_int)
4043 value -=
4044 (*(this->fe_interface->fe_face_values_neighbor))[extractor].gradient(
4045 dof_pair[1], q_point);
4046
4047 return value;
4048 }
4049
4050
4051
4052 template <int dim, int spacedim>
4054 Vector<dim, spacedim>::jump_gradient(const unsigned int interface_dof_index,
4055 const unsigned int q_point) const
4056 {
4057 return jump_in_gradients(interface_dof_index, q_point);
4058 }
4059
4060
4061
4062 template <int dim, int spacedim>
4065 const unsigned int interface_dof_index,
4066 const unsigned int q_point) const
4067 {
4068 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4069
4070 if (this->fe_interface->at_boundary())
4071 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
4072 dof_pair[0], q_point);
4073
4074 hessian_type value;
4075
4076 if (dof_pair[0] != numbers::invalid_unsigned_int)
4077 value +=
4078 0.5 *
4079 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
4080 q_point);
4081
4082 if (dof_pair[1] != numbers::invalid_unsigned_int)
4083 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
4084 .hessian(dof_pair[1], q_point);
4085
4086 return value;
4087 }
4088
4089
4090
4091 template <int dim, int spacedim>
4093 Vector<dim, spacedim>::average_hessian(const unsigned int interface_dof_index,
4094 const unsigned int q_point) const
4095 {
4096 return average_of_hessians(interface_dof_index, q_point);
4097 }
4098
4099
4100
4101 template <int dim, int spacedim>
4104 const unsigned int interface_dof_index,
4105 const unsigned int q_point) const
4106 {
4107 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4108
4109 if (this->fe_interface->at_boundary())
4110 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
4111 dof_pair[0], q_point);
4112
4113 hessian_type value;
4114
4115 if (dof_pair[0] != numbers::invalid_unsigned_int)
4116 value +=
4117 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
4118 q_point);
4119
4120 if (dof_pair[1] != numbers::invalid_unsigned_int)
4121 value -=
4122 (*(this->fe_interface->fe_face_values_neighbor))[extractor].hessian(
4123 dof_pair[1], q_point);
4124
4125 return value;
4126 }
4127
4128
4129
4130 template <int dim, int spacedim>
4132 Vector<dim, spacedim>::jump_hessian(const unsigned int interface_dof_index,
4133 const unsigned int q_point) const
4134 {
4135 return jump_in_hessians(interface_dof_index, q_point);
4136 }
4137
4138
4139
4140 template <int dim, int spacedim>
4143 const unsigned int interface_dof_index,
4144 const unsigned int q_point) const
4145 {
4146 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4147
4148 if (this->fe_interface->at_boundary())
4149 return (*(this->fe_interface->fe_face_values))[extractor]
4150 .third_derivative(dof_pair[0], q_point);
4151
4152 third_derivative_type value;
4153
4154 if (dof_pair[0] != numbers::invalid_unsigned_int)
4155 value +=
4156 (*(this->fe_interface->fe_face_values))[extractor].third_derivative(
4157 dof_pair[0], q_point);
4158
4159 if (dof_pair[1] != numbers::invalid_unsigned_int)
4160 value -= (*(this->fe_interface->fe_face_values_neighbor))[extractor]
4161 .third_derivative(dof_pair[1], q_point);
4162
4163 return value;
4164 }
4165
4166
4167
4168 template <int dim, int spacedim>
4169 template <class InputVector>
4170 void
4172 const bool here_or_there,
4173 const InputVector &local_dof_values,
4174 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4175 const
4176 {
4177 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
4178
4179 for (const auto dof_index : this->fe_interface->dof_indices())
4180 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4181 {
4182 if (dof_index == 0)
4183 values[q_index] = 0.;
4184
4185 values[q_index] += local_dof_values[dof_index] *
4186 value(here_or_there, dof_index, q_index);
4187 }
4188 }
4189
4190
4191
4192 template <int dim, int spacedim>
4193 template <class InputVector>
4194 void
4196 const bool here_or_there,
4197 const InputVector &fe_function,
4198 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4199 const
4200 {
4201 std::vector<typename InputVector::value_type> local_dof_values(
4202 this->fe_interface->n_current_interface_dofs());
4203 this->get_local_dof_values(fe_function, local_dof_values);
4204
4205 get_function_values_from_local_dof_values(here_or_there,
4206 local_dof_values,
4207 values);
4208 }
4209
4210
4211
4212 template <int dim, int spacedim>
4213 template <class InputVector>
4214 void
4216 const InputVector &local_dof_values,
4217 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4218 const
4219 {
4220 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
4221
4222 for (const auto dof_index : this->fe_interface->dof_indices())
4223 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4224 {
4225 if (dof_index == 0)
4226 values[q_index] = 0.;
4227
4228 values[q_index] +=
4229 local_dof_values[dof_index] * jump_in_values(dof_index, q_index);
4230 }
4231 }
4232
4233
4234
4235 template <int dim, int spacedim>
4236 template <class InputVector>
4237 void
4239 const InputVector &fe_function,
4240 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4241 const
4242 {
4243 std::vector<typename InputVector::value_type> local_dof_values(
4244 this->fe_interface->n_current_interface_dofs());
4245 this->get_local_dof_values(fe_function, local_dof_values);
4246
4247 get_jump_in_function_values_from_local_dof_values(local_dof_values, values);
4248 }
4249
4250
4251
4252 template <int dim, int spacedim>
4253 template <class InputVector>
4254 void
4256 const InputVector &local_dof_values,
4257 std::vector<solution_gradient_type<typename InputVector::value_type>>
4258 &gradients) const
4259 {
4260 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
4261
4262 for (const auto dof_index : this->fe_interface->dof_indices())
4263 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4264 {
4265 if (dof_index == 0)
4266 gradients[q_index] = 0.;
4267
4268 gradients[q_index] +=
4269 local_dof_values[dof_index] * jump_in_gradients(dof_index, q_index);
4270 }
4271 }
4272
4273
4274
4275 template <int dim, int spacedim>
4276 template <class InputVector>
4277 void
4279 const InputVector &fe_function,
4280 std::vector<solution_gradient_type<typename InputVector::value_type>>
4281 &gradients) const
4282 {
4283 std::vector<typename InputVector::value_type> local_dof_values(
4284 this->fe_interface->n_current_interface_dofs());
4285 this->get_local_dof_values(fe_function, local_dof_values);
4286
4287 get_jump_in_function_gradients_from_local_dof_values(local_dof_values,
4288 gradients);
4289 }
4290
4291
4292
4293 template <int dim, int spacedim>
4294 template <class InputVector>
4295 void
4297 const InputVector &local_dof_values,
4298 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4299 const
4300 {
4301 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
4302
4303 for (const auto dof_index : this->fe_interface->dof_indices())
4304 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4305 {
4306 if (dof_index == 0)
4307 values[q_index] = 0.;
4308
4309 values[q_index] +=
4310 local_dof_values[dof_index] * average_of_values(dof_index, q_index);
4311 }
4312 }
4313
4314
4315
4316 template <int dim, int spacedim>
4317 template <class InputVector>
4318 void
4320 const InputVector &fe_function,
4321 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4322 const
4323 {
4324 std::vector<typename InputVector::value_type> local_dof_values(
4325 this->fe_interface->n_current_interface_dofs());
4326 this->get_local_dof_values(fe_function, local_dof_values);
4327
4328 get_average_of_function_values_from_local_dof_values(local_dof_values,
4329 values);
4330 }
4331
4332
4333
4334 template <int dim, int spacedim>
4335 template <class InputVector>
4336 void
4339 const InputVector &local_dof_values,
4340 std::vector<solution_gradient_type<typename InputVector::value_type>>
4341 &gradients) const
4342 {
4343 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
4344
4345 for (const auto dof_index : this->fe_interface->dof_indices())
4346 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4347 {
4348 if (dof_index == 0)
4349 gradients[q_index] = 0.;
4350
4351 gradients[q_index] += local_dof_values[dof_index] *
4352 average_of_gradients(dof_index, q_index);
4353 }
4354 }
4355
4356
4357
4358 template <int dim, int spacedim>
4359 template <class InputVector>
4360 void
4362 const InputVector &fe_function,
4363 std::vector<solution_gradient_type<typename InputVector::value_type>>
4364 &gradients) const
4365 {
4366 std::vector<typename InputVector::value_type> local_dof_values(
4367 this->fe_interface->n_current_interface_dofs());
4368 this->get_local_dof_values(fe_function, local_dof_values);
4369
4370 get_average_of_function_gradients_from_local_dof_values(local_dof_values,
4371 gradients);
4372 }
4373
4374
4375
4376 template <int dim, int spacedim>
4377 template <class InputVector>
4378 void
4380 const InputVector &local_dof_values,
4381 std::vector<solution_hessian_type<typename InputVector::value_type>>
4382 &hessians) const
4383 {
4384 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
4385
4386 for (const auto dof_index : this->fe_interface->dof_indices())
4387 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4388 {
4389 if (dof_index == 0)
4390 hessians[q_index] = 0.;
4391
4392 hessians[q_index] +=
4393 local_dof_values[dof_index] * jump_in_hessians(dof_index, q_index);
4394 }
4395 }
4396
4397
4398
4399 template <int dim, int spacedim>
4400 template <class InputVector>
4401 void
4403 const InputVector &fe_function,
4404 std::vector<solution_hessian_type<typename InputVector::value_type>>
4405 &hessians) const
4406 {
4407 std::vector<typename InputVector::value_type> local_dof_values(
4408 this->fe_interface->n_current_interface_dofs());
4409 this->get_local_dof_values(fe_function, local_dof_values);
4410
4411 get_jump_in_function_hessians_from_local_dof_values(local_dof_values,
4412 hessians);
4413 }
4414
4415
4416
4417 template <int dim, int spacedim>
4418 template <class InputVector>
4419 void
4421 const InputVector &local_dof_values,
4422 std::vector<solution_hessian_type<typename InputVector::value_type>>
4423 &hessians) const
4424 {
4425 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
4426
4427 for (const unsigned int dof_index : this->fe_interface->dof_indices())
4428 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4429 {
4430 if (dof_index == 0)
4431 hessians[q_index] = 0.;
4432
4433 hessians[q_index] += local_dof_values[dof_index] *
4434 average_of_hessians(dof_index, q_index);
4435 }
4436 }
4437
4438
4439
4440 template <int dim, int spacedim>
4441 template <class InputVector>
4442 void
4444 const InputVector &fe_function,
4445 std::vector<solution_hessian_type<typename InputVector::value_type>>
4446 &hessians) const
4447 {
4448 std::vector<typename InputVector::value_type> local_dof_values(
4449 this->fe_interface->n_current_interface_dofs());
4450 this->get_local_dof_values(fe_function, local_dof_values);
4451
4452 get_average_of_function_hessians_from_local_dof_values(local_dof_values,
4453 hessians);
4454 }
4455
4456
4457
4458 template <int dim, int spacedim>
4459 template <class InputVector>
4460 void
4463 const InputVector &local_dof_values,
4464 std::vector<
4465 solution_third_derivative_type<typename InputVector::value_type>>
4466 &third_derivatives) const
4467 {
4468 AssertDimension(third_derivatives.size(),
4469 this->fe_interface->n_quadrature_points);
4470
4471 for (const unsigned int dof_index : this->fe_interface->dof_indices())
4472 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4473 {
4474 if (dof_index == 0)
4475 third_derivatives[q_index] = 0.;
4476
4477 third_derivatives[q_index] +=
4478 local_dof_values[dof_index] *
4479 jump_in_third_derivatives(dof_index, q_index);
4480 }
4481 }
4482
4483
4484
4485 template <int dim, int spacedim>
4486 template <class InputVector>
4487 void
4489 const InputVector &fe_function,
4490 std::vector<
4491 solution_third_derivative_type<typename InputVector::value_type>>
4492 &third_derivatives) const
4493 {
4494 std::vector<typename InputVector::value_type> local_dof_values(
4495 this->fe_interface->n_current_interface_dofs());
4496 this->get_local_dof_values(fe_function, local_dof_values);
4497
4498 get_jump_in_function_third_derivatives_from_local_dof_values(
4499 local_dof_values, third_derivatives);
4500 }
4501} // namespace FEInterfaceViews
4502
4503#endif // DOXYGEN
4504
4506
4507#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:318
Definition point.h:111
friend class Vector
Definition vector.h:1124
Number value_type
Definition vector.h:137
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:498
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:499
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)
Definition exceptions.h:489
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:294
std::size_t size
Definition mpi.cc:734
typename ::internal::FEInterfaceViews::ViewType< dim, spacedim, Extractor >::type View
Definition hp.h:117
const types::fe_index invalid_fe_index
Definition types.h:243
static const unsigned int invalid_unsigned_int
Definition types.h:220
const types::global_dof_index invalid_dof_index
Definition types.h:252
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