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
scratch_data.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2019 - 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
17
18#include <memory>
19
21
22namespace MeshWorker
23{
24 template <int dim, int spacedim>
26 const Mapping<dim, spacedim> & mapping,
28 const Quadrature<dim> & quadrature,
29 const UpdateFlags & update_flags,
30 const Quadrature<dim - 1> & face_quadrature,
31 const UpdateFlags & face_update_flags)
32 : mapping(&mapping)
33 , fe(&fe)
34 , cell_quadrature(quadrature)
35 , face_quadrature(face_quadrature)
36 , hp_capability_enabled(false)
37 , cell_update_flags(update_flags)
38 , neighbor_cell_update_flags(update_flags)
39 , face_update_flags(face_update_flags)
40 , neighbor_face_update_flags(face_update_flags)
41 , local_dof_indices(fe.n_dofs_per_cell())
42 , neighbor_dof_indices(fe.n_dofs_per_cell())
43 {}
44
45
46
47 template <int dim, int spacedim>
49 const Mapping<dim, spacedim> & mapping,
51 const Quadrature<dim> & quadrature,
52 const UpdateFlags & update_flags,
53 const UpdateFlags & neighbor_update_flags,
54 const Quadrature<dim - 1> & face_quadrature,
55 const UpdateFlags & face_update_flags,
56 const UpdateFlags & neighbor_face_update_flags)
57 : mapping(&mapping)
58 , fe(&fe)
59 , cell_quadrature(quadrature)
60 , face_quadrature(face_quadrature)
61 , hp_capability_enabled(false)
62 , cell_update_flags(update_flags)
63 , neighbor_cell_update_flags(neighbor_update_flags)
64 , face_update_flags(face_update_flags)
65 , neighbor_face_update_flags(neighbor_face_update_flags)
66 , local_dof_indices(fe.n_dofs_per_cell())
67 , neighbor_dof_indices(fe.n_dofs_per_cell())
68 {}
69
70
71
72 template <int dim, int spacedim>
75 const Quadrature<dim> & quadrature,
76 const UpdateFlags & update_flags,
77 const Quadrature<dim - 1> & face_quadrature,
78 const UpdateFlags & face_update_flags)
79 : ScratchData(fe.reference_cell()
80 .template get_default_linear_mapping<dim, spacedim>(),
81 fe,
82 quadrature,
83 update_flags,
84 face_quadrature,
85 face_update_flags)
86 {}
87
88
89
90 template <int dim, int spacedim>
93 const Quadrature<dim> & quadrature,
94 const UpdateFlags & update_flags,
95 const UpdateFlags & neighbor_update_flags,
96 const Quadrature<dim - 1> & face_quadrature,
97 const UpdateFlags & face_update_flags,
98 const UpdateFlags & neighbor_face_update_flags)
99 : ScratchData(fe.reference_cell()
100 .template get_default_linear_mapping<dim, spacedim>(),
101 fe,
102 quadrature,
103 update_flags,
104 neighbor_update_flags,
105 face_quadrature,
106 face_update_flags,
107 neighbor_face_update_flags)
108 {}
109
110
111
112 template <int dim, int spacedim>
114 const hp::MappingCollection<dim, spacedim> &mapping_collection,
115 const hp::FECollection<dim, spacedim> & fe_collection,
116 const hp::QCollection<dim> & cell_quadrature_collection,
117 const UpdateFlags & cell_update_flags,
118 const hp::QCollection<dim - 1> & face_quadrature_collection,
119 const UpdateFlags & face_update_flags)
120 : mapping_collection(&mapping_collection)
121 , fe_collection(&fe_collection)
122 , cell_quadrature_collection(cell_quadrature_collection)
123 , face_quadrature_collection(face_quadrature_collection)
124 , hp_capability_enabled(true)
125 , cell_update_flags(cell_update_flags)
126 , neighbor_cell_update_flags(cell_update_flags)
127 , face_update_flags(face_update_flags)
128 , neighbor_face_update_flags(face_update_flags)
129 {
130 local_dof_indices.reserve(fe_collection.max_dofs_per_cell());
131 neighbor_dof_indices.reserve(fe_collection.max_dofs_per_cell());
132 }
133
134
135
136 template <int dim, int spacedim>
138 const hp::MappingCollection<dim, spacedim> &mapping_collection,
139 const hp::FECollection<dim, spacedim> & fe_collection,
140 const hp::QCollection<dim> & cell_quadrature_collection,
141 const UpdateFlags & cell_update_flags,
142 const UpdateFlags & neighbor_cell_update_flags,
143 const hp::QCollection<dim - 1> & face_quadrature_collection,
144 const UpdateFlags & face_update_flags,
145 const UpdateFlags & neighbor_face_update_flags)
146 : mapping_collection(&mapping_collection)
147 , fe_collection(&fe_collection)
148 , cell_quadrature_collection(cell_quadrature_collection)
149 , face_quadrature_collection(face_quadrature_collection)
150 , hp_capability_enabled(true)
151 , cell_update_flags(cell_update_flags)
152 , neighbor_cell_update_flags(neighbor_cell_update_flags)
153 , face_update_flags(face_update_flags)
154 , neighbor_face_update_flags(neighbor_face_update_flags)
155 {
156 local_dof_indices.reserve(fe_collection.max_dofs_per_cell());
157 neighbor_dof_indices.reserve(fe_collection.max_dofs_per_cell());
158 }
159
160
161
162 template <int dim, int spacedim>
164 const hp::FECollection<dim, spacedim> &fe_collection,
165 const hp::QCollection<dim> & cell_quadrature_collection,
166 const UpdateFlags & cell_update_flags,
167 const hp::QCollection<dim - 1> & face_quadrature_collection,
168 const UpdateFlags & face_update_flags)
169 : ScratchData(fe_collection.get_reference_cell_default_linear_mapping(),
170 fe_collection,
171 cell_quadrature_collection,
172 cell_update_flags,
173 face_quadrature_collection,
174 face_update_flags)
175 {}
176
177
178
179 template <int dim, int spacedim>
181 const hp::FECollection<dim, spacedim> &fe_collection,
182 const hp::QCollection<dim> & cell_quadrature_collection,
183 const UpdateFlags & cell_update_flags,
184 const UpdateFlags & neighbor_cell_update_flags,
185 const hp::QCollection<dim - 1> & face_quadrature_collection,
186 const UpdateFlags & face_update_flags,
187 const UpdateFlags & neighbor_face_update_flags)
188 : ScratchData(fe_collection.get_reference_cell_default_linear_mapping(),
189 fe_collection,
190 cell_quadrature_collection,
191 cell_update_flags,
192 neighbor_cell_update_flags,
193 face_quadrature_collection,
194 face_update_flags,
195 neighbor_face_update_flags)
196 {}
197
198
199
200 template <int dim, int spacedim>
202 const ScratchData<dim, spacedim> &scratch)
203 : mapping(scratch.mapping)
204 , fe(scratch.fe)
205 , cell_quadrature(scratch.cell_quadrature)
206 , face_quadrature(scratch.face_quadrature)
207 , mapping_collection(scratch.mapping_collection)
208 , fe_collection(scratch.fe_collection)
209 , cell_quadrature_collection(scratch.cell_quadrature_collection)
210 , face_quadrature_collection(scratch.face_quadrature_collection)
211 , hp_capability_enabled(scratch.hp_capability_enabled)
212 , cell_update_flags(scratch.cell_update_flags)
213 , neighbor_cell_update_flags(scratch.neighbor_cell_update_flags)
214 , face_update_flags(scratch.face_update_flags)
215 , neighbor_face_update_flags(scratch.neighbor_face_update_flags)
216 , local_dof_indices(scratch.local_dof_indices)
217 , neighbor_dof_indices(scratch.neighbor_dof_indices)
218 , user_data_storage(scratch.user_data_storage)
219 , internal_data_storage(scratch.internal_data_storage)
220 {}
221
222
223
224 template <int dim, int spacedim>
228 {
229 if (hp_capability_enabled == false)
230 {
231 if (!fe_values)
232 fe_values = std::make_unique<FEValues<dim, spacedim>>(
233 *mapping, *fe, cell_quadrature, cell_update_flags);
234
235 fe_values->reinit(cell);
236 local_dof_indices.resize(fe_values->dofs_per_cell);
237 cell->get_dof_indices(local_dof_indices);
238 current_fe_values = fe_values.get();
239 return *fe_values;
240 }
241 else
242 {
243 if (!hp_fe_values)
244 hp_fe_values = std::make_unique<hp::FEValues<dim, spacedim>>(
245 *mapping_collection,
246 *fe_collection,
247 cell_quadrature_collection,
248 cell_update_flags);
249
250 hp_fe_values->reinit(cell);
251 const auto &fe_values = hp_fe_values->get_present_fe_values();
252
254 (*fe_collection)[cell->active_fe_index()].n_dofs_per_cell(),
255 fe_values.dofs_per_cell);
256 local_dof_indices.resize(fe_values.dofs_per_cell);
257 cell->get_dof_indices(local_dof_indices);
258
259 current_fe_values = &fe_values;
260 return fe_values;
261 }
262 }
263
264
265
266 template <int dim, int spacedim>
270 const unsigned int face_no)
271 {
272 if (hp_capability_enabled == false)
273 {
274 if (!fe_face_values)
275 fe_face_values = std::make_unique<FEFaceValues<dim, spacedim>>(
276 *mapping, *fe, face_quadrature, face_update_flags);
277
278 fe_face_values->reinit(cell, face_no);
279 local_dof_indices.resize(fe->n_dofs_per_cell());
280 cell->get_dof_indices(local_dof_indices);
281 current_fe_values = fe_face_values.get();
282 return *fe_face_values;
283 }
284 else
285 {
286 return reinit(cell, cell, face_no);
287 }
288 }
289
290
291
292 template <int dim, int spacedim>
297 & neighbor_cell,
298 const unsigned int face_no)
299 {
300 Assert(hp_capability_enabled, ExcOnlyAvailableWithHP());
301
302 if (!hp_fe_face_values)
303 hp_fe_face_values = std::make_unique<hp::FEFaceValues<dim, spacedim>>(
304 *mapping_collection,
305 *fe_collection,
306 face_quadrature_collection,
307 face_update_flags);
308
309 if (neighbor_cell == cell)
310 {
311 hp_fe_face_values->reinit(cell, face_no);
312 }
313 else
314 {
315 // When we want to ensure some agreement between the cell face and its
316 // neighbor on the quadrature order and mapping to use on this face,
317 // then we defer to the dominance of one FE over another. This should
318 // ensure that the optimal integration order and mapping order are
319 // selected for this situation.
320 unsigned int dominated_fe_index = fe_collection->find_dominated_fe(
321 {cell->active_fe_index(), neighbor_cell->active_fe_index()});
322
323 // TODO: find_dominated_fe returns invalid_fe_index when no dominated FE
324 // has been found. We want to pass this value to FEFaceValues, but it
325 // expects an invalid_unsigned_int in this case. We need to match the
326 // interfaces in the future.
327 if (dominated_fe_index == numbers::invalid_fe_index)
328 dominated_fe_index = numbers::invalid_unsigned_int;
329
330 hp_fe_face_values->reinit(cell,
331 face_no,
332 dominated_fe_index,
333 dominated_fe_index);
334 }
335
336 const auto &fe_face_values = hp_fe_face_values->get_present_fe_values();
337 const auto &fe = (*fe_collection)[cell->active_fe_index()];
338
339 local_dof_indices.resize(fe.n_dofs_per_cell());
340 cell->get_dof_indices(local_dof_indices);
341
342 current_fe_values = &fe_face_values;
343 return fe_face_values;
344 }
345
346
347
348 template <int dim, int spacedim>
352 const unsigned int face_no,
353 const unsigned int subface_no)
354 {
355 if (subface_no != numbers::invalid_unsigned_int)
356 {
357 if (hp_capability_enabled == false)
358 {
359 if (!fe_subface_values)
360 fe_subface_values =
361 std::make_unique<FESubfaceValues<dim, spacedim>>(
362 *mapping, *fe, face_quadrature, face_update_flags);
363 fe_subface_values->reinit(cell, face_no, subface_no);
364 local_dof_indices.resize(fe->n_dofs_per_cell());
365 cell->get_dof_indices(local_dof_indices);
366
367 current_fe_values = fe_subface_values.get();
368 return *fe_subface_values;
369 }
370 else
371 {
372 return reinit(cell, cell, face_no, subface_no);
373 }
374 }
375 else
376 return reinit(cell, face_no);
377 }
378
379
380
381 template <int dim, int spacedim>
386 & neighbor_cell,
387 const unsigned int face_no,
388 const unsigned int subface_no)
389 {
390 Assert(hp_capability_enabled, ExcOnlyAvailableWithHP());
391
392 if (subface_no != numbers::invalid_unsigned_int)
393 {
394 if (!hp_fe_subface_values)
395 hp_fe_subface_values =
396 std::make_unique<hp::FESubfaceValues<dim, spacedim>>(
397 *mapping_collection,
398 *fe_collection,
399 face_quadrature_collection,
400 face_update_flags);
401
402 if (neighbor_cell == cell)
403 {
404 hp_fe_subface_values->reinit(cell, face_no, subface_no);
405 }
406 else
407 {
408 // When we want to ensure some agreement between the cell face and
409 // its neighbor on the quadrature order and mapping to use on this
410 // face, then we defer to the dominance of one FE over another. This
411 // should ensure that the optimal integration order and mapping
412 // order are selected for this situation.
413 unsigned int dominated_fe_index = fe_collection->find_dominated_fe(
414 {cell->active_fe_index(), neighbor_cell->active_fe_index()});
415
416 // TODO: find_dominated_fe returns invalid_fe_index when no
417 // dominated FE has been found. We want to pass this value to
418 // FEFaceValues, but it expects an invalid_unsigned_int in this
419 // case. We need to match the interfaces in the future.
420 if (dominated_fe_index == numbers::invalid_fe_index)
421 dominated_fe_index = numbers::invalid_unsigned_int;
422
423 hp_fe_subface_values->reinit(cell,
424 face_no,
425 subface_no,
426 dominated_fe_index,
427 dominated_fe_index);
428 }
429
430 const auto &fe_subface_values =
431 hp_fe_subface_values->get_present_fe_values();
432 const auto &fe = (*fe_collection)[cell->active_fe_index()];
433
434 local_dof_indices.resize(fe.n_dofs_per_cell());
435 cell->get_dof_indices(local_dof_indices);
436
437 current_fe_values = &fe_subface_values;
438 return fe_subface_values;
439 }
440 else
441 {
442 return reinit(cell, neighbor_cell, face_no);
443 }
444 }
445
446
447
448 template <int dim, int spacedim>
452 const unsigned int face_no,
454 & cell_neighbor,
455 const unsigned int face_no_neighbor)
456 {
457 return reinit(cell,
458 face_no,
460 cell_neighbor,
461 face_no_neighbor,
463 }
464
465
466
467 template <int dim, int spacedim>
471 const unsigned int face_no,
472 const unsigned int sub_face_no,
474 & cell_neighbor,
475 const unsigned int face_no_neighbor,
476 const unsigned int sub_face_no_neighbor)
477 {
478 if (hp_capability_enabled == false)
479 {
480 if (!interface_fe_values)
481 interface_fe_values =
482 std::make_unique<FEInterfaceValues<dim, spacedim>>(
483 *mapping, *fe, face_quadrature, face_update_flags);
484
485 interface_fe_values->reinit(cell,
486 face_no,
487 sub_face_no,
488 cell_neighbor,
489 face_no_neighbor,
490 sub_face_no_neighbor);
491 }
492 else
493 {
494 if (!interface_fe_values)
495 interface_fe_values =
496 std::make_unique<FEInterfaceValues<dim, spacedim>>(
497 *mapping_collection,
498 *fe_collection,
499 face_quadrature_collection,
500 face_update_flags);
501
502 // When we want to ensure some agreement between the cell face and
503 // its neighbor on the quadrature order and mapping to use on this
504 // face, then we defer to the dominance of one FE over another. This
505 // should ensure that the optimal integration order and mapping
506 // order are selected for this situation.
507 unsigned int dominated_fe_index = fe_collection->find_dominated_fe(
508 {cell->active_fe_index(), cell_neighbor->active_fe_index()});
509
510 // TODO: find_dominated_fe returns invalid_fe_index when no dominated FE
511 // has been found. We want to pass this value to FEFaceValues, but it
512 // expects an invalid_unsigned_int in this case. We need to match the
513 // interfaces in the future.
514 if (dominated_fe_index == numbers::invalid_fe_index)
515 dominated_fe_index = numbers::invalid_unsigned_int;
516
517 interface_fe_values->reinit(cell,
518 face_no,
519 sub_face_no,
520 cell_neighbor,
521 face_no_neighbor,
522 sub_face_no_neighbor,
523 dominated_fe_index,
524 dominated_fe_index);
525 }
526
527 current_fe_values = &interface_fe_values->get_fe_face_values(0);
528 current_neighbor_fe_values = &interface_fe_values->get_fe_face_values(1);
529
530 local_dof_indices = interface_fe_values->get_interface_dof_indices();
531 return *interface_fe_values;
532 }
533
534
535
536 template <int dim, int spacedim>
540 {
541 if (hp_capability_enabled == false)
542 {
543 if (!neighbor_fe_values)
544 neighbor_fe_values = std::make_unique<FEValues<dim, spacedim>>(
545 *mapping, *fe, cell_quadrature, neighbor_cell_update_flags);
546
547 neighbor_fe_values->reinit(cell);
548 cell->get_dof_indices(neighbor_dof_indices);
549 current_neighbor_fe_values = neighbor_fe_values.get();
550 return *neighbor_fe_values;
551 }
552 else
553 {
554 if (!neighbor_hp_fe_values)
555 neighbor_hp_fe_values = std::make_unique<hp::FEValues<dim, spacedim>>(
556 *mapping_collection,
557 *fe_collection,
558 cell_quadrature_collection,
559 neighbor_cell_update_flags);
560
561 neighbor_hp_fe_values->reinit(cell);
562 const auto &neighbor_fe_values =
563 neighbor_hp_fe_values->get_present_fe_values();
564
566 (*fe_collection)[cell->active_fe_index()].n_dofs_per_cell(),
567 neighbor_fe_values.dofs_per_cell);
568 neighbor_dof_indices.resize(neighbor_fe_values.dofs_per_cell);
569 cell->get_dof_indices(neighbor_dof_indices);
570
571 current_neighbor_fe_values = &neighbor_fe_values;
572 return neighbor_fe_values;
573 }
574 }
575
576
577
578 template <int dim, int spacedim>
582 const unsigned int face_no)
583 {
584 if (hp_capability_enabled == false)
585 {
586 if (!neighbor_fe_face_values)
587 neighbor_fe_face_values =
588 std::make_unique<FEFaceValues<dim, spacedim>>(
589 *mapping, *fe, face_quadrature, neighbor_face_update_flags);
590 neighbor_fe_face_values->reinit(cell, face_no);
591 cell->get_dof_indices(neighbor_dof_indices);
592 current_neighbor_fe_values = neighbor_fe_face_values.get();
593 return *neighbor_fe_face_values;
594 }
595 else
596 {
597 return reinit_neighbor(cell, cell, face_no);
598 }
599 }
600
601
602
603 template <int dim, int spacedim>
608 & neighbor_cell,
609 const unsigned int face_no)
610 {
611 Assert(hp_capability_enabled, ExcOnlyAvailableWithHP());
612
613 if (!neighbor_hp_fe_face_values)
614 neighbor_hp_fe_face_values =
615 std::make_unique<hp::FEFaceValues<dim, spacedim>>(
616 *mapping_collection,
617 *fe_collection,
618 face_quadrature_collection,
619 neighbor_face_update_flags);
620
621 if (neighbor_cell == cell)
622 {
623 neighbor_hp_fe_face_values->reinit(neighbor_cell, face_no);
624 }
625 else
626 {
627 // When we want to ensure some agreement between the cell face and its
628 // neighbor on the quadrature order and mapping to use on this face,
629 // then we defer to the dominance of one FE over another. This should
630 // ensure that the optimal integration order and mapping order are
631 // selected for this situation.
632 unsigned int dominated_fe_index = fe_collection->find_dominated_fe(
633 {cell->active_fe_index(), neighbor_cell->active_fe_index()});
634
635 // TODO: find_dominated_fe returns invalid_fe_index when no dominated FE
636 // has been found. We want to pass this value to FEFaceValues, but it
637 // expects an invalid_unsigned_int in this case. We need to match the
638 // interfaces in the future.
639 if (dominated_fe_index == numbers::invalid_fe_index)
640 dominated_fe_index = numbers::invalid_unsigned_int;
641
642 neighbor_hp_fe_face_values->reinit(neighbor_cell,
643 face_no,
644 dominated_fe_index,
645 dominated_fe_index);
646 }
647
648 const auto &neighbor_fe_face_values =
649 neighbor_hp_fe_face_values->get_present_fe_values();
650 const auto &neighbor_fe =
651 (*fe_collection)[neighbor_cell->active_fe_index()];
652
653 neighbor_dof_indices.resize(neighbor_fe.n_dofs_per_cell());
654 neighbor_cell->get_dof_indices(neighbor_dof_indices);
655
656 current_neighbor_fe_values = &neighbor_fe_face_values;
657 return neighbor_fe_face_values;
658 }
659
660
661
662 template <int dim, int spacedim>
666 const unsigned int face_no,
667 const unsigned int subface_no)
668 {
669 if (subface_no != numbers::invalid_unsigned_int)
670 {
671 if (hp_capability_enabled == false)
672 {
673 if (!neighbor_fe_subface_values)
674 neighbor_fe_subface_values =
675 std::make_unique<FESubfaceValues<dim, spacedim>>(
676 *mapping, *fe, face_quadrature, neighbor_face_update_flags);
677 neighbor_fe_subface_values->reinit(cell, face_no, subface_no);
678 cell->get_dof_indices(neighbor_dof_indices);
679 current_neighbor_fe_values = neighbor_fe_subface_values.get();
680 return *neighbor_fe_subface_values;
681 }
682 else
683 {
684 return reinit_neighbor(cell, cell, face_no, subface_no);
685 }
686 }
687 else
688 return reinit_neighbor(cell, face_no);
689 }
690
691
692
693 template <int dim, int spacedim>
698 & neighbor_cell,
699 const unsigned int face_no,
700 const unsigned int subface_no)
701 {
702 Assert(hp_capability_enabled, ExcOnlyAvailableWithHP());
703
704 if (subface_no != numbers::invalid_unsigned_int)
705 {
706 if (!neighbor_hp_fe_subface_values)
707 neighbor_hp_fe_subface_values =
708 std::make_unique<hp::FESubfaceValues<dim, spacedim>>(
709 *mapping_collection,
710 *fe_collection,
711 face_quadrature_collection,
712 neighbor_face_update_flags);
713
714 if (neighbor_cell == cell)
715 {
716 neighbor_hp_fe_subface_values->reinit(neighbor_cell,
717 face_no,
718 subface_no);
719 }
720 else
721 {
722 // When we want to ensure some agreement between the cell face and
723 // its neighbor on the quadrature order and mapping to use on this
724 // face, then we defer to the dominance of one FE over another. This
725 // should ensure that the optimal integration order and mapping
726 // order are selected for this situation.
727 unsigned int dominated_fe_index = fe_collection->find_dominated_fe(
728 {cell->active_fe_index(), neighbor_cell->active_fe_index()});
729
730 // TODO: find_dominated_fe returns invalid_fe_index when no
731 // dominated FE has been found. We want to pass this value to
732 // FEFaceValues, but it expects an invalid_unsigned_int in this
733 // case. We need to match the interfaces in the future.
734 if (dominated_fe_index == numbers::invalid_fe_index)
735 dominated_fe_index = numbers::invalid_unsigned_int;
736
737 neighbor_hp_fe_subface_values->reinit(neighbor_cell,
738 face_no,
739 subface_no,
740 dominated_fe_index,
741 dominated_fe_index);
742 }
743
744 const auto &neighbor_fe_subface_values =
745 neighbor_hp_fe_subface_values->get_present_fe_values();
746 const auto &neighbor_fe =
747 (*fe_collection)[neighbor_cell->active_fe_index()];
748
749 neighbor_dof_indices.resize(neighbor_fe.n_dofs_per_cell());
750 neighbor_cell->get_dof_indices(neighbor_dof_indices);
751
752 current_neighbor_fe_values = &neighbor_fe_subface_values;
753 return neighbor_fe_subface_values;
754 }
755 else
756 {
757 return reinit_neighbor(cell, neighbor_cell, face_no);
758 }
759 }
760
761
762
763 template <int dim, int spacedim>
766 {
767 Assert(current_fe_values != nullptr,
768 ExcMessage("You have to initialize the cache using one of the "
769 "reinit functions first!"));
770 return *current_fe_values;
771 }
772
773
774
775 template <int dim, int spacedim>
778 {
779 Assert(interface_fe_values != nullptr,
780 ExcMessage("You have to initialize the cache using one of the "
781 "reinit functions first!"));
782 return *interface_fe_values;
783 }
784
785
786
787 template <int dim, int spacedim>
790 {
791 Assert(current_neighbor_fe_values != nullptr,
792 ExcMessage("You have to initialize the cache using one of the "
793 "reinit functions first!"));
794 return *current_neighbor_fe_values;
795 }
796
797
798
799 template <int dim, int spacedim>
800 const std::vector<Point<spacedim>> &
802 {
803 return get_current_fe_values().get_quadrature_points();
804 }
805
806
807
808 template <int dim, int spacedim>
809 const std::vector<double> &
811 {
812 return get_current_fe_values().get_JxW_values();
813 }
814
815
816
817 template <int dim, int spacedim>
818 const std::vector<double> &
820 {
821 return get_current_neighbor_fe_values().get_JxW_values();
822 }
823
824
825
826 template <int dim, int spacedim>
827 const std::vector<Tensor<1, spacedim>> &
829 {
830 return get_current_fe_values().get_normal_vectors();
831 }
832
833
834
835 template <int dim, int spacedim>
836 const std::vector<Tensor<1, spacedim>> &
838 {
839 return get_current_neighbor_fe_values().get_normal_vectors();
840 }
841
842
843
844 template <int dim, int spacedim>
845 const std::vector<types::global_dof_index> &
847 {
848 return local_dof_indices;
849 }
850
851
852
853 template <int dim, int spacedim>
854 unsigned int
856 {
857 return local_dof_indices.size();
858 }
859
860
861
862 template <int dim, int spacedim>
863 const std::vector<types::global_dof_index> &
865 {
866 return neighbor_dof_indices;
867 }
868
869
870
871 template <int dim, int spacedim>
872 unsigned int
874 {
875 return neighbor_dof_indices.size();
876 }
877
878
879
880 template <int dim, int spacedim>
883 {
884 return user_data_storage;
885 }
886
887
888
889 template <int dim, int spacedim>
890 const GeneralDataStorage &
892 {
893 return user_data_storage;
894 }
895
896
897
898 template <int dim, int spacedim>
901 {
902 Assert(hp_capability_enabled == false, ExcOnlyAvailableWithoutHP());
903 Assert(mapping, ExcNotInitialized());
904 return *mapping;
905 }
906
907
908
909 template <int dim, int spacedim>
912 {
913 Assert(hp_capability_enabled == false, ExcOnlyAvailableWithoutHP());
915 return *fe;
916 }
917
918
919
920 template <int dim, int spacedim>
921 const Quadrature<dim> &
923 {
924 Assert(hp_capability_enabled == false, ExcOnlyAvailableWithoutHP());
925 return cell_quadrature;
926 }
927
928
929
930 template <int dim, int spacedim>
931 const Quadrature<dim - 1> &
933 {
934 Assert(hp_capability_enabled == false, ExcOnlyAvailableWithoutHP());
935 return face_quadrature;
936 }
937
938
939
940 template <int dim, int spacedim>
943 {
944 Assert(hp_capability_enabled == true, ExcOnlyAvailableWithHP());
945 Assert(mapping_collection, ExcNotInitialized());
946 return *mapping_collection;
947 }
948
949
950
951 template <int dim, int spacedim>
954 {
955 Assert(hp_capability_enabled == true, ExcOnlyAvailableWithHP());
956 Assert(fe_collection, ExcNotInitialized());
957 return *fe_collection;
958 }
959
960
961
962 template <int dim, int spacedim>
965 {
966 Assert(hp_capability_enabled == true, ExcOnlyAvailableWithHP());
967 return cell_quadrature_collection;
968 }
969
970
971
972 template <int dim, int spacedim>
973 const hp::QCollection<dim - 1> &
975 {
976 Assert(hp_capability_enabled == true, ExcOnlyAvailableWithHP());
977 return face_quadrature_collection;
978 }
979
980
981
982 template <int dim, int spacedim>
983 bool
985 {
986 return hp_capability_enabled;
987 }
988
989
990
991 template <int dim, int spacedim>
994 {
995 return cell_update_flags;
996 }
997
998
999
1000 template <int dim, int spacedim>
1003 {
1004 return neighbor_cell_update_flags;
1005 }
1006
1007
1008
1009 template <int dim, int spacedim>
1012 {
1013 return face_update_flags;
1014 }
1015
1016
1017
1018 template <int dim, int spacedim>
1021 {
1022 return neighbor_face_update_flags;
1023 }
1024
1025} // namespace MeshWorker
1027
1028// Explicit instantiations
1030namespace MeshWorker
1031{
1032#include "scratch_data.inst"
1033}
void reinit(const Triangulation< dim, spacedim > &tria)
Abstract base class for mapping classes.
Definition mapping.h:317
bool has_hp_capabilities() const
const hp::MappingCollection< dim, spacedim > & get_mapping_collection() const
const FEValuesBase< dim, spacedim > & get_current_fe_values() const
GeneralDataStorage & get_general_data_storage()
const Quadrature< dim - 1 > & get_face_quadrature() const
const std::vector< Tensor< 1, spacedim > > & get_neighbor_normal_vectors()
unsigned int n_dofs_per_cell() const
unsigned int n_neighbor_dofs_per_cell() const
const std::vector< double > & get_JxW_values() const
const Mapping< dim, spacedim > & get_mapping() const
const std::vector< double > & get_neighbor_JxW_values() const
const FiniteElement< dim, spacedim > & get_fe() const
SmartPointer< const hp::FECollection< dim, spacedim > > fe_collection
const hp::QCollection< dim > & get_cell_quadrature_collection() const
std::vector< types::global_dof_index > neighbor_dof_indices
const std::vector< types::global_dof_index > & get_neighbor_dof_indices() const
UpdateFlags get_face_update_flags() const
UpdateFlags get_neighbor_cell_update_flags() const
const hp::FECollection< dim, spacedim > & get_fe_collection() const
const std::vector< Tensor< 1, spacedim > > & get_normal_vectors() const
const FEValues< dim, spacedim > & reinit_neighbor(const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell)
std::vector< types::global_dof_index > local_dof_indices
const std::vector< Point< spacedim > > & get_quadrature_points() const
const FEValuesBase< dim, spacedim > & get_current_neighbor_fe_values() const
const std::vector< types::global_dof_index > & get_local_dof_indices() const
const FEInterfaceValues< dim, spacedim > & get_current_interface_fe_values() const
UpdateFlags get_cell_update_flags() const
const Quadrature< dim > & get_cell_quadrature() const
const hp::QCollection< dim - 1 > & get_face_quadrature_collection() const
ScratchData(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, const UpdateFlags &update_flags, const Quadrature< dim - 1 > &face_quadrature=Quadrature< dim - 1 >(), const UpdateFlags &face_update_flags=update_default)
const FEValues< dim, spacedim > & reinit(const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell)
UpdateFlags get_neighbor_face_update_flags() const
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, lda > > &cell, 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)
Definition fe_values.cc:294
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcNotInitialized()
static ::ExceptionBase & ExcMessage(std::string arg1)
typename ActiveSelector::active_cell_iterator active_cell_iterator
UpdateFlags
const Mapping< dim, spacedim > & get_default_linear_mapping(const Triangulation< dim, spacedim > &triangulation)
Definition mapping.cc:285
const types::fe_index invalid_fe_index
Definition types.h:228
static const unsigned int invalid_unsigned_int
Definition types.h:213