Reference documentation for deal.II version GIT acb749f4f5 2024-02-21 16:10:02+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// Copyright (C) 2018 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16#ifndef dealii_fe_interface_values_h
17#define dealii_fe_interface_values_h
18
19#include <deal.II/base/config.h>
20
22
24#include <deal.II/fe/mapping.h>
25
29
31
32#ifndef DOXYGEN
33template <int dim, int spacedim>
35#endif
36
42{
46 template <int dim, int spacedim = dim>
47 class Base
48 {
49 public:
54
55 protected:
60
65 template <class InputVector, class OutputVector>
66 void
67 get_local_dof_values(const InputVector &dof_values,
68 OutputVector &local_dof_values) const;
69 };
70
71
72
76 template <int dim, int spacedim = dim>
77 class Scalar : public Base<dim, spacedim>
78 {
79 public:
83 using value_type = double;
84
91
98
105
112 template <typename Number>
114
121 template <typename Number>
124
131 template <typename Number>
134
141 template <typename Number>
144
149 const unsigned int component);
150
176 value(const bool here_or_there,
177 const unsigned int interface_dof_index,
178 const unsigned int q_point) const;
179
198 jump_in_values(const unsigned int interface_dof_index,
199 const unsigned int q_point) const;
200
212 jump_in_gradients(const unsigned int interface_dof_index,
213 const unsigned int q_point) const;
214
226 jump_in_hessians(const unsigned int interface_dof_index,
227 const unsigned int q_point) const;
228
240 jump_in_third_derivatives(const unsigned int interface_dof_index,
241 const unsigned int q_point) const;
242
261 average_of_values(const unsigned int interface_dof_index,
262 const unsigned int q_point) const;
263
275 average_of_gradients(const unsigned int interface_dof_index,
276 const unsigned int q_point) const;
277
290 average_of_hessians(const unsigned int interface_dof_index,
291 const unsigned int q_point) const;
292
319 template <class InputVector>
320 void
322 const bool here_or_there,
323 const InputVector &fe_function,
325 &values) const;
326
349 template <class InputVector>
350 void
352 const bool here_or_there,
353 const InputVector &local_dof_values,
355 &values) const;
356
377 template <class InputVector>
378 void
380 const InputVector &fe_function,
382 &values) const;
383
391 template <class InputVector>
392 void
394 const InputVector &local_dof_values,
396 &values) const;
397
411 template <class InputVector>
412 void
414 const InputVector &fe_function,
416 &gradients) const;
417
425 template <class InputVector>
426 void
428 const InputVector &local_dof_values,
430 &gradients) const;
431
445 template <class InputVector>
446 void
448 const InputVector &fe_function,
450 &hessians) const;
451
459 template <class InputVector>
460 void
462 const InputVector &local_dof_values,
464 &hessians) const;
465
480 template <class InputVector>
481 void
483 const InputVector &fe_function,
484 std::vector<
486 &third_derivatives) const;
487
495 template <class InputVector>
496 void
498 const InputVector &local_dof_values,
499 std::vector<
501 &third_derivatives) const;
502
523 template <class InputVector>
524 void
526 const InputVector &fe_function,
528 &values) const;
529
537 template <class InputVector>
538 void
540 const InputVector &local_dof_values,
542 &values) const;
543
557 template <class InputVector>
558 void
560 const InputVector &fe_function,
562 &gradients) const;
563
571 template <class InputVector>
572 void
574 const InputVector &local_dof_values,
576 &gradients) const;
577
591 template <class InputVector>
592 void
594 const InputVector &fe_function,
596 &hessians) const;
597
605 template <class InputVector>
606 void
608 const InputVector &local_dof_values,
610 &hessians) const;
611
614 private:
619 };
620
621
622
626 template <int dim, int spacedim = dim>
627 class Vector : public Base<dim, spacedim>
628 {
629 public:
635
642
650
658
665 template <typename Number>
667
674 template <typename Number>
677
684 template <typename Number>
687
694 template <typename Number>
697
702 const unsigned int first_vector_component);
703
729 value(const bool here_or_there,
730 const unsigned int interface_dof_index,
731 const unsigned int q_point) const;
732
750 jump_in_values(const unsigned int interface_dof_index,
751 const unsigned int q_point) const;
752
763 jump_in_gradients(const unsigned int interface_dof_index,
764 const unsigned int q_point) const;
765
772 jump_gradient(const unsigned int interface_dof_index,
773 const unsigned int q_point) const;
774
786 jump_in_hessians(const unsigned int interface_dof_index,
787 const unsigned int q_point) const;
788
795 jump_hessian(const unsigned int interface_dof_index,
796 const unsigned int q_point) const;
797
809 jump_in_third_derivatives(const unsigned int interface_dof_index,
810 const unsigned int q_point) const;
811
829 average_of_values(const unsigned int interface_dof_index,
830 const unsigned int q_point) const;
831
842 average_of_gradients(const unsigned int interface_dof_index,
843 const unsigned int q_point) const;
844
857 average_of_hessians(const unsigned int interface_dof_index,
858 const unsigned int q_point) const;
859
866 average_hessian(const unsigned int interface_dof_index,
867 const unsigned int q_point) const;
868
895 template <class InputVector>
896 void
898 const bool here_or_there,
899 const InputVector &fe_function,
901 &values) const;
902
925 template <class InputVector>
926 void
928 const bool here_or_there,
929 const InputVector &local_dof_values,
931 &values) const;
932
953 template <class InputVector>
954 void
956 const InputVector &fe_function,
958 &values) const;
959
967 template <class InputVector>
968 void
970 const InputVector &local_dof_values,
972 &values) const;
973
987 template <class InputVector>
988 void
990 const InputVector &fe_function,
992 &gradients) const;
993
1001 template <class InputVector>
1002 void
1004 const InputVector &local_dof_values,
1006 &gradients) const;
1007
1021 template <class InputVector>
1022 void
1024 const InputVector &fe_function,
1026 &hessians) const;
1027
1035 template <class InputVector>
1036 void
1038 const InputVector &local_dof_values,
1040 &hessians) const;
1041
1056 template <class InputVector>
1057 void
1059 const InputVector &fe_function,
1060 std::vector<
1062 &third_derivatives) const;
1063
1071 template <class InputVector>
1072 void
1074 const InputVector &local_dof_values,
1075 std::vector<
1077 &third_derivatives) const;
1078
1099 template <class InputVector>
1100 void
1102 const InputVector &fe_function,
1104 &values) const;
1105
1113 template <class InputVector>
1114 void
1116 const InputVector &local_dof_values,
1118 &values) const;
1119
1133 template <class InputVector>
1134 void
1136 const InputVector &fe_function,
1138 &gradients) const;
1139
1147 template <class InputVector>
1148 void
1150 const InputVector &local_dof_values,
1152 &gradients) const;
1153
1167 template <class InputVector>
1168 void
1170 const InputVector &fe_function,
1172 &hessians) const;
1173
1181 template <class InputVector>
1182 void
1184 const InputVector &local_dof_values,
1186 &hessians) const;
1187
1190 private:
1195 };
1196} // namespace FEInterfaceViews
1197
1198
1199namespace internal
1200{
1202 {
1207 template <int dim, int spacedim, typename Extractor>
1209 {};
1210
1218 template <int dim, int spacedim>
1219 struct ViewType<dim, spacedim, FEValuesExtractors::Scalar>
1220 {
1221 using type = typename ::FEInterfaceViews::Scalar<dim, spacedim>;
1222 };
1223
1231 template <int dim, int spacedim>
1232 struct ViewType<dim, spacedim, FEValuesExtractors::Vector>
1233 {
1234 using type = typename ::FEInterfaceViews::Vector<dim, spacedim>;
1235 };
1236 } // namespace FEInterfaceViews
1237} // namespace internal
1238
1239namespace FEInterfaceViews
1240{
1245 template <int dim, int spacedim, typename Extractor>
1246 using View = typename ::internal::FEInterfaceViews::
1247 ViewType<dim, spacedim, Extractor>::type;
1248} // namespace FEInterfaceViews
1249
1250
1251
1276template <int dim, int spacedim = dim>
1278{
1279public:
1283 const unsigned int n_quadrature_points;
1284
1292 const Quadrature<dim - 1> &quadrature,
1293 const UpdateFlags update_flags);
1294
1302 const hp::QCollection<dim - 1> &quadrature,
1303 const UpdateFlags update_flags);
1304
1312 const Quadrature<dim - 1> &quadrature,
1313 const UpdateFlags update_flags);
1314
1320 const hp::MappingCollection<dim, spacedim> &mapping_collection,
1321 const hp::FECollection<dim, spacedim> &fe_collection,
1322 const hp::QCollection<dim - 1> &quadrature_collection,
1323 const UpdateFlags update_flags);
1324
1329 const hp::QCollection<dim - 1> &quadrature_collection,
1330 const UpdateFlags update_flags);
1331
1446 template <typename CellIteratorType, typename CellNeighborIteratorType>
1447 void
1448 reinit(const CellIteratorType &cell,
1449 const unsigned int face_no,
1450 const unsigned int sub_face_no,
1451 const CellNeighborIteratorType &cell_neighbor,
1452 const unsigned int face_no_neighbor,
1453 const unsigned int sub_face_no_neighbor,
1454 const unsigned int q_index = numbers::invalid_unsigned_int,
1455 const unsigned int mapping_index = numbers::invalid_unsigned_int,
1456 const unsigned int fe_index = numbers::invalid_unsigned_int);
1457
1484 template <typename CellIteratorType>
1485 void
1486 reinit(const CellIteratorType &cell,
1487 const unsigned int face_no,
1488 const unsigned int q_index = numbers::invalid_unsigned_int,
1489 const unsigned int mapping_index = numbers::invalid_unsigned_int,
1490 const unsigned int fe_index = numbers::invalid_unsigned_int);
1491
1500 get_fe_face_values(const unsigned int cell_index) const;
1501
1507
1512 get_fe() const;
1513
1517 const Quadrature<dim - 1> &
1519
1525
1531
1535 const hp::QCollection<dim - 1> &
1537
1542 bool
1544
1550
1558 get_cell(const unsigned int cell_index) const;
1559
1567 unsigned int
1568 get_face_number(const unsigned int cell_index) const;
1569
1581 bool
1583
1595 double
1596 JxW(const unsigned int quadrature_point) const;
1597
1603 const std::vector<double> &
1605
1615 const std::vector<Tensor<1, spacedim>> &
1617
1627
1633 const std::vector<Point<spacedim>> &
1635
1646 unsigned
1648
1678
1689 std::vector<types::global_dof_index>
1691
1704 std::array<unsigned int, 2>
1705 interface_dof_to_dof_indices(const unsigned int interface_dof_index) const;
1706
1716 normal(const unsigned int q_point_index) const;
1717
1748 double
1749 shape_value(const bool here_or_there,
1750 const unsigned int interface_dof_index,
1751 const unsigned int q_point,
1752 const unsigned int component = 0) const;
1753
1781 double
1782 jump_in_shape_values(const unsigned int interface_dof_index,
1783 const unsigned int q_point,
1784 const unsigned int component = 0) const;
1785
1800 jump_in_shape_gradients(const unsigned int interface_dof_index,
1801 const unsigned int q_point,
1802 const unsigned int component = 0) const;
1803
1819 jump_in_shape_hessians(const unsigned int interface_dof_index,
1820 const unsigned int q_point,
1821 const unsigned int component = 0) const;
1822
1837 jump_in_shape_3rd_derivatives(const unsigned int interface_dof_index,
1838 const unsigned int q_point,
1839 const unsigned int component = 0) const;
1840
1863 double
1864 average_of_shape_values(const unsigned int interface_dof_index,
1865 const unsigned int q_point,
1866 const unsigned int component = 0) const;
1867
1882 average_of_shape_gradients(const unsigned int interface_dof_index,
1883 const unsigned int q_point,
1884 const unsigned int component = 0) const;
1885
1901 average_of_shape_hessians(const unsigned int interface_dof_index,
1902 const unsigned int q_point,
1903 const unsigned int component = 0) const;
1904
1924 template <class InputVector>
1925 void
1927 const InputVector &fe_function,
1928 std::vector<typename InputVector::value_type> &values) const;
1929
1938 template <class InputVector>
1939 void
1941 const InputVector &fe_function,
1943 &gradients) const;
1944
1952 template <class InputVector>
1953 void
1955 const InputVector &fe_function,
1957 &hessians) const;
1958
1967 template <class InputVector>
1968 void
1970 const InputVector &fe_function,
1972 &third_derivatives) const;
1973
1989 template <class InputVector>
1990 void
1992 const InputVector &fe_function,
1993 std::vector<typename InputVector::value_type> &values) const;
1994
2002 template <class InputVector>
2003 void
2005 const InputVector &fe_function,
2007 &gradients) const;
2008
2016 template <class InputVector>
2017 void
2019 const InputVector &fe_function,
2021 &hessians) const;
2022
2042
2051
2056private:
2060 std::vector<types::global_dof_index> interface_dof_indices;
2061
2067 std::vector<std::array<unsigned int, 2>> dofmap;
2068
2074
2081
// non-hp data
2086
2090 std::unique_ptr<FEFaceValues<dim, spacedim>> internal_fe_face_values;
2091
2095 std::unique_ptr<FESubfaceValues<dim, spacedim>> internal_fe_subface_values;
2096
2100 std::unique_ptr<FEFaceValues<dim, spacedim>> internal_fe_face_values_neighbor;
2101
2105 std::unique_ptr<FESubfaceValues<dim, spacedim>>
2107
// non-hp data
2109
// hp data
2114
2119 std::unique_ptr<hp::FEFaceValues<dim, spacedim>> internal_hp_fe_face_values;
2120
2125 std::unique_ptr<hp::FESubfaceValues<dim, spacedim>>
2127
2132 std::unique_ptr<hp::FEFaceValues<dim, spacedim>>
2134
2139 std::unique_ptr<hp::FESubfaceValues<dim, spacedim>>
2141
2149 "The current function doesn't make sense when used with a "
2150 "FEInterfaceValues object with hp-capabilities.");
2151
2159 "The current function doesn't make sense when used with a "
2160 "FEInterfaceValues object without hp-capabilities.");
2161
// hp data
2163
2164 /*
2165 * Make the view classes friends of this class, since they
2166 * access internal data.
2167 */
2168 template <int, int>
2170 template <int, int>
2172};
2173
2174
2175
2176#ifndef DOXYGEN
2177
2178/*---------------------- Inline functions ---------------------*/
2179
2180template <int dim, int spacedim>
2182 const Mapping<dim, spacedim> &mapping,
2184 const Quadrature<dim - 1> &quadrature,
2185 const UpdateFlags update_flags)
2186 : n_quadrature_points(quadrature.size())
2187 , fe_face_values(nullptr)
2188 , fe_face_values_neighbor(nullptr)
2189 , internal_fe_face_values(
2190 std::make_unique<FEFaceValues<dim, spacedim>>(mapping,
2191 fe,
2192 quadrature,
2193 update_flags))
2194 , internal_fe_subface_values(
2195 std::make_unique<FESubfaceValues<dim, spacedim>>(mapping,
2196 fe,
2197 quadrature,
2198 update_flags))
2199 , internal_fe_face_values_neighbor(
2200 std::make_unique<FEFaceValues<dim, spacedim>>(mapping,
2201 fe,
2202 quadrature,
2203 update_flags))
2204 , internal_fe_subface_values_neighbor(
2205 std::make_unique<FESubfaceValues<dim, spacedim>>(mapping,
2206 fe,
2207 quadrature,
2208 update_flags))
2209{}
2210
2211
2212
2213template <int dim, int spacedim>
2216 const Quadrature<dim - 1> &quadrature,
2217 const UpdateFlags update_flags)
2219 fe.reference_cell().template get_default_linear_mapping<dim, spacedim>(),
2220 fe,
2221 quadrature,
2222 update_flags)
2223{}
2224
2225
2226
2227template <int dim, int spacedim>
2229 const Mapping<dim, spacedim> &mapping,
2231 const hp::QCollection<dim - 1> &quadrature,
2232 const UpdateFlags update_flags)
2233 : n_quadrature_points(quadrature.max_n_quadrature_points())
2234 , fe_face_values(nullptr)
2235 , fe_face_values_neighbor(nullptr)
2236 , internal_fe_face_values(
2237 std::make_unique<FEFaceValues<dim, spacedim>>(mapping,
2238 fe,
2239 quadrature,
2240 update_flags))
2241 , internal_fe_subface_values(
2242 std::make_unique<FESubfaceValues<dim, spacedim>>(mapping,
2243 fe,
2244 quadrature,
2245 update_flags))
2246 , internal_fe_face_values_neighbor(
2247 std::make_unique<FEFaceValues<dim, spacedim>>(mapping,
2248 fe,
2249 quadrature[0],
2250 update_flags))
2251 , internal_fe_subface_values_neighbor(
2252 std::make_unique<FESubfaceValues<dim, spacedim>>(mapping,
2253 fe,
2254 quadrature[0],
2255 update_flags))
2256{}
2257
2258
2259
2260template <int dim, int spacedim>
2262 const hp::MappingCollection<dim, spacedim> &mapping_collection,
2263 const hp::FECollection<dim, spacedim> &fe_collection,
2264 const hp::QCollection<dim - 1> &quadrature_collection,
2265 const UpdateFlags update_flags)
2266 : n_quadrature_points(quadrature_collection.max_n_quadrature_points())
2267 , fe_face_values(nullptr)
2268 , fe_face_values_neighbor(nullptr)
2269 , internal_hp_fe_face_values(
2270 std::make_unique<hp::FEFaceValues<dim, spacedim>>(mapping_collection,
2271 fe_collection,
2272 quadrature_collection,
2273 update_flags))
2274 , internal_hp_fe_subface_values(
2275 std::make_unique<hp::FESubfaceValues<dim, spacedim>>(
2276 mapping_collection,
2277 fe_collection,
2278 quadrature_collection,
2279 update_flags))
2280 , internal_hp_fe_face_values_neighbor(
2281 std::make_unique<hp::FEFaceValues<dim, spacedim>>(mapping_collection,
2282 fe_collection,
2283 quadrature_collection,
2284 update_flags))
2285 , internal_hp_fe_subface_values_neighbor(
2286 std::make_unique<hp::FESubfaceValues<dim, spacedim>>(
2287 mapping_collection,
2288 fe_collection,
2289 quadrature_collection,
2290 update_flags))
2291{
2292 AssertDimension(dim, spacedim);
2293}
2294
2295
2296
2297template <int dim, int spacedim>
2299 const hp::FECollection<dim, spacedim> &fe_collection,
2300 const hp::QCollection<dim - 1> &quadrature_collection,
2301 const UpdateFlags update_flags)
2302 : FEInterfaceValues(fe_collection.get_reference_cell_default_linear_mapping(),
2303 fe_collection,
2304 quadrature_collection,
2305 update_flags)
2306{}
2307
2308
2309
2310template <int dim, int spacedim>
2311template <typename CellIteratorType, typename CellNeighborIteratorType>
2312void
2314 const CellIteratorType &cell,
2315 const unsigned int face_no,
2316 const unsigned int sub_face_no,
2317 const CellNeighborIteratorType &cell_neighbor,
2318 const unsigned int face_no_neighbor,
2319 const unsigned int sub_face_no_neighbor,
2320 const unsigned int q_index,
2321 const unsigned int mapping_index,
2322 const unsigned int fe_index)
2323{
2324 Assert(internal_fe_face_values || internal_hp_fe_face_values,
2326
2327 if (internal_fe_face_values)
2328 {
2329 if (sub_face_no == numbers::invalid_unsigned_int)
2330 {
2331 internal_fe_face_values->reinit(cell, face_no);
2332 fe_face_values = internal_fe_face_values.get();
2333 }
2334 else
2335 {
2336 internal_fe_subface_values->reinit(cell, face_no, sub_face_no);
2337 fe_face_values = internal_fe_subface_values.get();
2338 }
2339 if (sub_face_no_neighbor == numbers::invalid_unsigned_int)
2340 {
2341 internal_fe_face_values_neighbor->reinit(cell_neighbor,
2342 face_no_neighbor);
2343 fe_face_values_neighbor = internal_fe_face_values_neighbor.get();
2344 }
2345 else
2346 {
2347 internal_fe_subface_values_neighbor->reinit(cell_neighbor,
2348 face_no_neighbor,
2349 sub_face_no_neighbor);
2350 fe_face_values_neighbor = internal_fe_subface_values_neighbor.get();
2351 }
2352
2353 AssertDimension(fe_face_values->n_quadrature_points,
2354 fe_face_values_neighbor->n_quadrature_points);
2355
2356 const_cast<unsigned int &>(this->n_quadrature_points) =
2357 fe_face_values->n_quadrature_points;
2358 }
2359 else if (internal_hp_fe_face_values)
2360 {
2361 unsigned int used_q_index = q_index;
2362 unsigned int used_mapping_index = mapping_index;
2363
2364 // First check. If there is only one element in a collection, and if none
2365 // had been specified explicitly, then that's clearly the one to take:
2366 if (used_q_index == numbers::invalid_unsigned_int)
2367 if (internal_hp_fe_face_values->get_quadrature_collection().size() == 1)
2368 used_q_index = 0;
2369
2370 if (used_mapping_index == numbers::invalid_unsigned_int)
2371 if (internal_hp_fe_face_values->get_mapping_collection().size() == 1)
2372 used_mapping_index = 0;
2373
2374 // Second check: See if the two quadrature objects are the same, because
2375 // in that case it does not matter which one we use. Unfortunately, we
2376 // currently have no way of testing that two mapping objects are the
2377 // same :-(
2378 if (used_q_index == numbers::invalid_unsigned_int)
2379 if (internal_hp_fe_face_values
2380 ->get_quadrature_collection()[cell->active_fe_index()] ==
2381 internal_hp_fe_face_values
2382 ->get_quadrature_collection()[cell_neighbor->active_fe_index()])
2383 used_q_index = cell->active_fe_index();
2384
2385 // Third check, if the above did not already suffice. We see if we
2386 // can get somewhere via the dominated's finite element index.
2387 const unsigned int dominated_fe_index =
2388 ((used_q_index == numbers::invalid_unsigned_int) ||
2389 (used_mapping_index == numbers::invalid_unsigned_int) ?
2390 internal_hp_fe_face_values->get_fe_collection().find_dominated_fe(
2391 {cell->active_fe_index(), cell_neighbor->active_fe_index()}) :
2392 numbers::invalid_unsigned_int);
2393
2394 if (used_q_index == numbers::invalid_unsigned_int)
2395 {
2396 Assert(dominated_fe_index != numbers::invalid_fe_index,
2397 ExcMessage(
2398 "You called this function with 'q_index' left at its "
2399 "default value, but this can only work if one of "
2400 "the two finite elements adjacent to this face "
2401 "dominates the other. See the documentation "
2402 "of this function for more information of how "
2403 "to deal with this situation."));
2404 used_q_index = dominated_fe_index;
2405 }
2406
2407 if (used_mapping_index == numbers::invalid_unsigned_int)
2408 {
2409 Assert(dominated_fe_index != numbers::invalid_fe_index,
2410 ExcMessage(
2411 "You called this function with 'mapping_index' left "
2412 "at its default value, but this can only work if one "
2413 "of the two finite elements adjacent to this face "
2414 "dominates the other. See the documentation "
2415 "of this function for more information of how "
2416 "to deal with this situation."));
2417 used_mapping_index = dominated_fe_index;
2418 }
2419
2420 // Same as if above, but when hp is enabled.
2421 if (sub_face_no == numbers::invalid_unsigned_int)
2422 {
2423 internal_hp_fe_face_values->reinit(
2424 cell, face_no, used_q_index, used_mapping_index, fe_index);
2425 fe_face_values = &const_cast<FEFaceValues<dim, spacedim> &>(
2426 internal_hp_fe_face_values->get_present_fe_values());
2427 }
2428 else
2429 {
2430 internal_hp_fe_subface_values->reinit(
2431 cell, face_no, sub_face_no, used_q_index, used_mapping_index);
2432
2433 fe_face_values = &const_cast<FESubfaceValues<dim, spacedim> &>(
2434 internal_hp_fe_subface_values->get_present_fe_values());
2435 }
2436 if (sub_face_no_neighbor == numbers::invalid_unsigned_int)
2437 {
2438 internal_hp_fe_face_values_neighbor->reinit(cell_neighbor,
2439 face_no_neighbor,
2440 used_q_index,
2441 used_mapping_index);
2442
2443 fe_face_values_neighbor = &const_cast<FEFaceValues<dim, spacedim> &>(
2444 internal_hp_fe_face_values_neighbor->get_present_fe_values());
2445 }
2446 else
2447 {
2448 internal_hp_fe_subface_values_neighbor->reinit(cell_neighbor,
2449 face_no_neighbor,
2450 sub_face_no_neighbor,
2451 used_q_index,
2452 used_mapping_index);
2453
2454 fe_face_values_neighbor =
2455 &const_cast<FESubfaceValues<dim, spacedim> &>(
2456 internal_hp_fe_subface_values_neighbor->get_present_fe_values());
2457 }
2458
2459 AssertDimension(fe_face_values->n_quadrature_points,
2460 fe_face_values_neighbor->n_quadrature_points);
2461
2462 const_cast<unsigned int &>(this->n_quadrature_points) =
2463 fe_face_values->n_quadrature_points;
2464 }
2465
2466 // Set up dof mapping and remove duplicates (for continuous elements).
2467 {
2468 // Get dof indices first:
2469 std::vector<types::global_dof_index> v(
2470 fe_face_values->get_fe().n_dofs_per_cell());
2471 cell->get_active_or_mg_dof_indices(v);
2472 std::vector<types::global_dof_index> v2(
2473 fe_face_values_neighbor->get_fe().n_dofs_per_cell());
2474 cell_neighbor->get_active_or_mg_dof_indices(v2);
2475
2476 // Fill a map from the global dof index to the left and right
2477 // local index.
2478 std::map<types::global_dof_index, std::pair<unsigned int, unsigned int>>
2479 tempmap;
2480 std::pair<unsigned int, unsigned int> invalid_entry(
2482
2483 for (unsigned int i = 0; i < v.size(); ++i)
2484 {
2485 // If not already existing, add an invalid entry:
2486 auto result = tempmap.insert(std::make_pair(v[i], invalid_entry));
2487 result.first->second.first = i;
2488 }
2489
2490 for (unsigned int i = 0; i < v2.size(); ++i)
2491 {
2492 // If not already existing, add an invalid entry:
2493 auto result = tempmap.insert(std::make_pair(v2[i], invalid_entry));
2494 result.first->second.second = i;
2495 }
2496
2497 // Transfer from the map to the sorted std::vectors.
2498 dofmap.resize(tempmap.size());
2499 interface_dof_indices.resize(tempmap.size());
2500 unsigned int idx = 0;
2501 for (auto &x : tempmap)
2502 {
2503 interface_dof_indices[idx] = x.first;
2504 dofmap[idx] = {{x.second.first, x.second.second}};
2505 ++idx;
2506 }
2507 }
2508}
2509
2510
2511
2512template <int dim, int spacedim>
2513template <typename CellIteratorType>
2514void
2515FEInterfaceValues<dim, spacedim>::reinit(const CellIteratorType &cell,
2516 const unsigned int face_no,
2517 const unsigned int q_index,
2518 const unsigned int mapping_index,
2519 const unsigned int fe_index)
2520{
2521 Assert(internal_fe_face_values || internal_hp_fe_face_values,
2523
2524 if (internal_fe_face_values)
2525 {
2526 internal_fe_face_values->reinit(cell, face_no);
2527 fe_face_values = internal_fe_face_values.get();
2528 fe_face_values_neighbor = nullptr;
2529
2530 interface_dof_indices.resize(fe_face_values->get_fe().n_dofs_per_cell());
2531 cell->get_active_or_mg_dof_indices(interface_dof_indices);
2532 }
2533 else if (internal_hp_fe_face_values)
2534 {
2535 internal_hp_fe_face_values->reinit(
2536 cell, face_no, q_index, mapping_index, fe_index);
2537 fe_face_values = &const_cast<FEFaceValues<dim> &>(
2538 internal_hp_fe_face_values->get_present_fe_values());
2539 fe_face_values_neighbor = nullptr;
2540
2541 interface_dof_indices.resize(fe_face_values->get_fe().n_dofs_per_cell());
2542 cell->get_active_or_mg_dof_indices(interface_dof_indices);
2543 }
2544
2545 dofmap.resize(interface_dof_indices.size());
2546 for (unsigned int i = 0; i < interface_dof_indices.size(); ++i)
2547 {
2548 dofmap[i] = {{i, numbers::invalid_unsigned_int}};
2549 }
2550}
2551
2552
2553
2554template <int dim, int spacedim>
2555inline double
2556FEInterfaceValues<dim, spacedim>::JxW(const unsigned int q) const
2557{
2558 Assert(fe_face_values != nullptr,
2559 ExcMessage("This call requires a call to reinit() first."));
2560 return fe_face_values->JxW(q);
2561}
2562
2563
2564
2565template <int dim, int spacedim>
2566const std::vector<double> &
2568{
2569 Assert(fe_face_values != nullptr,
2570 ExcMessage("This call requires a call to reinit() first."));
2571 return fe_face_values->get_JxW_values();
2572}
2573
2574
2575
2576template <int dim, int spacedim>
2577const std::vector<Tensor<1, spacedim>> &
2579{
2580 Assert(fe_face_values != nullptr,
2581 ExcMessage("This call requires a call to reinit() first."));
2582 return fe_face_values->get_normal_vectors();
2583}
2584
2585
2586
2587template <int dim, int spacedim>
2590{
2591 Assert(!has_hp_capabilities(), ExcOnlyAvailableWithoutHP());
2592 return internal_fe_face_values->get_mapping();
2593}
2594
2595
2596
2597template <int dim, int spacedim>
2600{
2601 Assert(!has_hp_capabilities(), ExcOnlyAvailableWithoutHP());
2602 return internal_fe_face_values->get_fe();
2603}
2604
2605
2606
2607template <int dim, int spacedim>
2608const Quadrature<dim - 1> &
2610{
2611 Assert(!has_hp_capabilities(), ExcOnlyAvailableWithoutHP());
2612 return internal_fe_face_values->get_quadrature();
2613}
2614
2615
2616
2617template <int dim, int spacedim>
2620{
2621 Assert(has_hp_capabilities(), ExcOnlyAvailableWithHP());
2622 return internal_hp_fe_face_values->get_mapping_collection();
2623}
2624
2625
2626
2627template <int dim, int spacedim>
2630{
2631 Assert(has_hp_capabilities(), ExcOnlyAvailableWithHP());
2632 return internal_hp_fe_face_values->get_fe_collection();
2633}
2634
2635
2636
2637template <int dim, int spacedim>
2638const hp::QCollection<dim - 1> &
2640{
2641 Assert(has_hp_capabilities(), ExcOnlyAvailableWithHP());
2642 return internal_hp_fe_face_values->get_quadrature_collection();
2643}
2644
2645
2646
2647template <int dim, int spacedim>
2648bool
2650{
2651 if (internal_hp_fe_face_values || internal_hp_fe_subface_values ||
2652 internal_hp_fe_face_values_neighbor ||
2653 internal_hp_fe_subface_values_neighbor)
2654 {
2655 Assert(!internal_fe_face_values, ExcInternalError());
2656 Assert(!internal_fe_subface_values, ExcInternalError());
2657 Assert(!internal_fe_face_values_neighbor, ExcInternalError());
2658 Assert(!internal_fe_subface_values_neighbor, ExcInternalError());
2659
2660 return true;
2661 }
2662
2663 Assert(internal_fe_face_values || internal_fe_subface_values ||
2664 internal_fe_face_values_neighbor ||
2665 internal_fe_subface_values_neighbor,
2667 Assert(!internal_hp_fe_face_values, ExcInternalError());
2668 Assert(!internal_hp_fe_subface_values, ExcInternalError());
2669 Assert(!internal_hp_fe_face_values_neighbor, ExcInternalError());
2670 Assert(!internal_hp_fe_subface_values_neighbor, ExcInternalError());
2671
2672 return false;
2673}
2674
2675
2676
2677template <int dim, int spacedim>
2680{
2681 return {0U, n_quadrature_points};
2682}
2683
2684
2685
2686template <int dim, int spacedim>
2687const std::vector<Point<spacedim>> &
2689{
2690 Assert(fe_face_values != nullptr,
2691 ExcMessage("This call requires a call to reinit() first."));
2692 return fe_face_values->get_quadrature_points();
2693}
2694
2695
2696
2697template <int dim, int spacedim>
2700{
2701 if (has_hp_capabilities())
2702 return internal_hp_fe_face_values->get_update_flags();
2703 else
2704 return internal_fe_face_values->get_update_flags();
2705}
2706
2707
2708
2709template <int dim, int spacedim>
2712{
2713 return get_fe_face_values(cell_index).get_cell();
2714}
2715
2716
2717
2718template <int dim, int spacedim>
2719inline unsigned int
2721 const unsigned int cell_index) const
2722{
2723 return get_fe_face_values(cell_index).get_face_number();
2724}
2725
2726
2727
2728template <int dim, int spacedim>
2729unsigned
2731{
2732 Assert(
2733 interface_dof_indices.size() > 0,
2734 ExcMessage(
2735 "n_current_interface_dofs() is only available after a call to reinit()."));
2736 return interface_dof_indices.size();
2737}
2738
2739
2740
2741template <int dim, int spacedim>
2744{
2745 return {0U, n_current_interface_dofs()};
2746}
2747
2748
2749
2750template <int dim, int spacedim>
2751bool
2753{
2754 return fe_face_values_neighbor == nullptr;
2755}
2756
2757
2758
2759template <int dim, int spacedim>
2760std::vector<types::global_dof_index>
2762{
2763 return interface_dof_indices;
2764}
2765
2766
2767
2768template <int dim, int spacedim>
2769std::array<unsigned int, 2>
2771 const unsigned int interface_dof_index) const
2772{
2773 AssertIndexRange(interface_dof_index, n_current_interface_dofs());
2774 return dofmap[interface_dof_index];
2775}
2776
2777
2778
2779template <int dim, int spacedim>
2782 const unsigned int cell_index) const
2783{
2785 Assert(
2786 cell_index == 0 || !at_boundary(),
2787 ExcMessage(
2788 "You are on a boundary, so you can only ask for the first FEFaceValues object."));
2789
2790 return (cell_index == 0) ? *fe_face_values : *fe_face_values_neighbor;
2791}
2792
2793
2794
2795template <int dim, int spacedim>
2797FEInterfaceValues<dim, spacedim>::normal(const unsigned int q_point_index) const
2798{
2799 return fe_face_values->normal_vector(q_point_index);
2800}
2801
2802
2803
2804template <int dim, int spacedim>
2805double
2807 const bool here_or_there,
2808 const unsigned int interface_dof_index,
2809 const unsigned int q_point,
2810 const unsigned int component) const
2811{
2812 const auto dof_pair = dofmap[interface_dof_index];
2813
2814 if (here_or_there && dof_pair[0] != numbers::invalid_unsigned_int)
2815 return get_fe_face_values(0).shape_value_component(dof_pair[0],
2816 q_point,
2817 component);
2818 if (!here_or_there && dof_pair[1] != numbers::invalid_unsigned_int)
2819 return get_fe_face_values(1).shape_value_component(dof_pair[1],
2820 q_point,
2821 component);
2822
2823 return 0.0;
2824}
2825
2826
2827
2828template <int dim, int spacedim>
2829double
2831 const unsigned int interface_dof_index,
2832 const unsigned int q_point,
2833 const unsigned int component) const
2834{
2835 const auto dof_pair = dofmap[interface_dof_index];
2836
2837 double value = 0.0;
2838
2839 if (dof_pair[0] != numbers::invalid_unsigned_int)
2840 value += get_fe_face_values(0).shape_value_component(dof_pair[0],
2841 q_point,
2842 component);
2843 if (dof_pair[1] != numbers::invalid_unsigned_int)
2844 value -= get_fe_face_values(1).shape_value_component(dof_pair[1],
2845 q_point,
2846 component);
2847 return value;
2848}
2849
2850
2851
2852template <int dim, int spacedim>
2853double
2855 const unsigned int interface_dof_index,
2856 const unsigned int q_point,
2857 const unsigned int component) const
2858{
2859 const auto dof_pair = dofmap[interface_dof_index];
2860
2861 if (at_boundary())
2862 return get_fe_face_values(0).shape_value_component(dof_pair[0],
2863 q_point,
2864 component);
2865
2866 double value = 0.0;
2867
2868 if (dof_pair[0] != numbers::invalid_unsigned_int)
2869 value += 0.5 * get_fe_face_values(0).shape_value_component(dof_pair[0],
2870 q_point,
2871 component);
2872 if (dof_pair[1] != numbers::invalid_unsigned_int)
2873 value += 0.5 * get_fe_face_values(1).shape_value_component(dof_pair[1],
2874 q_point,
2875 component);
2876
2877 return value;
2878}
2879
2880
2881
2882template <int dim, int spacedim>
2885 const unsigned int interface_dof_index,
2886 const unsigned int q_point,
2887 const unsigned int component) const
2888{
2889 const auto dof_pair = dofmap[interface_dof_index];
2890
2891 if (at_boundary())
2892 return get_fe_face_values(0).shape_grad_component(dof_pair[0],
2893 q_point,
2894 component);
2895
2897
2898 if (dof_pair[0] != numbers::invalid_unsigned_int)
2899 value += 0.5 * get_fe_face_values(0).shape_grad_component(dof_pair[0],
2900 q_point,
2901 component);
2902 if (dof_pair[1] != numbers::invalid_unsigned_int)
2903 value += 0.5 * get_fe_face_values(1).shape_grad_component(dof_pair[1],
2904 q_point,
2905 component);
2906
2907 return value;
2908}
2909
2910
2911
2912template <int dim, int spacedim>
2915 const unsigned int interface_dof_index,
2916 const unsigned int q_point,
2917 const unsigned int component) const
2918{
2919 const auto dof_pair = dofmap[interface_dof_index];
2920
2921 if (at_boundary())
2922 return get_fe_face_values(0).shape_hessian_component(dof_pair[0],
2923 q_point,
2924 component);
2925
2927
2928 if (dof_pair[0] != numbers::invalid_unsigned_int)
2929 value += 0.5 * get_fe_face_values(0).shape_hessian_component(dof_pair[0],
2930 q_point,
2931 component);
2932 if (dof_pair[1] != numbers::invalid_unsigned_int)
2933 value += 0.5 * get_fe_face_values(1).shape_hessian_component(dof_pair[1],
2934 q_point,
2935 component);
2936
2937 return value;
2938}
2939
2940
2941
2942template <int dim, int spacedim>
2945 const unsigned int interface_dof_index,
2946 const unsigned int q_point,
2947 const unsigned int component) const
2948{
2949 const auto dof_pair = dofmap[interface_dof_index];
2950
2951 if (at_boundary())
2952 return get_fe_face_values(0).shape_grad_component(dof_pair[0],
2953 q_point,
2954 component);
2955
2957
2958 if (dof_pair[0] != numbers::invalid_unsigned_int)
2959 value += get_fe_face_values(0).shape_grad_component(dof_pair[0],
2960 q_point,
2961 component);
2962 if (dof_pair[1] != numbers::invalid_unsigned_int)
2963 value -= get_fe_face_values(1).shape_grad_component(dof_pair[1],
2964 q_point,
2965 component);
2966
2967 return value;
2968}
2969
2970
2971
2972template <int dim, int spacedim>
2975 const unsigned int interface_dof_index,
2976 const unsigned int q_point,
2977 const unsigned int component) const
2978{
2979 const auto dof_pair = dofmap[interface_dof_index];
2980
2981 if (at_boundary())
2982 return get_fe_face_values(0).shape_hessian_component(dof_pair[0],
2983 q_point,
2984 component);
2985
2987
2988 if (dof_pair[0] != numbers::invalid_unsigned_int)
2989 value += get_fe_face_values(0).shape_hessian_component(dof_pair[0],
2990 q_point,
2991 component);
2992 if (dof_pair[1] != numbers::invalid_unsigned_int)
2993 value -= get_fe_face_values(1).shape_hessian_component(dof_pair[1],
2994 q_point,
2995 component);
2996
2997 return value;
2998}
2999
3000
3001
3002template <int dim, int spacedim>
3005 const unsigned int interface_dof_index,
3006 const unsigned int q_point,
3007 const unsigned int component) const
3008{
3009 const auto dof_pair = dofmap[interface_dof_index];
3010
3011 if (at_boundary())
3012 return get_fe_face_values(0).shape_3rd_derivative_component(dof_pair[0],
3013 q_point,
3014 component);
3015
3017
3018 if (dof_pair[0] != numbers::invalid_unsigned_int)
3019 value += get_fe_face_values(0).shape_3rd_derivative_component(dof_pair[0],
3020 q_point,
3021 component);
3022 if (dof_pair[1] != numbers::invalid_unsigned_int)
3023 value -= get_fe_face_values(1).shape_3rd_derivative_component(dof_pair[1],
3024 q_point,
3025 component);
3026
3027 return value;
3028}
3029
3030
3031
3032template <int dim, int spacedim>
3033template <class InputVector>
3034void
3036 const InputVector &fe_function,
3037 std::vector<typename InputVector::value_type> &values) const
3038{
3039 AssertDimension(values.size(), n_quadrature_points);
3040
3042 this->operator[](scalar).get_jump_in_function_values(fe_function, values);
3043}
3044
3045
3046
3047template <int dim, int spacedim>
3048template <class InputVector>
3049void
3051 const InputVector &fe_function,
3053 const
3054{
3055 AssertDimension(gradients.size(), n_quadrature_points);
3056
3058 this->operator[](scalar).get_jump_in_function_gradients(fe_function,
3059 gradients);
3060}
3061
3062
3063
3064template <int dim, int spacedim>
3065template <class InputVector>
3066void
3068 const InputVector &fe_function,
3070 const
3071{
3072 AssertDimension(hessians.size(), n_quadrature_points);
3073
3075 this->operator[](scalar).get_jump_in_function_hessians(fe_function, hessians);
3076}
3077
3078
3079
3080template <int dim, int spacedim>
3081template <class InputVector>
3082void
3084 const InputVector &fe_function,
3086 &third_derivatives) const
3087{
3088 AssertDimension(third_derivatives.size(), n_quadrature_points);
3089
3091 this->operator[](scalar).get_jump_in_function_third_derivatives(
3092 fe_function, third_derivatives);
3093}
3094
3095
3096
3097template <int dim, int spacedim>
3098template <class InputVector>
3099void
3101 const InputVector &fe_function,
3102 std::vector<typename InputVector::value_type> &values) const
3103{
3104 AssertDimension(values.size(), n_quadrature_points);
3105
3107 this->operator[](scalar).get_average_of_function_values(fe_function, values);
3108}
3109
3110
3111
3112template <int dim, int spacedim>
3113template <class InputVector>
3114void
3116 const InputVector &fe_function,
3118 const
3119{
3120 AssertDimension(gradients.size(), n_quadrature_points);
3121
3123 this->operator[](scalar).get_average_of_function_gradients(fe_function,
3124 gradients);
3125}
3126
3127
3128
3129template <int dim, int spacedim>
3130template <class InputVector>
3131void
3133 const InputVector &fe_function,
3135 const
3136{
3137 AssertDimension(hessians.size(), n_quadrature_points);
3138
3140 this->operator[](scalar).get_average_of_function_hessians(fe_function,
3141 hessians);
3142}
3143
3144
3145
3146/*------------ Inline functions: FEInterfaceValues------------*/
3147template <int dim, int spacedim>
3150 const FEValuesExtractors::Scalar &scalar) const
3151{
3152 const unsigned int n_components =
3153 (this->has_hp_capabilities() ? this->get_fe_collection().n_components() :
3154 this->get_fe().n_components());
3155 (void)n_components;
3156 AssertIndexRange(scalar.component, n_components);
3157 return FEInterfaceViews::Scalar<dim, spacedim>(*this, scalar.component);
3158}
3159
3160
3161
3162template <int dim, int spacedim>
3165 const FEValuesExtractors::Vector &vector) const
3166{
3167 const unsigned int n_components =
3168 (this->has_hp_capabilities() ? this->get_fe_collection().n_components() :
3169 this->get_fe().n_components());
3170 const unsigned int n_vectors =
3173 0);
3174 (void)n_components;
3175 (void)n_vectors;
3176 AssertIndexRange(vector.first_vector_component, n_vectors);
3178 vector.first_vector_component);
3179}
3180
3181
3182
3183namespace FEInterfaceViews
3184{
3185 template <int dim, int spacedim>
3187 const FEInterfaceValues<dim, spacedim> &fe_interface)
3188 : fe_interface(&fe_interface)
3189 {}
3190
3191
3192
3193 template <int dim, int spacedim>
3194 Scalar<dim, spacedim>::Scalar(
3195 const FEInterfaceValues<dim, spacedim> &fe_interface,
3196 const unsigned int component)
3197 : Base<dim, spacedim>(fe_interface)
3198 , extractor(component)
3199 {}
3200
3201
3202
3203 template <int dim, int spacedim>
3204 template <class InputVector, class OutputVector>
3205 void
3206 Base<dim, spacedim>::get_local_dof_values(
3207 const InputVector &dof_values,
3208 OutputVector &local_dof_values) const
3209 {
3210 const auto &interface_dof_indices =
3211 this->fe_interface->get_interface_dof_indices();
3212
3213 AssertDimension(interface_dof_indices.size(), local_dof_values.size());
3214
3215 for (const unsigned int i : this->fe_interface->dof_indices())
3216 local_dof_values[i] = dof_values(interface_dof_indices[i]);
3217 }
3218
3219
3220
3221 template <int dim, int spacedim>
3222 typename Scalar<dim, spacedim>::value_type
3223 Scalar<dim, spacedim>::value(const bool here_or_there,
3224 const unsigned int interface_dof_index,
3225 const unsigned int q_point) const
3226 {
3227 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3228
3229 if (here_or_there && dof_pair[0] != numbers::invalid_unsigned_int)
3230 return (*(this->fe_interface->fe_face_values))[extractor].value(
3231 dof_pair[0], q_point);
3232
3233 if (!here_or_there && dof_pair[1] != numbers::invalid_unsigned_int)
3234 return (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3235 dof_pair[1], q_point);
3236
3237 return 0.0;
3238 }
3239
3240
3241
3242 template <int dim, int spacedim>
3243 typename Scalar<dim, spacedim>::value_type
3244 Scalar<dim, spacedim>::jump_in_values(const unsigned int interface_dof_index,
3245 const unsigned int q_point) const
3246 {
3247 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3248
3249 value_type value = 0.0;
3250
3251 if (dof_pair[0] != numbers::invalid_unsigned_int)
3252 value +=
3253 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3254 q_point);
3255
3256 if (dof_pair[1] != numbers::invalid_unsigned_int)
3257 value -=
3258 (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3259 dof_pair[1], q_point);
3260
3261 return value;
3262 }
3263
3264
3265
3266 template <int dim, int spacedim>
3267 typename Scalar<dim, spacedim>::value_type
3268 Scalar<dim, spacedim>::average_of_values(
3269 const unsigned int interface_dof_index,
3270 const unsigned int q_point) const
3271 {
3272 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3273
3274 if (this->fe_interface->at_boundary())
3275 return (*(this->fe_interface->fe_face_values))[extractor].value(
3276 dof_pair[0], q_point);
3277
3278 value_type value = 0.0;
3279
3280 if (dof_pair[0] != numbers::invalid_unsigned_int)
3281 value +=
3282 0.5 *
3283 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3284 q_point);
3285
3286 if (dof_pair[1] != numbers::invalid_unsigned_int)
3287 value +=
3288 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3289 dof_pair[1], q_point);
3290
3291 return value;
3292 }
3293
3294
3295
3296 template <int dim, int spacedim>
3297 typename Scalar<dim, spacedim>::gradient_type
3298 Scalar<dim, spacedim>::average_of_gradients(
3299 const unsigned int interface_dof_index,
3300 const unsigned int q_point) const
3301 {
3302 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3303
3304 if (this->fe_interface->at_boundary())
3305 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
3306 dof_pair[0], q_point);
3307
3308 gradient_type value;
3309
3310 if (dof_pair[0] != numbers::invalid_unsigned_int)
3311 value +=
3312 0.5 *
3313 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
3314 q_point);
3315
3316 if (dof_pair[1] != numbers::invalid_unsigned_int)
3317 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3318 .gradient(dof_pair[1], q_point);
3319
3320 return value;
3321 }
3322
3323
3324
3325 template <int dim, int spacedim>
3326 typename Scalar<dim, spacedim>::gradient_type
3327 Scalar<dim, spacedim>::jump_in_gradients(
3328 const unsigned int interface_dof_index,
3329 const unsigned int q_point) const
3330 {
3331 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3332
3333 if (this->fe_interface->at_boundary())
3334 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
3335 dof_pair[0], q_point);
3336
3337 gradient_type value;
3338
3339 if (dof_pair[0] != numbers::invalid_unsigned_int)
3340 value +=
3341 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
3342 q_point);
3343
3344 if (dof_pair[1] != numbers::invalid_unsigned_int)
3345 value -=
3346 (*(this->fe_interface->fe_face_values_neighbor))[extractor].gradient(
3347 dof_pair[1], q_point);
3348
3349 return value;
3350 }
3351
3352
3353
3354 template <int dim, int spacedim>
3355 typename Scalar<dim, spacedim>::hessian_type
3356 Scalar<dim, spacedim>::average_of_hessians(
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 (this->fe_interface->at_boundary())
3363 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
3364 dof_pair[0], q_point);
3365
3366 hessian_type value;
3367
3368 if (dof_pair[0] != numbers::invalid_unsigned_int)
3369 value +=
3370 0.5 *
3371 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
3372 q_point);
3373
3374 if (dof_pair[1] != numbers::invalid_unsigned_int)
3375 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3376 .hessian(dof_pair[1], q_point);
3377
3378 return value;
3379 }
3380
3381
3382
3383 template <int dim, int spacedim>
3384 typename Scalar<dim, spacedim>::third_derivative_type
3385 Scalar<dim, spacedim>::jump_in_third_derivatives(
3386 const unsigned int interface_dof_index,
3387 const unsigned int q_point) const
3388 {
3389 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3390
3391 if (this->fe_interface->at_boundary())
3392 return (*(this->fe_interface->fe_face_values))[extractor]
3393 .third_derivative(dof_pair[0], q_point);
3394
3395 third_derivative_type value;
3396
3397 if (dof_pair[0] != numbers::invalid_unsigned_int)
3398 value +=
3399 (*(this->fe_interface->fe_face_values))[extractor].third_derivative(
3400 dof_pair[0], q_point);
3401
3402 if (dof_pair[1] != numbers::invalid_unsigned_int)
3403 value -= (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3404 .third_derivative(dof_pair[1], q_point);
3405
3406 return value;
3407 }
3408
3409
3410
3411 template <int dim, int spacedim>
3412 typename Scalar<dim, spacedim>::hessian_type
3413 Scalar<dim, spacedim>::jump_in_hessians(
3414 const unsigned int interface_dof_index,
3415 const unsigned int q_point) const
3416 {
3417 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3418
3419 if (this->fe_interface->at_boundary())
3420 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
3421 dof_pair[0], q_point);
3422
3423 hessian_type value;
3424
3425 if (dof_pair[0] != numbers::invalid_unsigned_int)
3426 value +=
3427 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
3428 q_point);
3429
3430 if (dof_pair[1] != numbers::invalid_unsigned_int)
3431 value -=
3432 (*(this->fe_interface->fe_face_values_neighbor))[extractor].hessian(
3433 dof_pair[1], q_point);
3434
3435 return value;
3436 }
3437
3438
3439
3440 template <int dim, int spacedim>
3441 template <class InputVector>
3442 void
3443 Scalar<dim, spacedim>::get_function_values_from_local_dof_values(
3444 const bool here_or_there,
3445 const InputVector &local_dof_values,
3446 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3447 const
3448 {
3449 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
3450
3451 for (const auto dof_index : this->fe_interface->dof_indices())
3452 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3453 {
3454 if (dof_index == 0)
3455 values[q_index] = 0.;
3456
3457 values[q_index] += local_dof_values[dof_index] *
3458 value(here_or_there, dof_index, q_index);
3459 }
3460 }
3461
3462
3463
3464 template <int dim, int spacedim>
3465 template <class InputVector>
3466 void
3467 Scalar<dim, spacedim>::get_function_values(
3468 const bool here_or_there,
3469 const InputVector &fe_function,
3470 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3471 const
3472 {
3473 std::vector<typename InputVector::value_type> local_dof_values(
3474 this->fe_interface->n_current_interface_dofs());
3475 this->get_local_dof_values(fe_function, local_dof_values);
3476
3477 get_function_values_from_local_dof_values(here_or_there,
3478 local_dof_values,
3479 values);
3480 }
3481
3482
3483
3484 template <int dim, int spacedim>
3485 template <class InputVector>
3486 void
3487 Scalar<dim, spacedim>::get_jump_in_function_values_from_local_dof_values(
3488 const InputVector &local_dof_values,
3489 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3490 const
3491 {
3492 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
3493
3494 for (const auto dof_index : this->fe_interface->dof_indices())
3495 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3496 {
3497 if (dof_index == 0)
3498 values[q_index] = 0.;
3499
3500 values[q_index] +=
3501 local_dof_values[dof_index] * jump_in_values(dof_index, q_index);
3502 }
3503 }
3504
3505
3506
3507 template <int dim, int spacedim>
3508 template <class InputVector>
3509 void
3510 Scalar<dim, spacedim>::get_jump_in_function_values(
3511 const InputVector &fe_function,
3512 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3513 const
3514 {
3515 std::vector<typename InputVector::value_type> local_dof_values(
3516 this->fe_interface->n_current_interface_dofs());
3517 this->get_local_dof_values(fe_function, local_dof_values);
3518
3519 get_jump_in_function_values_from_local_dof_values(local_dof_values, values);
3520 }
3521
3522
3523
3524 template <int dim, int spacedim>
3525 template <class InputVector>
3526 void
3527 Scalar<dim, spacedim>::get_jump_in_function_gradients_from_local_dof_values(
3528 const InputVector &local_dof_values,
3529 std::vector<solution_gradient_type<typename InputVector::value_type>>
3530 &gradients) const
3531 {
3532 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
3533
3534 for (const auto dof_index : this->fe_interface->dof_indices())
3535 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3536 {
3537 if (dof_index == 0)
3538 gradients[q_index] = 0.;
3539
3540 gradients[q_index] +=
3541 local_dof_values[dof_index] * jump_in_gradients(dof_index, q_index);
3542 }
3543 }
3544
3545
3546
3547 template <int dim, int spacedim>
3548 template <class InputVector>
3549 void
3550 Scalar<dim, spacedim>::get_jump_in_function_gradients(
3551 const InputVector &fe_function,
3552 std::vector<solution_gradient_type<typename InputVector::value_type>>
3553 &gradients) const
3554 {
3555 std::vector<typename InputVector::value_type> local_dof_values(
3556 this->fe_interface->n_current_interface_dofs());
3557 this->get_local_dof_values(fe_function, local_dof_values);
3558
3559 get_jump_in_function_gradients_from_local_dof_values(local_dof_values,
3560 gradients);
3561 }
3562
3563
3564
3565 template <int dim, int spacedim>
3566 template <class InputVector>
3567 void
3568 Scalar<dim, spacedim>::get_average_of_function_values_from_local_dof_values(
3569 const InputVector &local_dof_values,
3570 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3571 const
3572 {
3573 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
3574
3575 for (const auto dof_index : this->fe_interface->dof_indices())
3576 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3577 {
3578 if (dof_index == 0)
3579 values[q_index] = 0.;
3580
3581 values[q_index] +=
3582 local_dof_values[dof_index] * average_of_values(dof_index, q_index);
3583 }
3584 }
3585
3586
3587
3588 template <int dim, int spacedim>
3589 template <class InputVector>
3590 void
3591 Scalar<dim, spacedim>::get_average_of_function_values(
3592 const InputVector &fe_function,
3593 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3594 const
3595 {
3596 std::vector<typename InputVector::value_type> local_dof_values(
3597 this->fe_interface->n_current_interface_dofs());
3598 this->get_local_dof_values(fe_function, local_dof_values);
3599
3600 get_average_of_function_values_from_local_dof_values(local_dof_values,
3601 values);
3602 }
3603
3604
3605
3606 template <int dim, int spacedim>
3607 template <class InputVector>
3608 void
3609 Scalar<dim, spacedim>::
3610 get_average_of_function_gradients_from_local_dof_values(
3611 const InputVector &local_dof_values,
3612 std::vector<solution_gradient_type<typename InputVector::value_type>>
3613 &gradients) const
3614 {
3615 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
3616
3617 for (const auto dof_index : this->fe_interface->dof_indices())
3618 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3619 {
3620 if (dof_index == 0)
3621 gradients[q_index] = 0.;
3622
3623 gradients[q_index] += local_dof_values[dof_index] *
3624 average_of_gradients(dof_index, q_index);
3625 }
3626 }
3627
3628
3629
3630 template <int dim, int spacedim>
3631 template <class InputVector>
3632 void
3633 Scalar<dim, spacedim>::get_average_of_function_gradients(
3634 const InputVector &fe_function,
3635 std::vector<solution_gradient_type<typename InputVector::value_type>>
3636 &gradients) const
3637 {
3638 std::vector<typename InputVector::value_type> local_dof_values(
3639 this->fe_interface->n_current_interface_dofs());
3640 this->get_local_dof_values(fe_function, local_dof_values);
3641
3642 get_average_of_function_gradients_from_local_dof_values(local_dof_values,
3643 gradients);
3644 }
3645
3646
3647
3648 template <int dim, int spacedim>
3649 template <class InputVector>
3650 void
3651 Scalar<dim, spacedim>::get_jump_in_function_hessians_from_local_dof_values(
3652 const InputVector &local_dof_values,
3653 std::vector<solution_hessian_type<typename InputVector::value_type>>
3654 &hessians) const
3655 {
3656 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
3657
3658 for (const auto dof_index : this->fe_interface->dof_indices())
3659 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3660 {
3661 if (dof_index == 0)
3662 hessians[q_index] = 0.;
3663
3664 hessians[q_index] +=
3665 local_dof_values[dof_index] * jump_in_hessians(dof_index, q_index);
3666 }
3667 }
3668
3669
3670
3671 template <int dim, int spacedim>
3672 template <class InputVector>
3673 void
3674 Scalar<dim, spacedim>::get_jump_in_function_hessians(
3675 const InputVector &fe_function,
3676 std::vector<solution_hessian_type<typename InputVector::value_type>>
3677 &hessians) const
3678 {
3679 std::vector<typename InputVector::value_type> local_dof_values(
3680 this->fe_interface->n_current_interface_dofs());
3681 this->get_local_dof_values(fe_function, local_dof_values);
3682
3683 get_jump_in_function_hessians_from_local_dof_values(local_dof_values,
3684 hessians);
3685 }
3686
3687
3688
3689 template <int dim, int spacedim>
3690 template <class InputVector>
3691 void
3692 Scalar<dim, spacedim>::get_average_of_function_hessians_from_local_dof_values(
3693 const InputVector &local_dof_values,
3694 std::vector<solution_hessian_type<typename InputVector::value_type>>
3695 &hessians) const
3696 {
3697 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
3698
3699 for (const unsigned int dof_index : this->fe_interface->dof_indices())
3700 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3701 {
3702 if (dof_index == 0)
3703 hessians[q_index] = 0.;
3704
3705 hessians[q_index] += local_dof_values[dof_index] *
3706 average_of_hessians(dof_index, q_index);
3707 }
3708 }
3709
3710
3711
3712 template <int dim, int spacedim>
3713 template <class InputVector>
3714 void
3715 Scalar<dim, spacedim>::get_average_of_function_hessians(
3716 const InputVector &fe_function,
3717 std::vector<solution_hessian_type<typename InputVector::value_type>>
3718 &hessians) const
3719 {
3720 std::vector<typename InputVector::value_type> local_dof_values(
3721 this->fe_interface->n_current_interface_dofs());
3722 this->get_local_dof_values(fe_function, local_dof_values);
3723
3724 get_average_of_function_hessians_from_local_dof_values(local_dof_values,
3725 hessians);
3726 }
3727
3728
3729
3730 template <int dim, int spacedim>
3731 template <class InputVector>
3732 void
3733 Scalar<dim, spacedim>::
3734 get_jump_in_function_third_derivatives_from_local_dof_values(
3735 const InputVector &local_dof_values,
3736 std::vector<
3737 solution_third_derivative_type<typename InputVector::value_type>>
3738 &third_derivatives) const
3739 {
3740 AssertDimension(third_derivatives.size(),
3741 this->fe_interface->n_quadrature_points);
3742
3743 for (const unsigned int dof_index : this->fe_interface->dof_indices())
3744 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3745 {
3746 if (dof_index == 0)
3747 third_derivatives[q_index] = 0.;
3748
3749 third_derivatives[q_index] +=
3750 local_dof_values[dof_index] *
3751 jump_in_third_derivatives(dof_index, q_index);
3752 }
3753 }
3754
3755
3756
3757 template <int dim, int spacedim>
3758 template <class InputVector>
3759 void
3760 Scalar<dim, spacedim>::get_jump_in_function_third_derivatives(
3761 const InputVector &fe_function,
3762 std::vector<
3763 solution_third_derivative_type<typename InputVector::value_type>>
3764 &third_derivatives) const
3765 {
3766 std::vector<typename InputVector::value_type> local_dof_values(
3767 this->fe_interface->n_current_interface_dofs());
3768 this->get_local_dof_values(fe_function, local_dof_values);
3769
3770 get_jump_in_function_third_derivatives_from_local_dof_values(
3771 local_dof_values, third_derivatives);
3772 }
3773
3774
3775
3776 template <int dim, int spacedim>
3778 const FEInterfaceValues<dim, spacedim> &fe_interface,
3779 const unsigned int first_vector_component)
3780 : Base<dim, spacedim>(fe_interface)
3781 , extractor(first_vector_component)
3782 {}
3783
3784
3785
3786 template <int dim, int spacedim>
3788 Vector<dim, spacedim>::value(const bool here_or_there,
3789 const unsigned int interface_dof_index,
3790 const unsigned int q_point) const
3791 {
3792 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3793
3794 if (here_or_there && dof_pair[0] != numbers::invalid_unsigned_int)
3795 return (*(this->fe_interface->fe_face_values))[extractor].value(
3796 dof_pair[0], q_point);
3797
3798 if (!here_or_there && dof_pair[1] != numbers::invalid_unsigned_int)
3799 return (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3800 dof_pair[1], q_point);
3801
3802 return value_type();
3803 }
3804
3805
3806
3807 template <int dim, int spacedim>
3809 Vector<dim, spacedim>::jump_in_values(const unsigned int interface_dof_index,
3810 const unsigned int q_point) const
3811 {
3812 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3813
3814 value_type value;
3815
3816 if (dof_pair[0] != numbers::invalid_unsigned_int)
3817 value +=
3818 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3819 q_point);
3820
3821 if (dof_pair[1] != numbers::invalid_unsigned_int)
3822 value -=
3823 (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3824 dof_pair[1], q_point);
3825
3826 return value;
3827 }
3828
3829
3830
3831 template <int dim, int spacedim>
3834 const unsigned int interface_dof_index,
3835 const unsigned int q_point) const
3836 {
3837 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3838
3839 if (this->fe_interface->at_boundary())
3840 return (*(this->fe_interface->fe_face_values))[extractor].value(
3841 dof_pair[0], q_point);
3842
3843 value_type value;
3844
3845 if (dof_pair[0] != numbers::invalid_unsigned_int)
3846 value +=
3847 0.5 *
3848 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3849 q_point);
3850
3851 if (dof_pair[1] != numbers::invalid_unsigned_int)
3852 value +=
3853 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3854 dof_pair[1], q_point);
3855
3856 return value;
3857 }
3858
3859
3860
3861 template <int dim, int spacedim>
3864 const unsigned int interface_dof_index,
3865 const unsigned int q_point) const
3866 {
3867 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3868
3869 if (this->fe_interface->at_boundary())
3870 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
3871 dof_pair[0], q_point);
3872
3873 gradient_type value;
3874
3875 if (dof_pair[0] != numbers::invalid_unsigned_int)
3876 value +=
3877 0.5 *
3878 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
3879 q_point);
3880
3881 if (dof_pair[1] != numbers::invalid_unsigned_int)
3882 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3883 .gradient(dof_pair[1], q_point);
3884
3885 return value;
3886 }
3887
3888
3889
3890 template <int dim, int spacedim>
3893 const unsigned int interface_dof_index,
3894 const unsigned int q_point) const
3895 {
3896 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3897
3898 if (this->fe_interface->at_boundary())
3899 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
3900 dof_pair[0], q_point);
3901
3902 gradient_type value;
3903
3904 if (dof_pair[0] != numbers::invalid_unsigned_int)
3905 value +=
3906 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
3907 q_point);
3908
3909 if (dof_pair[1] != numbers::invalid_unsigned_int)
3910 value -=
3911 (*(this->fe_interface->fe_face_values_neighbor))[extractor].gradient(
3912 dof_pair[1], q_point);
3913
3914 return value;
3915 }
3916
3917
3918
3919 template <int dim, int spacedim>
3921 Vector<dim, spacedim>::jump_gradient(const unsigned int interface_dof_index,
3922 const unsigned int q_point) const
3923 {
3924 return jump_in_gradients(interface_dof_index, q_point);
3925 }
3926
3927
3928
3929 template <int dim, int spacedim>
3932 const unsigned int interface_dof_index,
3933 const unsigned int q_point) const
3934 {
3935 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3936
3937 if (this->fe_interface->at_boundary())
3938 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
3939 dof_pair[0], q_point);
3940
3941 hessian_type value;
3942
3943 if (dof_pair[0] != numbers::invalid_unsigned_int)
3944 value +=
3945 0.5 *
3946 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
3947 q_point);
3948
3949 if (dof_pair[1] != numbers::invalid_unsigned_int)
3950 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3951 .hessian(dof_pair[1], q_point);
3952
3953 return value;
3954 }
3955
3956
3957
3958 template <int dim, int spacedim>
3960 Vector<dim, spacedim>::average_hessian(const unsigned int interface_dof_index,
3961 const unsigned int q_point) const
3962 {
3963 return average_of_hessians(interface_dof_index, q_point);
3964 }
3965
3966
3967
3968 template <int dim, int spacedim>
3971 const unsigned int interface_dof_index,
3972 const unsigned int q_point) const
3973 {
3974 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3975
3976 if (this->fe_interface->at_boundary())
3977 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
3978 dof_pair[0], q_point);
3979
3980 hessian_type value;
3981
3982 if (dof_pair[0] != numbers::invalid_unsigned_int)
3983 value +=
3984 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
3985 q_point);
3986
3987 if (dof_pair[1] != numbers::invalid_unsigned_int)
3988 value -=
3989 (*(this->fe_interface->fe_face_values_neighbor))[extractor].hessian(
3990 dof_pair[1], q_point);
3991
3992 return value;
3993 }
3994
3995
3996
3997 template <int dim, int spacedim>
3999 Vector<dim, spacedim>::jump_hessian(const unsigned int interface_dof_index,
4000 const unsigned int q_point) const
4001 {
4002 return jump_in_hessians(interface_dof_index, q_point);
4003 }
4004
4005
4006
4007 template <int dim, int spacedim>
4010 const unsigned int interface_dof_index,
4011 const unsigned int q_point) const
4012 {
4013 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4014
4015 if (this->fe_interface->at_boundary())
4016 return (*(this->fe_interface->fe_face_values))[extractor]
4017 .third_derivative(dof_pair[0], q_point);
4018
4019 third_derivative_type value;
4020
4021 if (dof_pair[0] != numbers::invalid_unsigned_int)
4022 value +=
4023 (*(this->fe_interface->fe_face_values))[extractor].third_derivative(
4024 dof_pair[0], q_point);
4025
4026 if (dof_pair[1] != numbers::invalid_unsigned_int)
4027 value -= (*(this->fe_interface->fe_face_values_neighbor))[extractor]
4028 .third_derivative(dof_pair[1], q_point);
4029
4030 return value;
4031 }
4032
4033
4034
4035 template <int dim, int spacedim>
4036 template <class InputVector>
4037 void
4039 const bool here_or_there,
4040 const InputVector &local_dof_values,
4041 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4042 const
4043 {
4044 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
4045
4046 for (const auto dof_index : this->fe_interface->dof_indices())
4047 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4048 {
4049 if (dof_index == 0)
4050 values[q_index] = 0.;
4051
4052 values[q_index] += local_dof_values[dof_index] *
4053 value(here_or_there, dof_index, q_index);
4054 }
4055 }
4056
4057
4058
4059 template <int dim, int spacedim>
4060 template <class InputVector>
4061 void
4063 const bool here_or_there,
4064 const InputVector &fe_function,
4065 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4066 const
4067 {
4068 std::vector<typename InputVector::value_type> local_dof_values(
4069 this->fe_interface->n_current_interface_dofs());
4070 this->get_local_dof_values(fe_function, local_dof_values);
4071
4072 get_function_values_from_local_dof_values(here_or_there,
4073 local_dof_values,
4074 values);
4075 }
4076
4077
4078
4079 template <int dim, int spacedim>
4080 template <class InputVector>
4081 void
4083 const InputVector &local_dof_values,
4084 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4085 const
4086 {
4087 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
4088
4089 for (const auto dof_index : this->fe_interface->dof_indices())
4090 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4091 {
4092 if (dof_index == 0)
4093 values[q_index] = 0.;
4094
4095 values[q_index] +=
4096 local_dof_values[dof_index] * jump_in_values(dof_index, q_index);
4097 }
4098 }
4099
4100
4101
4102 template <int dim, int spacedim>
4103 template <class InputVector>
4104 void
4106 const InputVector &fe_function,
4107 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4108 const
4109 {
4110 std::vector<typename InputVector::value_type> local_dof_values(
4111 this->fe_interface->n_current_interface_dofs());
4112 this->get_local_dof_values(fe_function, local_dof_values);
4113
4114 get_jump_in_function_values_from_local_dof_values(local_dof_values, values);
4115 }
4116
4117
4118
4119 template <int dim, int spacedim>
4120 template <class InputVector>
4121 void
4123 const InputVector &local_dof_values,
4124 std::vector<solution_gradient_type<typename InputVector::value_type>>
4125 &gradients) const
4126 {
4127 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
4128
4129 for (const auto dof_index : this->fe_interface->dof_indices())
4130 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4131 {
4132 if (dof_index == 0)
4133 gradients[q_index] = 0.;
4134
4135 gradients[q_index] +=
4136 local_dof_values[dof_index] * jump_in_gradients(dof_index, q_index);
4137 }
4138 }
4139
4140
4141
4142 template <int dim, int spacedim>
4143 template <class InputVector>
4144 void
4146 const InputVector &fe_function,
4147 std::vector<solution_gradient_type<typename InputVector::value_type>>
4148 &gradients) const
4149 {
4150 std::vector<typename InputVector::value_type> local_dof_values(
4151 this->fe_interface->n_current_interface_dofs());
4152 this->get_local_dof_values(fe_function, local_dof_values);
4153
4154 get_jump_in_function_gradients_from_local_dof_values(local_dof_values,
4155 gradients);
4156 }
4157
4158
4159
4160 template <int dim, int spacedim>
4161 template <class InputVector>
4162 void
4164 const InputVector &local_dof_values,
4165 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4166 const
4167 {
4168 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
4169
4170 for (const auto dof_index : this->fe_interface->dof_indices())
4171 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4172 {
4173 if (dof_index == 0)
4174 values[q_index] = 0.;
4175
4176 values[q_index] +=
4177 local_dof_values[dof_index] * average_of_values(dof_index, q_index);
4178 }
4179 }
4180
4181
4182
4183 template <int dim, int spacedim>
4184 template <class InputVector>
4185 void
4187 const InputVector &fe_function,
4188 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4189 const
4190 {
4191 std::vector<typename InputVector::value_type> local_dof_values(
4192 this->fe_interface->n_current_interface_dofs());
4193 this->get_local_dof_values(fe_function, local_dof_values);
4194
4195 get_average_of_function_values_from_local_dof_values(local_dof_values,
4196 values);
4197 }
4198
4199
4200
4201 template <int dim, int spacedim>
4202 template <class InputVector>
4203 void
4206 const InputVector &local_dof_values,
4207 std::vector<solution_gradient_type<typename InputVector::value_type>>
4208 &gradients) const
4209 {
4210 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
4211
4212 for (const auto dof_index : this->fe_interface->dof_indices())
4213 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4214 {
4215 if (dof_index == 0)
4216 gradients[q_index] = 0.;
4217
4218 gradients[q_index] += local_dof_values[dof_index] *
4219 average_of_gradients(dof_index, q_index);
4220 }
4221 }
4222
4223
4224
4225 template <int dim, int spacedim>
4226 template <class InputVector>
4227 void
4229 const InputVector &fe_function,
4230 std::vector<solution_gradient_type<typename InputVector::value_type>>
4231 &gradients) const
4232 {
4233 std::vector<typename InputVector::value_type> local_dof_values(
4234 this->fe_interface->n_current_interface_dofs());
4235 this->get_local_dof_values(fe_function, local_dof_values);
4236
4237 get_average_of_function_gradients_from_local_dof_values(local_dof_values,
4238 gradients);
4239 }
4240
4241
4242
4243 template <int dim, int spacedim>
4244 template <class InputVector>
4245 void
4247 const InputVector &local_dof_values,
4248 std::vector<solution_hessian_type<typename InputVector::value_type>>
4249 &hessians) const
4250 {
4251 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
4252
4253 for (const auto dof_index : this->fe_interface->dof_indices())
4254 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4255 {
4256 if (dof_index == 0)
4257 hessians[q_index] = 0.;
4258
4259 hessians[q_index] +=
4260 local_dof_values[dof_index] * jump_in_hessians(dof_index, q_index);
4261 }
4262 }
4263
4264
4265
4266 template <int dim, int spacedim>
4267 template <class InputVector>
4268 void
4270 const InputVector &fe_function,
4271 std::vector<solution_hessian_type<typename InputVector::value_type>>
4272 &hessians) const
4273 {
4274 std::vector<typename InputVector::value_type> local_dof_values(
4275 this->fe_interface->n_current_interface_dofs());
4276 this->get_local_dof_values(fe_function, local_dof_values);
4277
4278 get_jump_in_function_hessians_from_local_dof_values(local_dof_values,
4279 hessians);
4280 }
4281
4282
4283
4284 template <int dim, int spacedim>
4285 template <class InputVector>
4286 void
4288 const InputVector &local_dof_values,
4289 std::vector<solution_hessian_type<typename InputVector::value_type>>
4290 &hessians) const
4291 {
4292 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
4293
4294 for (const unsigned int dof_index : this->fe_interface->dof_indices())
4295 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4296 {
4297 if (dof_index == 0)
4298 hessians[q_index] = 0.;
4299
4300 hessians[q_index] += local_dof_values[dof_index] *
4301 average_of_hessians(dof_index, q_index);
4302 }
4303 }
4304
4305
4306
4307 template <int dim, int spacedim>
4308 template <class InputVector>
4309 void
4311 const InputVector &fe_function,
4312 std::vector<solution_hessian_type<typename InputVector::value_type>>
4313 &hessians) const
4314 {
4315 std::vector<typename InputVector::value_type> local_dof_values(
4316 this->fe_interface->n_current_interface_dofs());
4317 this->get_local_dof_values(fe_function, local_dof_values);
4318
4319 get_average_of_function_hessians_from_local_dof_values(local_dof_values,
4320 hessians);
4321 }
4322
4323
4324
4325 template <int dim, int spacedim>
4326 template <class InputVector>
4327 void
4330 const InputVector &local_dof_values,
4331 std::vector<
4332 solution_third_derivative_type<typename InputVector::value_type>>
4333 &third_derivatives) const
4334 {
4335 AssertDimension(third_derivatives.size(),
4336 this->fe_interface->n_quadrature_points);
4337
4338 for (const unsigned int dof_index : this->fe_interface->dof_indices())
4339 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4340 {
4341 if (dof_index == 0)
4342 third_derivatives[q_index] = 0.;
4343
4344 third_derivatives[q_index] +=
4345 local_dof_values[dof_index] *
4346 jump_in_third_derivatives(dof_index, q_index);
4347 }
4348 }
4349
4350
4351
4352 template <int dim, int spacedim>
4353 template <class InputVector>
4354 void
4356 const InputVector &fe_function,
4357 std::vector<
4358 solution_third_derivative_type<typename InputVector::value_type>>
4359 &third_derivatives) const
4360 {
4361 std::vector<typename InputVector::value_type> local_dof_values(
4362 this->fe_interface->n_current_interface_dofs());
4363 this->get_local_dof_values(fe_function, local_dof_values);
4364
4365 get_jump_in_function_third_derivatives_from_local_dof_values(
4366 local_dof_values, third_derivatives);
4367 }
4368} // namespace FEInterfaceViews
4369
4370#endif // DOXYGEN
4371
4373
4374#endif
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access > > &cell, const unsigned int face_no)
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
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
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
FEInterfaceViews::Vector< dim, spacedim > operator[](const FEValuesExtractors::Vector &vector) 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
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)
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 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:317
friend class Vector
Definition vector.h:1057
Number value_type
Definition vector.h:129
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:503
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:504
unsigned int cell_index
#define Assert(cond, exc)
static ::ExceptionBase & ExcOnlyAvailableWithHP()
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:495
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:285
typename ::internal::FEInterfaceViews::ViewType< dim, spacedim, Extractor >::type View
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
spacedim & mapping
for(unsigned int j=best_vertex+1;j< vertices.size();++j) if(vertices_to_use[j])
Definition hp.h:118
const types::fe_index invalid_fe_index
Definition types.h:244
static const unsigned int invalid_unsigned_int
Definition types.h:221
boost::integer_range< IncrementableType > iota_view
Definition iota_view.h:46
STL namespace.
typename internal::ProductTypeImpl< std::decay_t< T >, std::decay_t< U > >::type type