Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
fe_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
208 jump(const unsigned int interface_dof_index,
209 const unsigned int q_point) const;
210
222 jump_in_gradients(const unsigned int interface_dof_index,
223 const unsigned int q_point) const;
224
232 jump_gradient(const unsigned int interface_dof_index,
233 const unsigned int q_point) const;
234
246 jump_in_hessians(const unsigned int interface_dof_index,
247 const unsigned int q_point) const;
248
256 jump_hessian(const unsigned int interface_dof_index,
257 const unsigned int q_point) const;
258
270 jump_in_third_derivatives(const unsigned int interface_dof_index,
271 const unsigned int q_point) const;
272
280 jump_3rd_derivative(const unsigned int interface_dof_index,
281 const unsigned int q_point) const;
282
301 average_of_values(const unsigned int interface_dof_index,
302 const unsigned int q_point) const;
303
311 average_value(const unsigned int interface_dof_index,
312 const unsigned int q_point) const;
313
321 average(const unsigned int interface_dof_index,
322 const unsigned int q_point) const;
323
335 average_of_gradients(const unsigned int interface_dof_index,
336 const unsigned int q_point) const;
337
345 average_gradient(const unsigned int interface_dof_index,
346 const unsigned int q_point) const;
347
360 average_of_hessians(const unsigned int interface_dof_index,
361 const unsigned int q_point) const;
362
370 average_hessian(const unsigned int interface_dof_index,
371 const unsigned int q_point) const;
372
399 template <class InputVector>
400 void
402 const bool here_or_there,
403 const InputVector &fe_function,
405 &values) const;
406
429 template <class InputVector>
430 void
432 const bool here_or_there,
433 const InputVector &local_dof_values,
435 &values) const;
436
457 template <class InputVector>
458 void
460 const InputVector &fe_function,
462 &values) const;
463
471 template <class InputVector>
472 void
474 const InputVector &local_dof_values,
476 &values) const;
477
491 template <class InputVector>
492 void
494 const InputVector &fe_function,
496 &gradients) const;
497
505 template <class InputVector>
506 void
508 const InputVector &local_dof_values,
510 &gradients) const;
511
525 template <class InputVector>
526 void
528 const InputVector &fe_function,
530 &hessians) const;
531
539 template <class InputVector>
540 void
542 const InputVector &local_dof_values,
544 &hessians) const;
545
560 template <class InputVector>
561 void
563 const InputVector &fe_function,
564 std::vector<
566 &third_derivatives) const;
567
575 template <class InputVector>
576 void
578 const InputVector &local_dof_values,
579 std::vector<
581 &third_derivatives) const;
582
603 template <class InputVector>
604 void
606 const InputVector &fe_function,
608 &values) const;
609
617 template <class InputVector>
618 void
620 const InputVector &local_dof_values,
622 &values) const;
623
637 template <class InputVector>
638 void
640 const InputVector &fe_function,
642 &gradients) const;
643
651 template <class InputVector>
652 void
654 const InputVector &local_dof_values,
656 &gradients) const;
657
671 template <class InputVector>
672 void
674 const InputVector &fe_function,
676 &hessians) const;
677
685 template <class InputVector>
686 void
688 const InputVector &local_dof_values,
690 &hessians) const;
691
694 private:
699 };
700
701
702
706 template <int dim, int spacedim = dim>
707 class Vector : public Base<dim, spacedim>
708 {
709 public:
715
722
730
738
745 template <typename Number>
747
754 template <typename Number>
757
764 template <typename Number>
767
774 template <typename Number>
777
782 const unsigned int first_vector_component);
783
809 value(const bool here_or_there,
810 const unsigned int interface_dof_index,
811 const unsigned int q_point) const;
812
830 jump_in_values(const unsigned int interface_dof_index,
831 const unsigned int q_point) const;
832
840 jump(const unsigned int interface_dof_index,
841 const unsigned int q_point) const;
842
853 jump_in_gradients(const unsigned int interface_dof_index,
854 const unsigned int q_point) const;
855
862 jump_gradient(const unsigned int interface_dof_index,
863 const unsigned int q_point) const;
864
876 jump_in_hessians(const unsigned int interface_dof_index,
877 const unsigned int q_point) const;
878
885 jump_hessian(const unsigned int interface_dof_index,
886 const unsigned int q_point) const;
887
899 jump_in_third_derivatives(const unsigned int interface_dof_index,
900 const unsigned int q_point) const;
901
909 jump_3rd_derivative(const unsigned int interface_dof_index,
910 const unsigned int q_point) const;
911
929 average_of_values(const unsigned int interface_dof_index,
930 const unsigned int q_point) const;
931
939 average(const unsigned int interface_dof_index,
940 const unsigned int q_point) const;
941
952 average_of_gradients(const unsigned int interface_dof_index,
953 const unsigned int q_point) const;
954
962 average_gradient(const unsigned int interface_dof_index,
963 const unsigned int q_point) const;
964
977 average_of_hessians(const unsigned int interface_dof_index,
978 const unsigned int q_point) const;
979
986 average_hessian(const unsigned int interface_dof_index,
987 const unsigned int q_point) const;
988
1015 template <class InputVector>
1016 void
1018 const bool here_or_there,
1019 const InputVector &fe_function,
1021 &values) const;
1022
1045 template <class InputVector>
1046 void
1048 const bool here_or_there,
1049 const InputVector &local_dof_values,
1051 &values) const;
1052
1073 template <class InputVector>
1074 void
1076 const InputVector &fe_function,
1078 &values) const;
1079
1087 template <class InputVector>
1088 void
1090 const InputVector &local_dof_values,
1092 &values) const;
1093
1107 template <class InputVector>
1108 void
1110 const InputVector &fe_function,
1112 &gradients) const;
1113
1121 template <class InputVector>
1122 void
1124 const InputVector &local_dof_values,
1126 &gradients) const;
1127
1141 template <class InputVector>
1142 void
1144 const InputVector &fe_function,
1146 &hessians) const;
1147
1155 template <class InputVector>
1156 void
1158 const InputVector &local_dof_values,
1160 &hessians) const;
1161
1176 template <class InputVector>
1177 void
1179 const InputVector &fe_function,
1180 std::vector<
1182 &third_derivatives) const;
1183
1191 template <class InputVector>
1192 void
1194 const InputVector &local_dof_values,
1195 std::vector<
1197 &third_derivatives) const;
1198
1219 template <class InputVector>
1220 void
1222 const InputVector &fe_function,
1224 &values) const;
1225
1233 template <class InputVector>
1234 void
1236 const InputVector &local_dof_values,
1238 &values) const;
1239
1253 template <class InputVector>
1254 void
1256 const InputVector &fe_function,
1258 &gradients) const;
1259
1267 template <class InputVector>
1268 void
1270 const InputVector &local_dof_values,
1272 &gradients) const;
1273
1287 template <class InputVector>
1288 void
1290 const InputVector &fe_function,
1292 &hessians) const;
1293
1301 template <class InputVector>
1302 void
1304 const InputVector &local_dof_values,
1306 &hessians) const;
1307
1310 private:
1315 };
1316} // namespace FEInterfaceViews
1317
1318
1319namespace internal
1320{
1322 {
1327 template <int dim, int spacedim, typename Extractor>
1329 {};
1330
1338 template <int dim, int spacedim>
1339 struct ViewType<dim, spacedim, FEValuesExtractors::Scalar>
1340 {
1341 using type = typename ::FEInterfaceViews::Scalar<dim, spacedim>;
1342 };
1343
1351 template <int dim, int spacedim>
1352 struct ViewType<dim, spacedim, FEValuesExtractors::Vector>
1353 {
1354 using type = typename ::FEInterfaceViews::Vector<dim, spacedim>;
1355 };
1356 } // namespace FEInterfaceViews
1357} // namespace internal
1358
1359namespace FEInterfaceViews
1360{
1365 template <int dim, int spacedim, typename Extractor>
1366 using View = typename ::internal::FEInterfaceViews::
1367 ViewType<dim, spacedim, Extractor>::type;
1368} // namespace FEInterfaceViews
1369
1370
1371
1396template <int dim, int spacedim = dim>
1398{
1399public:
1403 const unsigned int n_quadrature_points;
1404
1412 const Quadrature<dim - 1> & quadrature,
1413 const UpdateFlags update_flags);
1414
1422 const hp::QCollection<dim - 1> & quadrature,
1423 const UpdateFlags update_flags);
1424
1432 const Quadrature<dim - 1> & quadrature,
1433 const UpdateFlags update_flags);
1434
1440 const hp::MappingCollection<dim, spacedim> &mapping_collection,
1441 const hp::FECollection<dim, spacedim> & fe_collection,
1442 const hp::QCollection<dim - 1> & quadrature_collection,
1443 const UpdateFlags update_flags);
1444
1449 const hp::QCollection<dim - 1> &quadrature_collection,
1450 const UpdateFlags update_flags);
1451
1566 template <class CellIteratorType, class CellNeighborIteratorType>
1567 void
1568 reinit(const CellIteratorType & cell,
1569 const unsigned int face_no,
1570 const unsigned int sub_face_no,
1571 const CellNeighborIteratorType &cell_neighbor,
1572 const unsigned int face_no_neighbor,
1573 const unsigned int sub_face_no_neighbor,
1574 const unsigned int q_index = numbers::invalid_unsigned_int,
1575 const unsigned int mapping_index = numbers::invalid_unsigned_int,
1576 const unsigned int fe_index = numbers::invalid_unsigned_int);
1577
1604 template <class CellIteratorType>
1605 void
1606 reinit(const CellIteratorType &cell,
1607 const unsigned int face_no,
1608 const unsigned int q_index = numbers::invalid_unsigned_int,
1609 const unsigned int mapping_index = numbers::invalid_unsigned_int,
1610 const unsigned int fe_index = numbers::invalid_unsigned_int);
1611
1620 get_fe_face_values(const unsigned int cell_index) const;
1621
1627
1632 get_fe() const;
1633
1637 const Quadrature<dim - 1> &
1639
1645
1651
1655 const hp::QCollection<dim - 1> &
1657
1662 bool
1664
1670
1678 get_cell(const unsigned int cell_index) const;
1679
1687 unsigned int
1688 get_face_number(const unsigned int cell_index) const;
1689
1701 bool
1703
1715 double
1716 JxW(const unsigned int quadrature_point) const;
1717
1723 const std::vector<double> &
1725
1735 const std::vector<Tensor<1, spacedim>> &
1737
1747
1753 const std::vector<Point<spacedim>> &
1755
1766 unsigned
1768
1798
1809 std::vector<types::global_dof_index>
1811
1824 std::array<unsigned int, 2>
1825 interface_dof_to_dof_indices(const unsigned int interface_dof_index) const;
1826
1836 normal(const unsigned int q_point_index) const;
1837
1868 double
1869 shape_value(const bool here_or_there,
1870 const unsigned int interface_dof_index,
1871 const unsigned int q_point,
1872 const unsigned int component = 0) const;
1873
1901 double
1902 jump_in_shape_values(const unsigned int interface_dof_index,
1903 const unsigned int q_point,
1904 const unsigned int component = 0) const;
1905
1912 double
1913 jump(const unsigned int interface_dof_index,
1914 const unsigned int q_point,
1915 const unsigned int component = 0) const;
1916
1931 jump_in_shape_gradients(const unsigned int interface_dof_index,
1932 const unsigned int q_point,
1933 const unsigned int component = 0) const;
1934
1942 jump_gradient(const unsigned int interface_dof_index,
1943 const unsigned int q_point,
1944 const unsigned int component = 0) const;
1945
1961 jump_in_shape_hessians(const unsigned int interface_dof_index,
1962 const unsigned int q_point,
1963 const unsigned int component = 0) const;
1964
1972 jump_hessian(const unsigned int interface_dof_index,
1973 const unsigned int q_point,
1974 const unsigned int component = 0) const;
1975
1990 jump_in_shape_3rd_derivatives(const unsigned int interface_dof_index,
1991 const unsigned int q_point,
1992 const unsigned int component = 0) const;
1993
2001 jump_3rd_derivative(const unsigned int interface_dof_index,
2002 const unsigned int q_point,
2003 const unsigned int component = 0) const;
2004
2027 double
2028 average_of_shape_values(const unsigned int interface_dof_index,
2029 const unsigned int q_point,
2030 const unsigned int component = 0) const;
2031
2038 double
2039 average(const unsigned int interface_dof_index,
2040 const unsigned int q_point,
2041 const unsigned int component = 0) const;
2042
2057 average_of_shape_gradients(const unsigned int interface_dof_index,
2058 const unsigned int q_point,
2059 const unsigned int component = 0) const;
2060
2068 average_gradient(const unsigned int interface_dof_index,
2069 const unsigned int q_point,
2070 const unsigned int component = 0) const;
2071
2087 average_of_shape_hessians(const unsigned int interface_dof_index,
2088 const unsigned int q_point,
2089 const unsigned int component = 0) const;
2090
2098 average_hessian(const unsigned int interface_dof_index,
2099 const unsigned int q_point,
2100 const unsigned int component = 0) const;
2101
2121 template <class InputVector>
2122 void
2124 const InputVector & fe_function,
2125 std::vector<typename InputVector::value_type> &values) const;
2126
2135 template <class InputVector>
2136 void
2138 const InputVector &fe_function,
2140 &gradients) const;
2141
2149 template <class InputVector>
2150 void
2152 const InputVector &fe_function,
2154 &hessians) const;
2155
2164 template <class InputVector>
2165 void
2167 const InputVector &fe_function,
2169 &third_derivatives) const;
2170
2186 template <class InputVector>
2187 void
2189 const InputVector & fe_function,
2190 std::vector<typename InputVector::value_type> &values) const;
2191
2199 template <class InputVector>
2200 void
2202 const InputVector &fe_function,
2204 &gradients) const;
2205
2213 template <class InputVector>
2214 void
2216 const InputVector &fe_function,
2218 &hessians) const;
2219
2239
2248
2253private:
2257 std::vector<types::global_dof_index> interface_dof_indices;
2258
2264 std::vector<std::array<unsigned int, 2>> dofmap;
2265
2271
2278 // non-hp data
2283
2287 std::unique_ptr<FEFaceValues<dim, spacedim>> internal_fe_face_values;
2288
2292 std::unique_ptr<FESubfaceValues<dim, spacedim>> internal_fe_subface_values;
2293
2297 std::unique_ptr<FEFaceValues<dim, spacedim>> internal_fe_face_values_neighbor;
2298
2302 std::unique_ptr<FESubfaceValues<dim, spacedim>>
2304 // non-hp data
2306 // hp data
2311
2316 std::unique_ptr<hp::FEFaceValues<dim, spacedim>> internal_hp_fe_face_values;
2317
2322 std::unique_ptr<hp::FESubfaceValues<dim, spacedim>>
2324
2329 std::unique_ptr<hp::FEFaceValues<dim, spacedim>>
2331
2336 std::unique_ptr<hp::FESubfaceValues<dim, spacedim>>
2338
2346 "The current function doesn't make sense when used with a "
2347 "FEInterfaceValues object with hp-capabilities.");
2348
2356 "The current function doesn't make sense when used with a "
2357 "FEInterfaceValues object without hp-capabilities.");
2358 // hp data
2360
2361 /*
2362 * Make the view classes friends of this class, since they
2363 * access internal data.
2364 */
2365 template <int, int>
2367 template <int, int>
2369};
2370
2371
2372
2373#ifndef DOXYGEN
2374
2375/*---------------------- Inline functions ---------------------*/
2376
2377template <int dim, int spacedim>
2379 const Mapping<dim, spacedim> & mapping,
2381 const Quadrature<dim - 1> & quadrature,
2382 const UpdateFlags update_flags)
2383 : n_quadrature_points(quadrature.size())
2384 , fe_face_values(nullptr)
2385 , fe_face_values_neighbor(nullptr)
2386 , internal_fe_face_values(
2387 std::make_unique<FEFaceValues<dim, spacedim>>(mapping,
2388 fe,
2389 quadrature,
2390 update_flags))
2391 , internal_fe_subface_values(
2392 std::make_unique<FESubfaceValues<dim, spacedim>>(mapping,
2393 fe,
2394 quadrature,
2395 update_flags))
2396 , internal_fe_face_values_neighbor(
2397 std::make_unique<FEFaceValues<dim, spacedim>>(mapping,
2398 fe,
2399 quadrature,
2400 update_flags))
2401 , internal_fe_subface_values_neighbor(
2402 std::make_unique<FESubfaceValues<dim, spacedim>>(mapping,
2403 fe,
2404 quadrature,
2405 update_flags))
2406{}
2407
2408
2409
2410template <int dim, int spacedim>
2413 const Quadrature<dim - 1> & quadrature,
2414 const UpdateFlags update_flags)
2416 fe.reference_cell().template get_default_linear_mapping<dim, spacedim>(),
2417 fe,
2418 quadrature,
2419 update_flags)
2420{}
2421
2422
2423
2424template <int dim, int spacedim>
2426 const Mapping<dim, spacedim> & mapping,
2428 const hp::QCollection<dim - 1> & quadrature,
2429 const UpdateFlags update_flags)
2430 : n_quadrature_points(quadrature.max_n_quadrature_points())
2431 , fe_face_values(nullptr)
2432 , fe_face_values_neighbor(nullptr)
2433 , internal_fe_face_values(
2434 std::make_unique<FEFaceValues<dim, spacedim>>(mapping,
2435 fe,
2436 quadrature,
2437 update_flags))
2438 , internal_fe_subface_values(
2439 std::make_unique<FESubfaceValues<dim, spacedim>>(mapping,
2440 fe,
2441 quadrature,
2442 update_flags))
2443 , internal_fe_face_values_neighbor(
2444 std::make_unique<FEFaceValues<dim, spacedim>>(mapping,
2445 fe,
2446 quadrature[0],
2447 update_flags))
2448 , internal_fe_subface_values_neighbor(
2449 std::make_unique<FESubfaceValues<dim, spacedim>>(mapping,
2450 fe,
2451 quadrature[0],
2452 update_flags))
2453{}
2454
2455
2456
2457template <int dim, int spacedim>
2459 const hp::MappingCollection<dim, spacedim> &mapping_collection,
2460 const hp::FECollection<dim, spacedim> & fe_collection,
2461 const hp::QCollection<dim - 1> & quadrature_collection,
2462 const UpdateFlags update_flags)
2463 : n_quadrature_points(quadrature_collection.max_n_quadrature_points())
2464 , fe_face_values(nullptr)
2465 , fe_face_values_neighbor(nullptr)
2466 , internal_hp_fe_face_values(
2467 std::make_unique<hp::FEFaceValues<dim, spacedim>>(mapping_collection,
2468 fe_collection,
2469 quadrature_collection,
2470 update_flags))
2471 , internal_hp_fe_subface_values(
2472 std::make_unique<hp::FESubfaceValues<dim, spacedim>>(
2473 mapping_collection,
2474 fe_collection,
2475 quadrature_collection,
2476 update_flags))
2477 , internal_hp_fe_face_values_neighbor(
2478 std::make_unique<hp::FEFaceValues<dim, spacedim>>(mapping_collection,
2479 fe_collection,
2480 quadrature_collection,
2481 update_flags))
2482 , internal_hp_fe_subface_values_neighbor(
2483 std::make_unique<hp::FESubfaceValues<dim, spacedim>>(
2484 mapping_collection,
2485 fe_collection,
2486 quadrature_collection,
2487 update_flags))
2488{
2489 AssertDimension(dim, spacedim);
2490}
2491
2492
2493
2494template <int dim, int spacedim>
2496 const hp::FECollection<dim, spacedim> &fe_collection,
2497 const hp::QCollection<dim - 1> & quadrature_collection,
2498 const UpdateFlags update_flags)
2499 : FEInterfaceValues(fe_collection.get_reference_cell_default_linear_mapping(),
2500 fe_collection,
2501 quadrature_collection,
2502 update_flags)
2503{}
2504
2505
2506
2507template <int dim, int spacedim>
2508template <class CellIteratorType, class CellNeighborIteratorType>
2509void
2511 const CellIteratorType & cell,
2512 const unsigned int face_no,
2513 const unsigned int sub_face_no,
2514 const CellNeighborIteratorType &cell_neighbor,
2515 const unsigned int face_no_neighbor,
2516 const unsigned int sub_face_no_neighbor,
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 if (sub_face_no == numbers::invalid_unsigned_int)
2527 {
2528 internal_fe_face_values->reinit(cell, face_no);
2529 fe_face_values = internal_fe_face_values.get();
2530 }
2531 else
2532 {
2533 internal_fe_subface_values->reinit(cell, face_no, sub_face_no);
2534 fe_face_values = internal_fe_subface_values.get();
2535 }
2536 if (sub_face_no_neighbor == numbers::invalid_unsigned_int)
2537 {
2538 internal_fe_face_values_neighbor->reinit(cell_neighbor,
2539 face_no_neighbor);
2540 fe_face_values_neighbor = internal_fe_face_values_neighbor.get();
2541 }
2542 else
2543 {
2544 internal_fe_subface_values_neighbor->reinit(cell_neighbor,
2545 face_no_neighbor,
2546 sub_face_no_neighbor);
2547 fe_face_values_neighbor = internal_fe_subface_values_neighbor.get();
2548 }
2549
2550 AssertDimension(fe_face_values->n_quadrature_points,
2551 fe_face_values_neighbor->n_quadrature_points);
2552
2553 const_cast<unsigned int &>(this->n_quadrature_points) =
2554 fe_face_values->n_quadrature_points;
2555 }
2556 else if (internal_hp_fe_face_values)
2557 {
2558 unsigned int used_q_index = q_index;
2559 unsigned int used_mapping_index = mapping_index;
2560
2561 // First check. If there is only one element in a collection, and if none
2562 // had been specified explicitly, then that's clearly the one to take:
2563 if (used_q_index == numbers::invalid_unsigned_int)
2564 if (internal_hp_fe_face_values->get_quadrature_collection().size() == 1)
2565 used_q_index = 0;
2566
2567 if (used_mapping_index == numbers::invalid_unsigned_int)
2568 if (internal_hp_fe_face_values->get_mapping_collection().size() == 1)
2569 used_mapping_index = 0;
2570
2571 // Second check: See if the two quadrature objects are the same, because
2572 // in that case it does not matter which one we use. Unfortunately, we
2573 // currently have no way of testing that two mapping objects are the
2574 // same :-(
2575 if (used_q_index == numbers::invalid_unsigned_int)
2576 if (internal_hp_fe_face_values
2577 ->get_quadrature_collection()[cell->active_fe_index()] ==
2578 internal_hp_fe_face_values
2579 ->get_quadrature_collection()[cell_neighbor->active_fe_index()])
2580 used_q_index = cell->active_fe_index();
2581
2582 // Third check, if the above did not already suffice. We see if we
2583 // can get somewhere via the dominated's finite element index.
2584 const unsigned int dominated_fe_index =
2585 ((used_q_index == numbers::invalid_unsigned_int) ||
2586 (used_mapping_index == numbers::invalid_unsigned_int) ?
2587 internal_hp_fe_face_values->get_fe_collection().find_dominated_fe(
2588 {cell->active_fe_index(), cell_neighbor->active_fe_index()}) :
2589 numbers::invalid_unsigned_int);
2590
2591 if (used_q_index == numbers::invalid_unsigned_int)
2592 {
2593 Assert(dominated_fe_index != numbers::invalid_fe_index,
2594 ExcMessage(
2595 "You called this function with 'q_index' left at its "
2596 "default value, but this can only work if one of "
2597 "the two finite elements adjacent to this face "
2598 "dominates the other. See the documentation "
2599 "of this function for more information of how "
2600 "to deal with this situation."));
2601 used_q_index = dominated_fe_index;
2602 }
2603
2604 if (used_mapping_index == numbers::invalid_unsigned_int)
2605 {
2606 Assert(dominated_fe_index != numbers::invalid_fe_index,
2607 ExcMessage(
2608 "You called this function with 'mapping_index' left "
2609 "at its default value, but this can only work if one "
2610 "of the two finite elements adjacent to this face "
2611 "dominates the other. See the documentation "
2612 "of this function for more information of how "
2613 "to deal with this situation."));
2614 used_mapping_index = dominated_fe_index;
2615 }
2616
2617 // Same as if above, but when hp is enabled.
2618 if (sub_face_no == numbers::invalid_unsigned_int)
2619 {
2620 internal_hp_fe_face_values->reinit(
2621 cell, face_no, used_q_index, used_mapping_index, fe_index);
2622 fe_face_values = &const_cast<FEFaceValues<dim, spacedim> &>(
2623 internal_hp_fe_face_values->get_present_fe_values());
2624 }
2625 else
2626 {
2627 internal_hp_fe_subface_values->reinit(
2628 cell, face_no, sub_face_no, used_q_index, used_mapping_index);
2629
2630 fe_face_values = &const_cast<FESubfaceValues<dim, spacedim> &>(
2631 internal_hp_fe_subface_values->get_present_fe_values());
2632 }
2633 if (sub_face_no_neighbor == numbers::invalid_unsigned_int)
2634 {
2635 internal_hp_fe_face_values_neighbor->reinit(cell_neighbor,
2636 face_no_neighbor,
2637 used_q_index,
2638 used_mapping_index);
2639
2640 fe_face_values_neighbor = &const_cast<FEFaceValues<dim, spacedim> &>(
2641 internal_hp_fe_face_values_neighbor->get_present_fe_values());
2642 }
2643 else
2644 {
2645 internal_hp_fe_subface_values_neighbor->reinit(cell_neighbor,
2646 face_no_neighbor,
2647 sub_face_no_neighbor,
2648 used_q_index,
2649 used_mapping_index);
2650
2651 fe_face_values_neighbor =
2652 &const_cast<FESubfaceValues<dim, spacedim> &>(
2653 internal_hp_fe_subface_values_neighbor->get_present_fe_values());
2654 }
2655
2656 AssertDimension(fe_face_values->n_quadrature_points,
2657 fe_face_values_neighbor->n_quadrature_points);
2658
2659 const_cast<unsigned int &>(this->n_quadrature_points) =
2660 fe_face_values->n_quadrature_points;
2661 }
2662
2663 // Set up dof mapping and remove duplicates (for continuous elements).
2664 {
2665 // Get dof indices first:
2666 std::vector<types::global_dof_index> v(
2667 fe_face_values->get_fe().n_dofs_per_cell());
2668 cell->get_active_or_mg_dof_indices(v);
2669 std::vector<types::global_dof_index> v2(
2670 fe_face_values_neighbor->get_fe().n_dofs_per_cell());
2671 cell_neighbor->get_active_or_mg_dof_indices(v2);
2672
2673 // Fill a map from the global dof index to the left and right
2674 // local index.
2675 std::map<types::global_dof_index, std::pair<unsigned int, unsigned int>>
2676 tempmap;
2677 std::pair<unsigned int, unsigned int> invalid_entry(
2679
2680 for (unsigned int i = 0; i < v.size(); ++i)
2681 {
2682 // If not already existing, add an invalid entry:
2683 auto result = tempmap.insert(std::make_pair(v[i], invalid_entry));
2684 result.first->second.first = i;
2685 }
2686
2687 for (unsigned int i = 0; i < v2.size(); ++i)
2688 {
2689 // If not already existing, add an invalid entry:
2690 auto result = tempmap.insert(std::make_pair(v2[i], invalid_entry));
2691 result.first->second.second = i;
2692 }
2693
2694 // Transfer from the map to the sorted std::vectors.
2695 dofmap.resize(tempmap.size());
2696 interface_dof_indices.resize(tempmap.size());
2697 unsigned int idx = 0;
2698 for (auto &x : tempmap)
2699 {
2700 interface_dof_indices[idx] = x.first;
2701 dofmap[idx] = {{x.second.first, x.second.second}};
2702 ++idx;
2703 }
2704 }
2705}
2706
2707
2708
2709template <int dim, int spacedim>
2710template <class CellIteratorType>
2711void
2712FEInterfaceValues<dim, spacedim>::reinit(const CellIteratorType &cell,
2713 const unsigned int face_no,
2714 const unsigned int q_index,
2715 const unsigned int mapping_index,
2716 const unsigned int fe_index)
2717{
2718 Assert(internal_fe_face_values || internal_hp_fe_face_values,
2720
2721 if (internal_fe_face_values)
2722 {
2723 internal_fe_face_values->reinit(cell, face_no);
2724 fe_face_values = internal_fe_face_values.get();
2725 fe_face_values_neighbor = nullptr;
2726
2727 interface_dof_indices.resize(fe_face_values->get_fe().n_dofs_per_cell());
2728 cell->get_active_or_mg_dof_indices(interface_dof_indices);
2729 }
2730 else if (internal_hp_fe_face_values)
2731 {
2732 internal_hp_fe_face_values->reinit(
2733 cell, face_no, q_index, mapping_index, fe_index);
2734 fe_face_values = &const_cast<FEFaceValues<dim> &>(
2735 internal_hp_fe_face_values->get_present_fe_values());
2736 fe_face_values_neighbor = nullptr;
2737
2738 interface_dof_indices.resize(fe_face_values->get_fe().n_dofs_per_cell());
2739 cell->get_active_or_mg_dof_indices(interface_dof_indices);
2740 }
2741
2742 dofmap.resize(interface_dof_indices.size());
2743 for (unsigned int i = 0; i < interface_dof_indices.size(); ++i)
2744 {
2745 dofmap[i] = {{i, numbers::invalid_unsigned_int}};
2746 }
2747}
2748
2749
2750
2751template <int dim, int spacedim>
2752inline double
2753FEInterfaceValues<dim, spacedim>::JxW(const unsigned int q) const
2754{
2755 Assert(fe_face_values != nullptr,
2756 ExcMessage("This call requires a call to reinit() first."));
2757 return fe_face_values->JxW(q);
2758}
2759
2760
2761
2762template <int dim, int spacedim>
2763const std::vector<double> &
2765{
2766 Assert(fe_face_values != nullptr,
2767 ExcMessage("This call requires a call to reinit() first."));
2768 return fe_face_values->get_JxW_values();
2769}
2770
2771
2772
2773template <int dim, int spacedim>
2774const std::vector<Tensor<1, spacedim>> &
2776{
2777 Assert(fe_face_values != nullptr,
2778 ExcMessage("This call requires a call to reinit() first."));
2779 return fe_face_values->get_normal_vectors();
2780}
2781
2782
2783
2784template <int dim, int spacedim>
2787{
2788 Assert(!has_hp_capabilities(), ExcOnlyAvailableWithoutHP());
2789 return internal_fe_face_values->get_mapping();
2790}
2791
2792
2793
2794template <int dim, int spacedim>
2797{
2798 Assert(!has_hp_capabilities(), ExcOnlyAvailableWithoutHP());
2799 return internal_fe_face_values->get_fe();
2800}
2801
2802
2803
2804template <int dim, int spacedim>
2805const Quadrature<dim - 1> &
2807{
2808 Assert(!has_hp_capabilities(), ExcOnlyAvailableWithoutHP());
2809 return internal_fe_face_values->get_quadrature();
2810}
2811
2812
2813
2814template <int dim, int spacedim>
2817{
2818 Assert(has_hp_capabilities(), ExcOnlyAvailableWithHP());
2819 return internal_hp_fe_face_values->get_mapping_collection();
2820}
2821
2822
2823
2824template <int dim, int spacedim>
2827{
2828 Assert(has_hp_capabilities(), ExcOnlyAvailableWithHP());
2829 return internal_hp_fe_face_values->get_fe_collection();
2830}
2831
2832
2833
2834template <int dim, int spacedim>
2835const hp::QCollection<dim - 1> &
2837{
2838 Assert(has_hp_capabilities(), ExcOnlyAvailableWithHP());
2839 return internal_hp_fe_face_values->get_quadrature_collection();
2840}
2841
2842
2843
2844template <int dim, int spacedim>
2845bool
2847{
2848 if (internal_hp_fe_face_values || internal_hp_fe_subface_values ||
2849 internal_hp_fe_face_values_neighbor ||
2850 internal_hp_fe_subface_values_neighbor)
2851 {
2852 Assert(!internal_fe_face_values, ExcInternalError());
2853 Assert(!internal_fe_subface_values, ExcInternalError());
2854 Assert(!internal_fe_face_values_neighbor, ExcInternalError());
2855 Assert(!internal_fe_subface_values_neighbor, ExcInternalError());
2856
2857 return true;
2858 }
2859
2860 Assert(internal_fe_face_values || internal_fe_subface_values ||
2861 internal_fe_face_values_neighbor ||
2862 internal_fe_subface_values_neighbor,
2864 Assert(!internal_hp_fe_face_values, ExcInternalError());
2865 Assert(!internal_hp_fe_subface_values, ExcInternalError());
2866 Assert(!internal_hp_fe_face_values_neighbor, ExcInternalError());
2867 Assert(!internal_hp_fe_subface_values_neighbor, ExcInternalError());
2868
2869 return false;
2870}
2871
2872
2873
2874template <int dim, int spacedim>
2877{
2878 return {0U, n_quadrature_points};
2879}
2880
2881
2882
2883template <int dim, int spacedim>
2884const std::vector<Point<spacedim>> &
2886{
2887 Assert(fe_face_values != nullptr,
2888 ExcMessage("This call requires a call to reinit() first."));
2889 return fe_face_values->get_quadrature_points();
2890}
2891
2892
2893
2894template <int dim, int spacedim>
2897{
2898 return internal_fe_face_values->get_update_flags();
2899}
2900
2901
2902
2903template <int dim, int spacedim>
2906{
2907 return get_fe_face_values(cell_index).get_cell();
2908}
2909
2910
2911
2912template <int dim, int spacedim>
2913inline unsigned int
2915 const unsigned int cell_index) const
2916{
2917 return get_fe_face_values(cell_index).get_face_number();
2918}
2919
2920
2921
2922template <int dim, int spacedim>
2923unsigned
2925{
2926 Assert(
2927 interface_dof_indices.size() > 0,
2928 ExcMessage(
2929 "n_current_interface_dofs() is only available after a call to reinit()."));
2930 return interface_dof_indices.size();
2931}
2932
2933
2934
2935template <int dim, int spacedim>
2938{
2939 return {0U, n_current_interface_dofs()};
2940}
2941
2942
2943
2944template <int dim, int spacedim>
2945bool
2947{
2948 return fe_face_values_neighbor == nullptr;
2949}
2950
2951
2952
2953template <int dim, int spacedim>
2954std::vector<types::global_dof_index>
2956{
2957 return interface_dof_indices;
2958}
2959
2960
2961
2962template <int dim, int spacedim>
2963std::array<unsigned int, 2>
2965 const unsigned int interface_dof_index) const
2966{
2967 AssertIndexRange(interface_dof_index, n_current_interface_dofs());
2968 return dofmap[interface_dof_index];
2969}
2970
2971
2972
2973template <int dim, int spacedim>
2976 const unsigned int cell_index) const
2977{
2979 Assert(
2980 cell_index == 0 || !at_boundary(),
2981 ExcMessage(
2982 "You are on a boundary, so you can only ask for the first FEFaceValues object."));
2983
2984 return (cell_index == 0) ? *fe_face_values : *fe_face_values_neighbor;
2985}
2986
2987
2988
2989template <int dim, int spacedim>
2991FEInterfaceValues<dim, spacedim>::normal(const unsigned int q_point_index) const
2992{
2993 return fe_face_values->normal_vector(q_point_index);
2994}
2995
2996
2997
2998template <int dim, int spacedim>
2999double
3001 const bool here_or_there,
3002 const unsigned int interface_dof_index,
3003 const unsigned int q_point,
3004 const unsigned int component) const
3005{
3006 const auto dof_pair = dofmap[interface_dof_index];
3007
3008 if (here_or_there && dof_pair[0] != numbers::invalid_unsigned_int)
3009 return get_fe_face_values(0).shape_value_component(dof_pair[0],
3010 q_point,
3011 component);
3012 if (!here_or_there && dof_pair[1] != numbers::invalid_unsigned_int)
3013 return get_fe_face_values(1).shape_value_component(dof_pair[1],
3014 q_point,
3015 component);
3016
3017 return 0.0;
3018}
3019
3020
3021
3022template <int dim, int spacedim>
3023double
3025 const unsigned int interface_dof_index,
3026 const unsigned int q_point,
3027 const unsigned int component) const
3028{
3029 const auto dof_pair = dofmap[interface_dof_index];
3030
3031 double value = 0.0;
3032
3033 if (dof_pair[0] != numbers::invalid_unsigned_int)
3034 value += get_fe_face_values(0).shape_value_component(dof_pair[0],
3035 q_point,
3036 component);
3037 if (dof_pair[1] != numbers::invalid_unsigned_int)
3038 value -= get_fe_face_values(1).shape_value_component(dof_pair[1],
3039 q_point,
3040 component);
3041 return value;
3042}
3043
3044
3045
3046template <int dim, int spacedim>
3047double
3048FEInterfaceValues<dim, spacedim>::jump(const unsigned int interface_dof_index,
3049 const unsigned int q_point,
3050 const unsigned int component) const
3051{
3052 return jump_in_shape_values(interface_dof_index, q_point, component);
3053}
3054
3055
3056
3057template <int dim, int spacedim>
3058double
3060 const unsigned int interface_dof_index,
3061 const unsigned int q_point,
3062 const unsigned int component) const
3063{
3064 const auto dof_pair = dofmap[interface_dof_index];
3065
3066 if (at_boundary())
3067 return get_fe_face_values(0).shape_value_component(dof_pair[0],
3068 q_point,
3069 component);
3070
3071 double value = 0.0;
3072
3073 if (dof_pair[0] != numbers::invalid_unsigned_int)
3074 value += 0.5 * get_fe_face_values(0).shape_value_component(dof_pair[0],
3075 q_point,
3076 component);
3077 if (dof_pair[1] != numbers::invalid_unsigned_int)
3078 value += 0.5 * get_fe_face_values(1).shape_value_component(dof_pair[1],
3079 q_point,
3080 component);
3081
3082 return value;
3083}
3084
3085
3086
3087template <int dim, int spacedim>
3088double
3090 const unsigned int interface_dof_index,
3091 const unsigned int q_point,
3092 const unsigned int component) const
3093{
3094 return average_of_shape_values(interface_dof_index, q_point, component);
3095}
3096
3097
3098
3099template <int dim, int spacedim>
3102 const unsigned int interface_dof_index,
3103 const unsigned int q_point,
3104 const unsigned int component) const
3105{
3106 const auto dof_pair = dofmap[interface_dof_index];
3107
3108 if (at_boundary())
3109 return get_fe_face_values(0).shape_grad_component(dof_pair[0],
3110 q_point,
3111 component);
3112
3114
3115 if (dof_pair[0] != numbers::invalid_unsigned_int)
3116 value += 0.5 * get_fe_face_values(0).shape_grad_component(dof_pair[0],
3117 q_point,
3118 component);
3119 if (dof_pair[1] != numbers::invalid_unsigned_int)
3120 value += 0.5 * get_fe_face_values(1).shape_grad_component(dof_pair[1],
3121 q_point,
3122 component);
3123
3124 return value;
3125}
3126
3127
3128
3129template <int dim, int spacedim>
3132 const unsigned int interface_dof_index,
3133 const unsigned int q_point,
3134 const unsigned int component) const
3135{
3136 return average_of_shape_gradients(interface_dof_index, q_point, component);
3137}
3138
3139
3140
3141template <int dim, int spacedim>
3144 const unsigned int interface_dof_index,
3145 const unsigned int q_point,
3146 const unsigned int component) const
3147{
3148 const auto dof_pair = dofmap[interface_dof_index];
3149
3150 if (at_boundary())
3151 return get_fe_face_values(0).shape_hessian_component(dof_pair[0],
3152 q_point,
3153 component);
3154
3156
3157 if (dof_pair[0] != numbers::invalid_unsigned_int)
3158 value += 0.5 * get_fe_face_values(0).shape_hessian_component(dof_pair[0],
3159 q_point,
3160 component);
3161 if (dof_pair[1] != numbers::invalid_unsigned_int)
3162 value += 0.5 * get_fe_face_values(1).shape_hessian_component(dof_pair[1],
3163 q_point,
3164 component);
3165
3166 return value;
3167}
3168
3169
3170
3171template <int dim, int spacedim>
3174 const unsigned int interface_dof_index,
3175 const unsigned int q_point,
3176 const unsigned int component) const
3177{
3178 return average_of_shape_hessians(interface_dof_index, q_point, component);
3179}
3180
3181
3182
3183template <int dim, int spacedim>
3186 const unsigned int interface_dof_index,
3187 const unsigned int q_point,
3188 const unsigned int component) const
3189{
3190 const auto dof_pair = dofmap[interface_dof_index];
3191
3192 if (at_boundary())
3193 return get_fe_face_values(0).shape_grad_component(dof_pair[0],
3194 q_point,
3195 component);
3196
3198
3199 if (dof_pair[0] != numbers::invalid_unsigned_int)
3200 value += get_fe_face_values(0).shape_grad_component(dof_pair[0],
3201 q_point,
3202 component);
3203 if (dof_pair[1] != numbers::invalid_unsigned_int)
3204 value -= get_fe_face_values(1).shape_grad_component(dof_pair[1],
3205 q_point,
3206 component);
3207
3208 return value;
3209}
3210
3211
3212
3213template <int dim, int spacedim>
3216 const unsigned int interface_dof_index,
3217 const unsigned int q_point,
3218 const unsigned int component) const
3219{
3220 return jump_in_shape_gradients(interface_dof_index, q_point, component);
3221}
3222
3223
3224
3225template <int dim, int spacedim>
3228 const unsigned int interface_dof_index,
3229 const unsigned int q_point,
3230 const unsigned int component) const
3231{
3232 const auto dof_pair = dofmap[interface_dof_index];
3233
3234 if (at_boundary())
3235 return get_fe_face_values(0).shape_hessian_component(dof_pair[0],
3236 q_point,
3237 component);
3238
3240
3241 if (dof_pair[0] != numbers::invalid_unsigned_int)
3242 value += get_fe_face_values(0).shape_hessian_component(dof_pair[0],
3243 q_point,
3244 component);
3245 if (dof_pair[1] != numbers::invalid_unsigned_int)
3246 value -= get_fe_face_values(1).shape_hessian_component(dof_pair[1],
3247 q_point,
3248 component);
3249
3250 return value;
3251}
3252
3253
3254
3255template <int dim, int spacedim>
3258 const unsigned int interface_dof_index,
3259 const unsigned int q_point,
3260 const unsigned int component) const
3261{
3262 return jump_in_shape_hessians(interface_dof_index, q_point, component);
3263}
3264
3265
3266
3267template <int dim, int spacedim>
3270 const unsigned int interface_dof_index,
3271 const unsigned int q_point,
3272 const unsigned int component) const
3273{
3274 const auto dof_pair = dofmap[interface_dof_index];
3275
3276 if (at_boundary())
3277 return get_fe_face_values(0).shape_3rd_derivative_component(dof_pair[0],
3278 q_point,
3279 component);
3280
3282
3283 if (dof_pair[0] != numbers::invalid_unsigned_int)
3284 value += get_fe_face_values(0).shape_3rd_derivative_component(dof_pair[0],
3285 q_point,
3286 component);
3287 if (dof_pair[1] != numbers::invalid_unsigned_int)
3288 value -= get_fe_face_values(1).shape_3rd_derivative_component(dof_pair[1],
3289 q_point,
3290 component);
3291
3292 return value;
3293}
3294
3295
3296
3297template <int dim, int spacedim>
3300 const unsigned int interface_dof_index,
3301 const unsigned int q_point,
3302 const unsigned int component) const
3303{
3304 return jump_in_shape_3rd_derivatives(interface_dof_index, q_point, component);
3305}
3306
3307
3308
3309template <int dim, int spacedim>
3310template <class InputVector>
3311void
3313 const InputVector & fe_function,
3314 std::vector<typename InputVector::value_type> &values) const
3315{
3316 AssertDimension(values.size(), n_quadrature_points);
3317
3319 this->operator[](scalar).get_jump_in_function_values(fe_function, values);
3320}
3321
3322
3323
3324template <int dim, int spacedim>
3325template <class InputVector>
3326void
3328 const InputVector &fe_function,
3330 const
3331{
3332 AssertDimension(gradients.size(), n_quadrature_points);
3333
3335 this->operator[](scalar).get_jump_in_function_gradients(fe_function,
3336 gradients);
3337}
3338
3339
3340
3341template <int dim, int spacedim>
3342template <class InputVector>
3343void
3345 const InputVector &fe_function,
3347 const
3348{
3349 AssertDimension(hessians.size(), n_quadrature_points);
3350
3352 this->operator[](scalar).get_jump_in_function_hessians(fe_function, hessians);
3353}
3354
3355
3356
3357template <int dim, int spacedim>
3358template <class InputVector>
3359void
3361 const InputVector &fe_function,
3363 &third_derivatives) const
3364{
3365 AssertDimension(third_derivatives.size(), n_quadrature_points);
3366
3368 this->operator[](scalar).get_jump_in_function_third_derivatives(
3369 fe_function, third_derivatives);
3370}
3371
3372
3373
3374template <int dim, int spacedim>
3375template <class InputVector>
3376void
3378 const InputVector & fe_function,
3379 std::vector<typename InputVector::value_type> &values) const
3380{
3381 AssertDimension(values.size(), n_quadrature_points);
3382
3384 this->operator[](scalar).get_average_of_function_values(fe_function, values);
3385}
3386
3387
3388
3389template <int dim, int spacedim>
3390template <class InputVector>
3391void
3393 const InputVector &fe_function,
3395 const
3396{
3397 AssertDimension(gradients.size(), n_quadrature_points);
3398
3400 this->operator[](scalar).get_average_of_function_gradients(fe_function,
3401 gradients);
3402}
3403
3404
3405
3406template <int dim, int spacedim>
3407template <class InputVector>
3408void
3410 const InputVector &fe_function,
3412 const
3413{
3414 AssertDimension(hessians.size(), n_quadrature_points);
3415
3417 this->operator[](scalar).get_average_of_function_hessians(fe_function,
3418 hessians);
3419}
3420
3421
3422
3423/*------------ Inline functions: FEInterfaceValues------------*/
3424template <int dim, int spacedim>
3427 const FEValuesExtractors::Scalar &scalar) const
3428{
3429 const unsigned int n_components =
3430 (this->has_hp_capabilities() ? this->get_fe_collection().n_components() :
3431 this->get_fe().n_components());
3432 (void)n_components;
3433 AssertIndexRange(scalar.component, n_components);
3434 return FEInterfaceViews::Scalar<dim, spacedim>(*this, scalar.component);
3435}
3436
3437
3438
3439template <int dim, int spacedim>
3442 const FEValuesExtractors::Vector &vector) const
3443{
3444 const unsigned int n_components =
3445 (this->has_hp_capabilities() ? this->get_fe_collection().n_components() :
3446 this->get_fe().n_components());
3447 const unsigned int n_vectors =
3450 0);
3451 (void)n_components;
3452 (void)n_vectors;
3453 AssertIndexRange(vector.first_vector_component, n_vectors);
3455 vector.first_vector_component);
3456}
3457
3458
3459
3460namespace FEInterfaceViews
3461{
3462 template <int dim, int spacedim>
3464 const FEInterfaceValues<dim, spacedim> &fe_interface)
3465 : fe_interface(&fe_interface)
3466 {}
3467
3468
3469
3470 template <int dim, int spacedim>
3471 Scalar<dim, spacedim>::Scalar(
3472 const FEInterfaceValues<dim, spacedim> &fe_interface,
3473 const unsigned int component)
3474 : Base<dim, spacedim>(fe_interface)
3475 , extractor(component)
3476 {}
3477
3478
3479
3480 template <int dim, int spacedim>
3481 template <class InputVector, class OutputVector>
3482 void
3483 Base<dim, spacedim>::get_local_dof_values(
3484 const InputVector &dof_values,
3485 OutputVector & local_dof_values) const
3486 {
3487 const auto &interface_dof_indices =
3488 this->fe_interface->get_interface_dof_indices();
3489
3490 AssertDimension(interface_dof_indices.size(), local_dof_values.size());
3491
3492 for (const unsigned int i : this->fe_interface->dof_indices())
3493 local_dof_values[i] = dof_values(interface_dof_indices[i]);
3494 }
3495
3496
3497
3498 template <int dim, int spacedim>
3499 typename Scalar<dim, spacedim>::value_type
3500 Scalar<dim, spacedim>::value(const bool here_or_there,
3501 const unsigned int interface_dof_index,
3502 const unsigned int q_point) const
3503 {
3504 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3505
3506 if (here_or_there && dof_pair[0] != numbers::invalid_unsigned_int)
3507 return (*(this->fe_interface->fe_face_values))[extractor].value(
3508 dof_pair[0], q_point);
3509
3510 if (!here_or_there && dof_pair[1] != numbers::invalid_unsigned_int)
3511 return (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3512 dof_pair[1], q_point);
3513
3514 return 0.0;
3515 }
3516
3517
3518
3519 template <int dim, int spacedim>
3520 typename Scalar<dim, spacedim>::value_type
3521 Scalar<dim, spacedim>::jump_in_values(const unsigned int interface_dof_index,
3522 const unsigned int q_point) const
3523 {
3524 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3525
3526 value_type value = 0.0;
3527
3528 if (dof_pair[0] != numbers::invalid_unsigned_int)
3529 value +=
3530 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3531 q_point);
3532
3533 if (dof_pair[1] != numbers::invalid_unsigned_int)
3534 value -=
3535 (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3536 dof_pair[1], q_point);
3537
3538 return value;
3539 }
3540
3541
3542
3543 template <int dim, int spacedim>
3544 typename Scalar<dim, spacedim>::value_type
3545 Scalar<dim, spacedim>::jump(const unsigned int interface_dof_index,
3546 const unsigned int q_point) const
3547 {
3548 return jump_in_values(interface_dof_index, q_point);
3549 }
3550
3551
3552
3553 template <int dim, int spacedim>
3554 typename Scalar<dim, spacedim>::value_type
3555 Scalar<dim, spacedim>::average_of_values(
3556 const unsigned int interface_dof_index,
3557 const unsigned int q_point) const
3558 {
3559 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3560
3561 if (this->fe_interface->at_boundary())
3562 return (*(this->fe_interface->fe_face_values))[extractor].value(
3563 dof_pair[0], q_point);
3564
3565 value_type value = 0.0;
3566
3567 if (dof_pair[0] != numbers::invalid_unsigned_int)
3568 value +=
3569 0.5 *
3570 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
3571 q_point);
3572
3573 if (dof_pair[1] != numbers::invalid_unsigned_int)
3574 value +=
3575 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
3576 dof_pair[1], q_point);
3577
3578 return value;
3579 }
3580
3581
3582
3583 template <int dim, int spacedim>
3584 typename Scalar<dim, spacedim>::value_type
3585 Scalar<dim, spacedim>::average(const unsigned int interface_dof_index,
3586 const unsigned int q_point) const
3587 {
3588 return average_of_values(interface_dof_index, q_point);
3589 }
3590
3591
3592
3593 template <int dim, int spacedim>
3594 typename Scalar<dim, spacedim>::gradient_type
3595 Scalar<dim, spacedim>::average_of_gradients(
3596 const unsigned int interface_dof_index,
3597 const unsigned int q_point) const
3598 {
3599 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3600
3601 if (this->fe_interface->at_boundary())
3602 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
3603 dof_pair[0], q_point);
3604
3605 gradient_type value;
3606
3607 if (dof_pair[0] != numbers::invalid_unsigned_int)
3608 value +=
3609 0.5 *
3610 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
3611 q_point);
3612
3613 if (dof_pair[1] != numbers::invalid_unsigned_int)
3614 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3615 .gradient(dof_pair[1], q_point);
3616
3617 return value;
3618 }
3619
3620
3621 template <int dim, int spacedim>
3622 typename Scalar<dim, spacedim>::gradient_type
3623 Scalar<dim, spacedim>::average_gradient(
3624 const unsigned int interface_dof_index,
3625 const unsigned int q_point) const
3626 {
3627 return average_of_gradients(interface_dof_index, q_point);
3628 }
3629
3630
3631
3632 template <int dim, int spacedim>
3633 typename Scalar<dim, spacedim>::gradient_type
3634 Scalar<dim, spacedim>::jump_in_gradients(
3635 const unsigned int interface_dof_index,
3636 const unsigned int q_point) const
3637 {
3638 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3639
3640 if (this->fe_interface->at_boundary())
3641 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
3642 dof_pair[0], q_point);
3643
3644 gradient_type value;
3645
3646 if (dof_pair[0] != numbers::invalid_unsigned_int)
3647 value +=
3648 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
3649 q_point);
3650
3651 if (dof_pair[1] != numbers::invalid_unsigned_int)
3652 value -=
3653 (*(this->fe_interface->fe_face_values_neighbor))[extractor].gradient(
3654 dof_pair[1], q_point);
3655
3656 return value;
3657 }
3658
3659
3660
3661 template <int dim, int spacedim>
3662 typename Scalar<dim, spacedim>::gradient_type
3663 Scalar<dim, spacedim>::jump_gradient(const unsigned int interface_dof_index,
3664 const unsigned int q_point) const
3665 {
3666 return jump_in_gradients(interface_dof_index, q_point);
3667 }
3668
3669
3670
3671 template <int dim, int spacedim>
3672 typename Scalar<dim, spacedim>::hessian_type
3673 Scalar<dim, spacedim>::average_of_hessians(
3674 const unsigned int interface_dof_index,
3675 const unsigned int q_point) const
3676 {
3677 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3678
3679 if (this->fe_interface->at_boundary())
3680 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
3681 dof_pair[0], q_point);
3682
3683 hessian_type value;
3684
3685 if (dof_pair[0] != numbers::invalid_unsigned_int)
3686 value +=
3687 0.5 *
3688 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
3689 q_point);
3690
3691 if (dof_pair[1] != numbers::invalid_unsigned_int)
3692 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3693 .hessian(dof_pair[1], q_point);
3694
3695 return value;
3696 }
3697
3698
3699
3700 template <int dim, int spacedim>
3701 typename Scalar<dim, spacedim>::hessian_type
3702 Scalar<dim, spacedim>::average_hessian(const unsigned int interface_dof_index,
3703 const unsigned int q_point) const
3704 {
3705 return average_of_hessians(interface_dof_index, q_point);
3706 }
3707
3708
3709
3710 template <int dim, int spacedim>
3711 typename Scalar<dim, spacedim>::third_derivative_type
3712 Scalar<dim, spacedim>::jump_in_third_derivatives(
3713 const unsigned int interface_dof_index,
3714 const unsigned int q_point) const
3715 {
3716 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3717
3718 if (this->fe_interface->at_boundary())
3719 return (*(this->fe_interface->fe_face_values))[extractor]
3720 .third_derivative(dof_pair[0], q_point);
3721
3722 third_derivative_type value;
3723
3724 if (dof_pair[0] != numbers::invalid_unsigned_int)
3725 value +=
3726 (*(this->fe_interface->fe_face_values))[extractor].third_derivative(
3727 dof_pair[0], q_point);
3728
3729 if (dof_pair[1] != numbers::invalid_unsigned_int)
3730 value -= (*(this->fe_interface->fe_face_values_neighbor))[extractor]
3731 .third_derivative(dof_pair[1], q_point);
3732
3733 return value;
3734 }
3735
3736
3737
3738 template <int dim, int spacedim>
3739 typename Scalar<dim, spacedim>::third_derivative_type
3740 Scalar<dim, spacedim>::jump_3rd_derivative(
3741 const unsigned int interface_dof_index,
3742 const unsigned int q_point) const
3743 {
3744 return jump_in_third_derivatives(interface_dof_index, q_point);
3745 }
3746
3747
3748
3749 template <int dim, int spacedim>
3750 typename Scalar<dim, spacedim>::hessian_type
3751 Scalar<dim, spacedim>::jump_in_hessians(
3752 const unsigned int interface_dof_index,
3753 const unsigned int q_point) const
3754 {
3755 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
3756
3757 if (this->fe_interface->at_boundary())
3758 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
3759 dof_pair[0], q_point);
3760
3761 hessian_type value;
3762
3763 if (dof_pair[0] != numbers::invalid_unsigned_int)
3764 value +=
3765 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
3766 q_point);
3767
3768 if (dof_pair[1] != numbers::invalid_unsigned_int)
3769 value -=
3770 (*(this->fe_interface->fe_face_values_neighbor))[extractor].hessian(
3771 dof_pair[1], q_point);
3772
3773 return value;
3774 }
3775
3776
3777
3778 template <int dim, int spacedim>
3779 typename Scalar<dim, spacedim>::hessian_type
3780 Scalar<dim, spacedim>::jump_hessian(const unsigned int interface_dof_index,
3781 const unsigned int q_point) const
3782 {
3783 return jump_in_hessians(interface_dof_index, q_point);
3784 }
3785
3786
3787
3788 template <int dim, int spacedim>
3789 template <class InputVector>
3790 void
3791 Scalar<dim, spacedim>::get_function_values_from_local_dof_values(
3792 const bool here_or_there,
3793 const InputVector &local_dof_values,
3794 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3795 const
3796 {
3797 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
3798
3799 for (const auto dof_index : this->fe_interface->dof_indices())
3800 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3801 {
3802 if (dof_index == 0)
3803 values[q_index] = 0.;
3804
3805 values[q_index] += local_dof_values[dof_index] *
3806 value(here_or_there, dof_index, q_index);
3807 }
3808 }
3809
3810
3811
3812 template <int dim, int spacedim>
3813 template <class InputVector>
3814 void
3815 Scalar<dim, spacedim>::get_function_values(
3816 const bool here_or_there,
3817 const InputVector &fe_function,
3818 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3819 const
3820 {
3821 std::vector<typename InputVector::value_type> local_dof_values(
3822 this->fe_interface->n_current_interface_dofs());
3823 this->get_local_dof_values(fe_function, local_dof_values);
3824
3825 get_function_values_from_local_dof_values(here_or_there,
3826 local_dof_values,
3827 values);
3828 }
3829
3830
3831
3832 template <int dim, int spacedim>
3833 template <class InputVector>
3834 void
3835 Scalar<dim, spacedim>::get_jump_in_function_values_from_local_dof_values(
3836 const InputVector &local_dof_values,
3837 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3838 const
3839 {
3840 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
3841
3842 for (const auto dof_index : this->fe_interface->dof_indices())
3843 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3844 {
3845 if (dof_index == 0)
3846 values[q_index] = 0.;
3847
3848 values[q_index] +=
3849 local_dof_values[dof_index] * jump_in_values(dof_index, q_index);
3850 }
3851 }
3852
3853
3854
3855 template <int dim, int spacedim>
3856 template <class InputVector>
3857 void
3858 Scalar<dim, spacedim>::get_jump_in_function_values(
3859 const InputVector &fe_function,
3860 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3861 const
3862 {
3863 std::vector<typename InputVector::value_type> local_dof_values(
3864 this->fe_interface->n_current_interface_dofs());
3865 this->get_local_dof_values(fe_function, local_dof_values);
3866
3867 get_jump_in_function_values_from_local_dof_values(local_dof_values, values);
3868 }
3869
3870
3871
3872 template <int dim, int spacedim>
3873 template <class InputVector>
3874 void
3875 Scalar<dim, spacedim>::get_jump_in_function_gradients_from_local_dof_values(
3876 const InputVector &local_dof_values,
3877 std::vector<solution_gradient_type<typename InputVector::value_type>>
3878 &gradients) const
3879 {
3880 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
3881
3882 for (const auto dof_index : this->fe_interface->dof_indices())
3883 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3884 {
3885 if (dof_index == 0)
3886 gradients[q_index] = 0.;
3887
3888 gradients[q_index] +=
3889 local_dof_values[dof_index] * jump_in_gradients(dof_index, q_index);
3890 }
3891 }
3892
3893
3894
3895 template <int dim, int spacedim>
3896 template <class InputVector>
3897 void
3898 Scalar<dim, spacedim>::get_jump_in_function_gradients(
3899 const InputVector &fe_function,
3900 std::vector<solution_gradient_type<typename InputVector::value_type>>
3901 &gradients) const
3902 {
3903 std::vector<typename InputVector::value_type> local_dof_values(
3904 this->fe_interface->n_current_interface_dofs());
3905 this->get_local_dof_values(fe_function, local_dof_values);
3906
3907 get_jump_in_function_gradients_from_local_dof_values(local_dof_values,
3908 gradients);
3909 }
3910
3911
3912
3913 template <int dim, int spacedim>
3914 template <class InputVector>
3915 void
3916 Scalar<dim, spacedim>::get_average_of_function_values_from_local_dof_values(
3917 const InputVector &local_dof_values,
3918 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3919 const
3920 {
3921 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
3922
3923 for (const auto dof_index : this->fe_interface->dof_indices())
3924 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3925 {
3926 if (dof_index == 0)
3927 values[q_index] = 0.;
3928
3929 values[q_index] +=
3930 local_dof_values[dof_index] * average_of_values(dof_index, q_index);
3931 }
3932 }
3933
3934
3935
3936 template <int dim, int spacedim>
3937 template <class InputVector>
3938 void
3939 Scalar<dim, spacedim>::get_average_of_function_values(
3940 const InputVector &fe_function,
3941 std::vector<solution_value_type<typename InputVector::value_type>> &values)
3942 const
3943 {
3944 std::vector<typename InputVector::value_type> local_dof_values(
3945 this->fe_interface->n_current_interface_dofs());
3946 this->get_local_dof_values(fe_function, local_dof_values);
3947
3948 get_average_of_function_values_from_local_dof_values(local_dof_values,
3949 values);
3950 }
3951
3952
3953
3954 template <int dim, int spacedim>
3955 template <class InputVector>
3956 void
3957 Scalar<dim, spacedim>::
3958 get_average_of_function_gradients_from_local_dof_values(
3959 const InputVector &local_dof_values,
3960 std::vector<solution_gradient_type<typename InputVector::value_type>>
3961 &gradients) const
3962 {
3963 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
3964
3965 for (const auto dof_index : this->fe_interface->dof_indices())
3966 for (const auto q_index : this->fe_interface->quadrature_point_indices())
3967 {
3968 if (dof_index == 0)
3969 gradients[q_index] = 0.;
3970
3971 gradients[q_index] += local_dof_values[dof_index] *
3972 average_of_gradients(dof_index, q_index);
3973 }
3974 }
3975
3976
3977
3978 template <int dim, int spacedim>
3979 template <class InputVector>
3980 void
3981 Scalar<dim, spacedim>::get_average_of_function_gradients(
3982 const InputVector &fe_function,
3983 std::vector<solution_gradient_type<typename InputVector::value_type>>
3984 &gradients) const
3985 {
3986 std::vector<typename InputVector::value_type> local_dof_values(
3987 this->fe_interface->n_current_interface_dofs());
3988 this->get_local_dof_values(fe_function, local_dof_values);
3989
3990 get_average_of_function_gradients_from_local_dof_values(local_dof_values,
3991 gradients);
3992 }
3993
3994
3995
3996 template <int dim, int spacedim>
3997 template <class InputVector>
3998 void
3999 Scalar<dim, spacedim>::get_jump_in_function_hessians_from_local_dof_values(
4000 const InputVector &local_dof_values,
4001 std::vector<solution_hessian_type<typename InputVector::value_type>>
4002 &hessians) const
4003 {
4004 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
4005
4006 for (const auto dof_index : this->fe_interface->dof_indices())
4007 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4008 {
4009 if (dof_index == 0)
4010 hessians[q_index] = 0.;
4011
4012 hessians[q_index] +=
4013 local_dof_values[dof_index] * jump_in_hessians(dof_index, q_index);
4014 }
4015 }
4016
4017
4018
4019 template <int dim, int spacedim>
4020 template <class InputVector>
4021 void
4022 Scalar<dim, spacedim>::get_jump_in_function_hessians(
4023 const InputVector &fe_function,
4024 std::vector<solution_hessian_type<typename InputVector::value_type>>
4025 &hessians) const
4026 {
4027 std::vector<typename InputVector::value_type> local_dof_values(
4028 this->fe_interface->n_current_interface_dofs());
4029 this->get_local_dof_values(fe_function, local_dof_values);
4030
4031 get_jump_in_function_hessians_from_local_dof_values(local_dof_values,
4032 hessians);
4033 }
4034
4035
4036
4037 template <int dim, int spacedim>
4038 template <class InputVector>
4039 void
4040 Scalar<dim, spacedim>::get_average_of_function_hessians_from_local_dof_values(
4041 const InputVector &local_dof_values,
4042 std::vector<solution_hessian_type<typename InputVector::value_type>>
4043 &hessians) const
4044 {
4045 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
4046
4047 for (const unsigned int dof_index : this->fe_interface->dof_indices())
4048 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4049 {
4050 if (dof_index == 0)
4051 hessians[q_index] = 0.;
4052
4053 hessians[q_index] += local_dof_values[dof_index] *
4054 average_of_hessians(dof_index, q_index);
4055 }
4056 }
4057
4058
4059
4060 template <int dim, int spacedim>
4061 template <class InputVector>
4062 void
4063 Scalar<dim, spacedim>::get_average_of_function_hessians(
4064 const InputVector &fe_function,
4065 std::vector<solution_hessian_type<typename InputVector::value_type>>
4066 &hessians) 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_average_of_function_hessians_from_local_dof_values(local_dof_values,
4073 hessians);
4074 }
4075
4076
4077
4078 template <int dim, int spacedim>
4079 template <class InputVector>
4080 void
4081 Scalar<dim, spacedim>::
4082 get_jump_in_function_third_derivatives_from_local_dof_values(
4083 const InputVector &local_dof_values,
4084 std::vector<
4085 solution_third_derivative_type<typename InputVector::value_type>>
4086 &third_derivatives) const
4087 {
4088 AssertDimension(third_derivatives.size(),
4089 this->fe_interface->n_quadrature_points);
4090
4091 for (const unsigned int dof_index : this->fe_interface->dof_indices())
4092 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4093 {
4094 if (dof_index == 0)
4095 third_derivatives[q_index] = 0.;
4096
4097 third_derivatives[q_index] +=
4098 local_dof_values[dof_index] *
4099 jump_in_third_derivatives(dof_index, q_index);
4100 }
4101 }
4102
4103
4104
4105 template <int dim, int spacedim>
4106 template <class InputVector>
4107 void
4108 Scalar<dim, spacedim>::get_jump_in_function_third_derivatives(
4109 const InputVector &fe_function,
4110 std::vector<
4111 solution_third_derivative_type<typename InputVector::value_type>>
4112 &third_derivatives) const
4113 {
4114 std::vector<typename InputVector::value_type> local_dof_values(
4115 this->fe_interface->n_current_interface_dofs());
4116 this->get_local_dof_values(fe_function, local_dof_values);
4117
4118 get_jump_in_function_third_derivatives_from_local_dof_values(
4119 local_dof_values, third_derivatives);
4120 }
4121
4122
4123
4124 template <int dim, int spacedim>
4126 const FEInterfaceValues<dim, spacedim> &fe_interface,
4127 const unsigned int first_vector_component)
4128 : Base<dim, spacedim>(fe_interface)
4129 , extractor(first_vector_component)
4130 {}
4131
4132
4133
4134 template <int dim, int spacedim>
4136 Vector<dim, spacedim>::value(const bool here_or_there,
4137 const unsigned int interface_dof_index,
4138 const unsigned int q_point) const
4139 {
4140 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4141
4142 if (here_or_there && dof_pair[0] != numbers::invalid_unsigned_int)
4143 return (*(this->fe_interface->fe_face_values))[extractor].value(
4144 dof_pair[0], q_point);
4145
4146 if (!here_or_there && dof_pair[1] != numbers::invalid_unsigned_int)
4147 return (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
4148 dof_pair[1], q_point);
4149
4150 return value_type();
4151 }
4152
4153
4154
4155 template <int dim, int spacedim>
4157 Vector<dim, spacedim>::jump_in_values(const unsigned int interface_dof_index,
4158 const unsigned int q_point) const
4159 {
4160 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4161
4162 value_type value;
4163
4164 if (dof_pair[0] != numbers::invalid_unsigned_int)
4165 value +=
4166 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
4167 q_point);
4168
4169 if (dof_pair[1] != numbers::invalid_unsigned_int)
4170 value -=
4171 (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
4172 dof_pair[1], q_point);
4173
4174 return value;
4175 }
4176
4177
4178
4179 template <int dim, int spacedim>
4181 Vector<dim, spacedim>::jump(const unsigned int interface_dof_index,
4182 const unsigned int q_point) const
4183 {
4184 return jump_in_values(interface_dof_index, q_point);
4185 }
4186
4187
4188
4189 template <int dim, int spacedim>
4192 const unsigned int interface_dof_index,
4193 const unsigned int q_point) const
4194 {
4195 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4196
4197 if (this->fe_interface->at_boundary())
4198 return (*(this->fe_interface->fe_face_values))[extractor].value(
4199 dof_pair[0], q_point);
4200
4201 value_type value;
4202
4203 if (dof_pair[0] != numbers::invalid_unsigned_int)
4204 value +=
4205 0.5 *
4206 (*(this->fe_interface->fe_face_values))[extractor].value(dof_pair[0],
4207 q_point);
4208
4209 if (dof_pair[1] != numbers::invalid_unsigned_int)
4210 value +=
4211 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor].value(
4212 dof_pair[1], q_point);
4213
4214 return value;
4215 }
4216
4217
4218
4219 template <int dim, int spacedim>
4221 Vector<dim, spacedim>::average(const unsigned int interface_dof_index,
4222 const unsigned int q_point) const
4223 {
4224 return average_of_values(interface_dof_index, q_point);
4225 }
4226
4227
4228
4229 template <int dim, int spacedim>
4232 const unsigned int interface_dof_index,
4233 const unsigned int q_point) const
4234 {
4235 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4236
4237 if (this->fe_interface->at_boundary())
4238 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
4239 dof_pair[0], q_point);
4240
4241 gradient_type value;
4242
4243 if (dof_pair[0] != numbers::invalid_unsigned_int)
4244 value +=
4245 0.5 *
4246 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
4247 q_point);
4248
4249 if (dof_pair[1] != numbers::invalid_unsigned_int)
4250 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
4251 .gradient(dof_pair[1], q_point);
4252
4253 return value;
4254 }
4255
4256
4257
4258 template <int dim, int spacedim>
4261 const unsigned int interface_dof_index,
4262 const unsigned int q_point) const
4263 {
4264 return average_of_gradients(interface_dof_index, q_point);
4265 }
4266
4267
4268
4269 template <int dim, int spacedim>
4272 const unsigned int interface_dof_index,
4273 const unsigned int q_point) const
4274 {
4275 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4276
4277 if (this->fe_interface->at_boundary())
4278 return (*(this->fe_interface->fe_face_values))[extractor].gradient(
4279 dof_pair[0], q_point);
4280
4281 gradient_type value;
4282
4283 if (dof_pair[0] != numbers::invalid_unsigned_int)
4284 value +=
4285 (*(this->fe_interface->fe_face_values))[extractor].gradient(dof_pair[0],
4286 q_point);
4287
4288 if (dof_pair[1] != numbers::invalid_unsigned_int)
4289 value -=
4290 (*(this->fe_interface->fe_face_values_neighbor))[extractor].gradient(
4291 dof_pair[1], q_point);
4292
4293 return value;
4294 }
4295
4296
4297
4298 template <int dim, int spacedim>
4300 Vector<dim, spacedim>::jump_gradient(const unsigned int interface_dof_index,
4301 const unsigned int q_point) const
4302 {
4303 return jump_in_gradients(interface_dof_index, q_point);
4304 }
4305
4306
4307
4308 template <int dim, int spacedim>
4311 const unsigned int interface_dof_index,
4312 const unsigned int q_point) const
4313 {
4314 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4315
4316 if (this->fe_interface->at_boundary())
4317 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
4318 dof_pair[0], q_point);
4319
4320 hessian_type value;
4321
4322 if (dof_pair[0] != numbers::invalid_unsigned_int)
4323 value +=
4324 0.5 *
4325 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
4326 q_point);
4327
4328 if (dof_pair[1] != numbers::invalid_unsigned_int)
4329 value += 0.5 * (*(this->fe_interface->fe_face_values_neighbor))[extractor]
4330 .hessian(dof_pair[1], q_point);
4331
4332 return value;
4333 }
4334
4335
4336
4337 template <int dim, int spacedim>
4339 Vector<dim, spacedim>::average_hessian(const unsigned int interface_dof_index,
4340 const unsigned int q_point) const
4341 {
4342 return average_of_hessians(interface_dof_index, q_point);
4343 }
4344
4345
4346
4347 template <int dim, int spacedim>
4350 const unsigned int interface_dof_index,
4351 const unsigned int q_point) const
4352 {
4353 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4354
4355 if (this->fe_interface->at_boundary())
4356 return (*(this->fe_interface->fe_face_values))[extractor].hessian(
4357 dof_pair[0], q_point);
4358
4359 hessian_type value;
4360
4361 if (dof_pair[0] != numbers::invalid_unsigned_int)
4362 value +=
4363 (*(this->fe_interface->fe_face_values))[extractor].hessian(dof_pair[0],
4364 q_point);
4365
4366 if (dof_pair[1] != numbers::invalid_unsigned_int)
4367 value -=
4368 (*(this->fe_interface->fe_face_values_neighbor))[extractor].hessian(
4369 dof_pair[1], q_point);
4370
4371 return value;
4372 }
4373
4374
4375
4376 template <int dim, int spacedim>
4378 Vector<dim, spacedim>::jump_hessian(const unsigned int interface_dof_index,
4379 const unsigned int q_point) const
4380 {
4381 return jump_in_hessians(interface_dof_index, q_point);
4382 }
4383
4384
4385
4386 template <int dim, int spacedim>
4389 const unsigned int interface_dof_index,
4390 const unsigned int q_point) const
4391 {
4392 const auto dof_pair = this->fe_interface->dofmap[interface_dof_index];
4393
4394 if (this->fe_interface->at_boundary())
4395 return (*(this->fe_interface->fe_face_values))[extractor]
4396 .third_derivative(dof_pair[0], q_point);
4397
4398 third_derivative_type value;
4399
4400 if (dof_pair[0] != numbers::invalid_unsigned_int)
4401 value +=
4402 (*(this->fe_interface->fe_face_values))[extractor].third_derivative(
4403 dof_pair[0], q_point);
4404
4405 if (dof_pair[1] != numbers::invalid_unsigned_int)
4406 value -= (*(this->fe_interface->fe_face_values_neighbor))[extractor]
4407 .third_derivative(dof_pair[1], q_point);
4408
4409 return value;
4410 }
4411
4412
4413
4414 template <int dim, int spacedim>
4417 const unsigned int interface_dof_index,
4418 const unsigned int q_point) const
4419 {
4420 return jump_in_third_derivatives(interface_dof_index, q_point);
4421 }
4422
4423
4424
4425 template <int dim, int spacedim>
4426 template <class InputVector>
4427 void
4429 const bool here_or_there,
4430 const InputVector &local_dof_values,
4431 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4432 const
4433 {
4434 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
4435
4436 for (const auto dof_index : this->fe_interface->dof_indices())
4437 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4438 {
4439 if (dof_index == 0)
4440 values[q_index] = 0.;
4441
4442 values[q_index] += local_dof_values[dof_index] *
4443 value(here_or_there, dof_index, q_index);
4444 }
4445 }
4446
4447
4448
4449 template <int dim, int spacedim>
4450 template <class InputVector>
4451 void
4453 const bool here_or_there,
4454 const InputVector &fe_function,
4455 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4456 const
4457 {
4458 std::vector<typename InputVector::value_type> local_dof_values(
4459 this->fe_interface->n_current_interface_dofs());
4460 this->get_local_dof_values(fe_function, local_dof_values);
4461
4462 get_function_values_from_local_dof_values(here_or_there,
4463 local_dof_values,
4464 values);
4465 }
4466
4467
4468
4469 template <int dim, int spacedim>
4470 template <class InputVector>
4471 void
4473 const InputVector &local_dof_values,
4474 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4475 const
4476 {
4477 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
4478
4479 for (const auto dof_index : this->fe_interface->dof_indices())
4480 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4481 {
4482 if (dof_index == 0)
4483 values[q_index] = 0.;
4484
4485 values[q_index] +=
4486 local_dof_values[dof_index] * jump_in_values(dof_index, q_index);
4487 }
4488 }
4489
4490
4491
4492 template <int dim, int spacedim>
4493 template <class InputVector>
4494 void
4496 const InputVector &fe_function,
4497 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4498 const
4499 {
4500 std::vector<typename InputVector::value_type> local_dof_values(
4501 this->fe_interface->n_current_interface_dofs());
4502 this->get_local_dof_values(fe_function, local_dof_values);
4503
4504 get_jump_in_function_values_from_local_dof_values(local_dof_values, values);
4505 }
4506
4507
4508
4509 template <int dim, int spacedim>
4510 template <class InputVector>
4511 void
4513 const InputVector &local_dof_values,
4514 std::vector<solution_gradient_type<typename InputVector::value_type>>
4515 &gradients) const
4516 {
4517 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
4518
4519 for (const auto dof_index : this->fe_interface->dof_indices())
4520 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4521 {
4522 if (dof_index == 0)
4523 gradients[q_index] = 0.;
4524
4525 gradients[q_index] +=
4526 local_dof_values[dof_index] * jump_in_gradients(dof_index, q_index);
4527 }
4528 }
4529
4530
4531
4532 template <int dim, int spacedim>
4533 template <class InputVector>
4534 void
4536 const InputVector &fe_function,
4537 std::vector<solution_gradient_type<typename InputVector::value_type>>
4538 &gradients) const
4539 {
4540 std::vector<typename InputVector::value_type> local_dof_values(
4541 this->fe_interface->n_current_interface_dofs());
4542 this->get_local_dof_values(fe_function, local_dof_values);
4543
4544 get_jump_in_function_gradients_from_local_dof_values(local_dof_values,
4545 gradients);
4546 }
4547
4548
4549
4550 template <int dim, int spacedim>
4551 template <class InputVector>
4552 void
4554 const InputVector &local_dof_values,
4555 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4556 const
4557 {
4558 AssertDimension(values.size(), this->fe_interface->n_quadrature_points);
4559
4560 for (const auto dof_index : this->fe_interface->dof_indices())
4561 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4562 {
4563 if (dof_index == 0)
4564 values[q_index] = 0.;
4565
4566 values[q_index] +=
4567 local_dof_values[dof_index] * average_of_values(dof_index, q_index);
4568 }
4569 }
4570
4571
4572
4573 template <int dim, int spacedim>
4574 template <class InputVector>
4575 void
4577 const InputVector &fe_function,
4578 std::vector<solution_value_type<typename InputVector::value_type>> &values)
4579 const
4580 {
4581 std::vector<typename InputVector::value_type> local_dof_values(
4582 this->fe_interface->n_current_interface_dofs());
4583 this->get_local_dof_values(fe_function, local_dof_values);
4584
4585 get_average_of_function_values_from_local_dof_values(local_dof_values,
4586 values);
4587 }
4588
4589
4590
4591 template <int dim, int spacedim>
4592 template <class InputVector>
4593 void
4596 const InputVector &local_dof_values,
4597 std::vector<solution_gradient_type<typename InputVector::value_type>>
4598 &gradients) const
4599 {
4600 AssertDimension(gradients.size(), this->fe_interface->n_quadrature_points);
4601
4602 for (const auto dof_index : this->fe_interface->dof_indices())
4603 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4604 {
4605 if (dof_index == 0)
4606 gradients[q_index] = 0.;
4607
4608 gradients[q_index] += local_dof_values[dof_index] *
4609 average_of_gradients(dof_index, q_index);
4610 }
4611 }
4612
4613
4614
4615 template <int dim, int spacedim>
4616 template <class InputVector>
4617 void
4619 const InputVector &fe_function,
4620 std::vector<solution_gradient_type<typename InputVector::value_type>>
4621 &gradients) const
4622 {
4623 std::vector<typename InputVector::value_type> local_dof_values(
4624 this->fe_interface->n_current_interface_dofs());
4625 this->get_local_dof_values(fe_function, local_dof_values);
4626
4627 get_average_of_function_gradients_from_local_dof_values(local_dof_values,
4628 gradients);
4629 }
4630
4631
4632
4633 template <int dim, int spacedim>
4634 template <class InputVector>
4635 void
4637 const InputVector &local_dof_values,
4638 std::vector<solution_hessian_type<typename InputVector::value_type>>
4639 &hessians) const
4640 {
4641 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
4642
4643 for (const auto dof_index : this->fe_interface->dof_indices())
4644 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4645 {
4646 if (dof_index == 0)
4647 hessians[q_index] = 0.;
4648
4649 hessians[q_index] +=
4650 local_dof_values[dof_index] * jump_in_hessians(dof_index, q_index);
4651 }
4652 }
4653
4654
4655
4656 template <int dim, int spacedim>
4657 template <class InputVector>
4658 void
4660 const InputVector &fe_function,
4661 std::vector<solution_hessian_type<typename InputVector::value_type>>
4662 &hessians) const
4663 {
4664 std::vector<typename InputVector::value_type> local_dof_values(
4665 this->fe_interface->n_current_interface_dofs());
4666 this->get_local_dof_values(fe_function, local_dof_values);
4667
4668 get_jump_in_function_hessians_from_local_dof_values(local_dof_values,
4669 hessians);
4670 }
4671
4672
4673
4674 template <int dim, int spacedim>
4675 template <class InputVector>
4676 void
4678 const InputVector &local_dof_values,
4679 std::vector<solution_hessian_type<typename InputVector::value_type>>
4680 &hessians) const
4681 {
4682 AssertDimension(hessians.size(), this->fe_interface->n_quadrature_points);
4683
4684 for (const unsigned int dof_index : this->fe_interface->dof_indices())
4685 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4686 {
4687 if (dof_index == 0)
4688 hessians[q_index] = 0.;
4689
4690 hessians[q_index] += local_dof_values[dof_index] *
4691 average_of_hessians(dof_index, q_index);
4692 }
4693 }
4694
4695
4696
4697 template <int dim, int spacedim>
4698 template <class InputVector>
4699 void
4701 const InputVector &fe_function,
4702 std::vector<solution_hessian_type<typename InputVector::value_type>>
4703 &hessians) const
4704 {
4705 std::vector<typename InputVector::value_type> local_dof_values(
4706 this->fe_interface->n_current_interface_dofs());
4707 this->get_local_dof_values(fe_function, local_dof_values);
4708
4709 get_average_of_function_hessians_from_local_dof_values(local_dof_values,
4710 hessians);
4711 }
4712
4713
4714
4715 template <int dim, int spacedim>
4716 template <class InputVector>
4717 void
4720 const InputVector &local_dof_values,
4721 std::vector<
4722 solution_third_derivative_type<typename InputVector::value_type>>
4723 &third_derivatives) const
4724 {
4725 AssertDimension(third_derivatives.size(),
4726 this->fe_interface->n_quadrature_points);
4727
4728 for (const unsigned int dof_index : this->fe_interface->dof_indices())
4729 for (const auto q_index : this->fe_interface->quadrature_point_indices())
4730 {
4731 if (dof_index == 0)
4732 third_derivatives[q_index] = 0.;
4733
4734 third_derivatives[q_index] +=
4735 local_dof_values[dof_index] *
4736 jump_in_third_derivatives(dof_index, q_index);
4737 }
4738 }
4739
4740
4741
4742 template <int dim, int spacedim>
4743 template <class InputVector>
4744 void
4746 const InputVector &fe_function,
4747 std::vector<
4748 solution_third_derivative_type<typename InputVector::value_type>>
4749 &third_derivatives) const
4750 {
4751 std::vector<typename InputVector::value_type> local_dof_values(
4752 this->fe_interface->n_current_interface_dofs());
4753 this->get_local_dof_values(fe_function, local_dof_values);
4754
4755 get_jump_in_function_third_derivatives_from_local_dof_values(
4756 local_dof_values, third_derivatives);
4757 }
4758} // namespace FEInterfaceViews
4759
4760#endif // DOXYGEN
4761
4763
4764#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
Tensor< 2, spacedim > average_hessian(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
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
const Triangulation< dim, spacedim >::cell_iterator get_cell(const unsigned int cell_index) 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
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
double average(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) 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)
double jump(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
const hp::QCollection< dim - 1 > & get_quadrature_collection() const
Tensor< 3, spacedim > jump_3rd_derivative(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
FEInterfaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const hp::QCollection< dim - 1 > &quadrature, const UpdateFlags update_flags)
FEFaceValuesBase< dim, spacedim > * fe_face_values
std::vector< types::global_dof_index > interface_dof_indices
Tensor< 1, spacedim > jump_gradient(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
std::unique_ptr< FESubfaceValues< dim, spacedim > > internal_fe_subface_values_neighbor
std::vector< std::array< unsigned int, 2 > > dofmap
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< 2, spacedim > jump_hessian(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) 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< 1, spacedim > average_gradient(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
Tensor< 2, spacedim > jump_in_shape_hessians(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
std::unique_ptr< hp::FEFaceValues< dim, spacedim > > internal_hp_fe_face_values_neighbor
std_cxx20::ranges::iota_view< unsigned int, unsigned int > quadrature_point_indices() const
double jump_in_shape_values(const unsigned int interface_dof_index, const unsigned int q_point, const unsigned int component=0) const
const FEInterfaceViews::Scalar< dim, spacedim > operator[](const FEValuesExtractors::Scalar &scalar) 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
const FEInterfaceViews::Vector< dim, spacedim > operator[](const FEValuesExtractors::Vector &vector) 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 average_gradient(const unsigned int interface_dof_index, const unsigned int q_point) const
gradient_type jump_in_gradients(const unsigned int interface_dof_index, const unsigned int q_point) const
hessian_type jump_hessian(const unsigned int interface_dof_index, const unsigned int q_point) const
gradient_type jump_gradient(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
value_type jump(const unsigned int interface_dof_index, const unsigned int q_point) 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
value_type average(const unsigned int interface_dof_index, const unsigned int q_point) 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
value_type average_value(const unsigned int interface_dof_index, const unsigned int q_point) const
void get_average_of_function_values(const InputVector &fe_function, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
hessian_type average_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_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
third_derivative_type jump_3rd_derivative(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
value_type average(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
third_derivative_type jump_3rd_derivative(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
value_type jump_in_values(const unsigned int interface_dof_index, const unsigned int q_point) const
value_type jump(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
gradient_type average_gradient(const unsigned int interface_dof_index, const unsigned int q_point) const
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:1049
Number value_type
Definition vector.h:128
#define DEAL_II_DEPRECATED
Definition config.h:172
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
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:488
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)
Definition hp.h:118
const types::fe_index invalid_fe_index
Definition types.h:228
static const unsigned int invalid_unsigned_int
Definition types.h:213
boost::integer_range< IncrementableType > iota_view
Definition iota_view.h:46
STL namespace.
typename internal::ProductTypeImpl< typename std::decay< T >::type, typename std::decay< U >::type >::type type