deal.II version 9.7.1
\(\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
face_setup_internal.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2018 - 2025 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15
16#ifndef dealii_face_setup_internal_h
17#define dealii_face_setup_internal_h
18
19#include <deal.II/base/config.h>
20
21#include <deal.II/base/memory_consumption.h>
22
23#include <deal.II/distributed/tria_base.h>
24
25#include <deal.II/grid/tria.h>
26#include <deal.II/grid/tria_accessor.h>
27
28#include <deal.II/matrix_free/face_info.h>
29#include <deal.II/matrix_free/shape_info.h>
30#include <deal.II/matrix_free/task_info.h>
31
32#include <fstream>
33#include <set>
34
35
36DEAL_II_NAMESPACE_OPEN
37
38
39namespace internal
40{
41 namespace MatrixFreeFunctions
42 {
58
59
60
69 template <int dim>
70 struct FaceSetup
71 {
72 FaceSetup();
73
80 void
81 initialize(
82 const ::Triangulation<dim> &triangulation,
83 const unsigned int mg_level,
84 const bool hold_all_faces_to_owned_cells,
85 const bool build_inner_faces,
86 std::vector<std::pair<unsigned int, unsigned int>> &cell_levels);
87
94 void
95 generate_faces(
96 const ::Triangulation<dim> &triangulation,
97 const std::vector<std::pair<unsigned int, unsigned int>> &cell_levels,
98 TaskInfo &task_info);
99
106 FaceToCellTopology<1>
107 create_face(
108 const unsigned int face_no,
109 const typename ::Triangulation<dim>::cell_iterator &cell,
110 const unsigned int number_cell_interior,
111 const typename ::Triangulation<dim>::cell_iterator &neighbor,
112 const unsigned int number_cell_exterior,
113 const bool is_mixed_mesh);
114
115 bool use_active_cells;
116
129
130 std::vector<FaceCategory> face_is_owned;
131 std::vector<bool> at_processor_boundary;
132 std::vector<FaceToCellTopology<1>> inner_faces;
133 std::vector<FaceToCellTopology<1>> boundary_faces;
134 std::vector<FaceToCellTopology<1>> inner_ghost_faces;
135 std::vector<FaceToCellTopology<1>> refinement_edge_faces;
136 };
137
138
139
143 template <int vectorization_width>
144 void
145 collect_faces_vectorization(
146 const std::vector<FaceToCellTopology<1>> &faces_in,
147 const std::vector<bool> &hard_vectorization_boundary,
148 std::vector<unsigned int> &face_partition_data,
149 std::vector<FaceToCellTopology<vectorization_width>> &faces_out);
150
151
152
153 /* -------------------------------------------------------------------- */
154
155#ifndef DOXYGEN
156
157 template <int dim>
158 FaceSetup<dim>::FaceSetup()
159 : use_active_cells(true)
160 {}
161
162
163
164 template <int dim>
165 void
166 FaceSetup<dim>::initialize(
167 const ::Triangulation<dim> &triangulation,
168 const unsigned int mg_level,
169 const bool hold_all_faces_to_owned_cells,
170 const bool build_inner_faces,
171 std::vector<std::pair<unsigned int, unsigned int>> &cell_levels)
172 {
173 use_active_cells = mg_level == numbers::invalid_unsigned_int;
174
175 if constexpr (running_in_debug_mode())
176 {
177 // safety check
178 if (use_active_cells)
179 for (const auto &cell_level : cell_levels)
180 {
181 typename ::Triangulation<dim>::cell_iterator dcell(
182 &triangulation, cell_level.first, cell_level.second);
183 Assert(dcell->is_active(), ExcInternalError());
184 }
185 }
186
187 // step 1: add ghost cells for those cells that we identify as
188 // interesting
189
190 at_processor_boundary.resize(cell_levels.size(), false);
191 face_is_owned.resize(dim > 1 ? triangulation.n_raw_faces() :
192 triangulation.n_vertices(),
193 FaceCategory::locally_active_done_elsewhere);
194
195 // go through the mesh and divide the faces on the processor
196 // boundaries as evenly as possible between the processors
197 std::map<types::subdomain_id, FaceIdentifier>
198 inner_faces_at_proc_boundary;
199 if (triangulation.locally_owned_subdomain() !=
200 numbers::invalid_subdomain_id)
201 {
202 const types::subdomain_id my_domain =
203 triangulation.locally_owned_subdomain();
204 for (unsigned int i = 0; i < cell_levels.size(); ++i)
205 {
206 if (i > 0 && cell_levels[i] == cell_levels[i - 1])
207 continue;
208 typename ::Triangulation<dim>::cell_iterator dcell(
209 &triangulation, cell_levels[i].first, cell_levels[i].second);
210 for (const unsigned int f : dcell->face_indices())
211 {
212 if (dcell->at_boundary(f) && !dcell->has_periodic_neighbor(f))
213 continue;
214 typename ::Triangulation<dim>::cell_iterator neighbor =
215 dcell->neighbor_or_periodic_neighbor(f);
216
217 // faces at hanging nodes are always treated by the processor
218 // who owns the element on the fine side. but we need to count
219 // the number of inner faces in order to balance the remaining
220 // faces properly
221 const CellId id_mine = dcell->id();
222 if (use_active_cells && neighbor->has_children())
223 for (unsigned int c = 0;
224 c < (dcell->has_periodic_neighbor(f) ?
225 dcell->periodic_neighbor(f)
226 ->face(dcell->periodic_neighbor_face_no(f))
227 ->n_children() :
228 dcell->face(f)->n_children());
229 ++c)
230 {
231 typename ::Triangulation<dim>::cell_iterator
232 neighbor_c =
233 dcell->at_boundary(f) ?
234 dcell->periodic_neighbor_child_on_subface(f, c) :
235 dcell->neighbor_child_on_subface(f, c);
236 const types::subdomain_id neigh_domain =
237 neighbor_c->subdomain_id();
238 if (my_domain < neigh_domain)
239 inner_faces_at_proc_boundary[neigh_domain]
240 .n_hanging_faces_larger_subdomain++;
241 else if (my_domain > neigh_domain)
242 inner_faces_at_proc_boundary[neigh_domain]
243 .n_hanging_faces_smaller_subdomain++;
244 }
245 else
246 {
247 const types::subdomain_id neigh_domain =
248 use_active_cells ? neighbor->subdomain_id() :
249 neighbor->level_subdomain_id();
250 if (neighbor->level() < dcell->level() &&
251 use_active_cells)
252 {
253 if (my_domain < neigh_domain)
254 inner_faces_at_proc_boundary[neigh_domain]
255 .n_hanging_faces_smaller_subdomain++;
256 else if (my_domain > neigh_domain)
257 inner_faces_at_proc_boundary[neigh_domain]
258 .n_hanging_faces_larger_subdomain++;
259 }
260 else if (neighbor->level() == dcell->level() &&
261 my_domain != neigh_domain)
262 {
263 // always list the cell whose owner has the lower
264 // subdomain id first. this applies to both processors
265 // involved, so both processors will generate the same
266 // list that we will later order
267 const CellId id_neigh = neighbor->id();
268 if (my_domain < neigh_domain)
269 inner_faces_at_proc_boundary[neigh_domain]
270 .shared_faces.emplace_back(id_mine, id_neigh);
271 else
272 inner_faces_at_proc_boundary[neigh_domain]
273 .shared_faces.emplace_back(id_neigh, id_mine);
274 }
275 }
276 }
277 }
278
279 // sort the cell ids related to each neighboring processor. This
280 // algorithm is symmetric so every processor combination should
281 // arrive here and no deadlock should be possible
282 for (auto &inner_face : inner_faces_at_proc_boundary)
283 {
284 Assert(inner_face.first != my_domain,
285 ExcInternalError("Should not send info to myself"));
286 std::sort(inner_face.second.shared_faces.begin(),
287 inner_face.second.shared_faces.end());
288 inner_face.second.shared_faces.erase(
289 std::unique(inner_face.second.shared_faces.begin(),
290 inner_face.second.shared_faces.end()),
291 inner_face.second.shared_faces.end());
292
293 // safety check: both involved processors should see the same list
294 // because the pattern of ghosting is symmetric. We test this by
295 // looking at the length of the lists of faces
296# if defined(DEAL_II_WITH_MPI) && defined(DEBUG)
297 MPI_Comm comm = MPI_COMM_SELF;
298 if (const ::parallel::TriangulationBase<dim> *ptria =
299 dynamic_cast<const ::parallel::TriangulationBase<dim>
300 *>(&triangulation))
301 comm = ptria->get_mpi_communicator();
302
303 MPI_Status status;
304 unsigned int mysize = inner_face.second.shared_faces.size();
305 unsigned int othersize = numbers::invalid_unsigned_int;
306
307 int ierr = MPI_Sendrecv(&mysize,
308 1,
309 MPI_UNSIGNED,
310 inner_face.first,
311 600 + my_domain,
312 &othersize,
313 1,
314 MPI_UNSIGNED,
315 inner_face.first,
316 600 + inner_face.first,
317 comm,
318 &status);
319 AssertThrowMPI(ierr);
320 AssertDimension(mysize, othersize);
321 mysize = inner_face.second.n_hanging_faces_smaller_subdomain;
322 ierr = MPI_Sendrecv(&mysize,
323 1,
324 MPI_UNSIGNED,
325 inner_face.first,
326 700 + my_domain,
327 &othersize,
328 1,
329 MPI_UNSIGNED,
330 inner_face.first,
331 700 + inner_face.first,
332 comm,
333 &status);
334 AssertThrowMPI(ierr);
335 AssertDimension(mysize, othersize);
336 mysize = inner_face.second.n_hanging_faces_larger_subdomain;
337 ierr = MPI_Sendrecv(&mysize,
338 1,
339 MPI_UNSIGNED,
340 inner_face.first,
341 800 + my_domain,
342 &othersize,
343 1,
344 MPI_UNSIGNED,
345 inner_face.first,
346 800 + inner_face.first,
347 comm,
348 &status);
349 AssertThrowMPI(ierr);
350 AssertDimension(mysize, othersize);
351# endif
352
353 // Arrange the face "ownership" such that cells that are access
354 // by more than one face (think of a cell in a corner) get
355 // ghosted. This arrangement has the advantage that we need to
356 // send less data because the same data is used twice. The
357 // strategy applied here is to ensure the same order of face
358 // pairs on both processors that share some faces, and make the
359 // same decision on both sides.
360
361 // Create a vector with cell ids sorted over the processor with
362 // the larger rank. In the code below we need to be able to
363 // identify the same cell once for the processor with higher
364 // rank and once for the processor with the lower rank. The
365 // format for the processor with the higher rank is already
366 // contained in `shared_faces`, whereas we need a copy that we
367 // sort differently for the other way around.
368 std::vector<std::tuple<CellId, CellId, unsigned int>> other_range(
369 inner_face.second.shared_faces.size());
370 for (unsigned int i = 0; i < other_range.size(); ++i)
371 other_range[i] =
372 std::make_tuple(inner_face.second.shared_faces[i].second,
373 inner_face.second.shared_faces[i].first,
374 i);
375 std::sort(other_range.begin(), other_range.end());
376
377 // the vector 'assignment' sets whether a particular cell
378 // appears more often and acts as a pre-selection of the rank. A
379 // value of 1 means that the process with the higher rank gets
380 // those faces, a value -1 means that the process with the lower
381 // rank gets it, whereas a value 0 means that the decision can
382 // be made in an arbitrary way.
383 unsigned int n_faces_lower_proc = 0, n_faces_higher_proc = 0;
384 std::vector<signed char> assignment(other_range.size(), 0);
385 if (inner_face.second.shared_faces.size() > 0)
386 {
387 // identify faces that go to the processor with the higher
388 // rank
389 unsigned int count = 0;
390 for (unsigned int i = 1;
391 i < inner_face.second.shared_faces.size();
392 ++i)
393 if (inner_face.second.shared_faces[i].first ==
394 inner_face.second.shared_faces[i - 1 - count].first)
395 ++count;
396 else
397 {
398 AssertThrow(count < 2 * dim, ExcInternalError());
399 if (count > 0)
400 {
401 for (unsigned int k = 0; k <= count; ++k)
402 assignment[i - 1 - k] = 1;
403 n_faces_higher_proc += count + 1;
404 }
405 count = 0;
406 }
407
408 // identify faces that definitely go to the processor with
409 // the lower rank - this must use the sorting of CellId
410 // variables from the processor with the higher rank, i.e.,
411 // other_range rather than `shared_faces`.
412 count = 0;
413 for (unsigned int i = 1; i < other_range.size(); ++i)
414 if (std::get<0>(other_range[i]) ==
415 std::get<0>(other_range[i - 1 - count]))
416 ++count;
417 else
418 {
419 AssertThrow(count < 2 * dim, ExcInternalError());
420 if (count > 0)
421 {
422 for (unsigned int k = 0; k <= count; ++k)
423 {
424 Assert(inner_face.second
425 .shared_faces[std::get<2>(
426 other_range[i - 1])]
427 .second ==
428 inner_face.second
429 .shared_faces[std::get<2>(
430 other_range[i - 1 - k])]
431 .second,
432 ExcInternalError());
433 // only assign to -1 if higher rank was not
434 // yet set
435 if (assignment[std::get<2>(
436 other_range[i - 1 - k])] == 0)
437 {
438 assignment[std::get<2>(
439 other_range[i - 1 - k])] = -1;
440 ++n_faces_lower_proc;
441 }
442 }
443 }
444 count = 0;
445 }
446 }
447
448
449 // divide the faces evenly between the two processors. the
450 // processor with small rank takes the first half, the processor
451 // with larger rank the second half. Adjust for the hanging
452 // faces that always get assigned to one side, and the faces we
453 // have already assigned due to the criterion above
454 n_faces_lower_proc +=
455 inner_face.second.n_hanging_faces_smaller_subdomain;
456 n_faces_higher_proc +=
457 inner_face.second.n_hanging_faces_larger_subdomain;
458 const unsigned int n_total_faces_at_proc_boundary =
459 (inner_face.second.shared_faces.size() +
460 inner_face.second.n_hanging_faces_smaller_subdomain +
461 inner_face.second.n_hanging_faces_larger_subdomain);
462 unsigned int split_index = n_total_faces_at_proc_boundary / 2;
463 if (split_index < n_faces_lower_proc)
464 split_index = 0;
465 else if (split_index <
466 n_total_faces_at_proc_boundary - n_faces_higher_proc)
467 split_index -= n_faces_lower_proc;
468 else
469 split_index = n_total_faces_at_proc_boundary -
470 n_faces_higher_proc - n_faces_lower_proc;
471
472 // make sure the splitting is consistent between both sides
473# if defined(DEAL_II_WITH_MPI) && defined(DEBUG)
474 ierr = MPI_Sendrecv(&split_index,
475 1,
476 MPI_UNSIGNED,
477 inner_face.first,
478 900 + my_domain,
479 &othersize,
480 1,
481 MPI_UNSIGNED,
482 inner_face.first,
483 900 + inner_face.first,
484 comm,
485 &status);
486 AssertThrowMPI(ierr);
487 AssertDimension(split_index, othersize);
488 ierr = MPI_Sendrecv(&n_faces_lower_proc,
489 1,
490 MPI_UNSIGNED,
491 inner_face.first,
492 1000 + my_domain,
493 &othersize,
494 1,
495 MPI_UNSIGNED,
496 inner_face.first,
497 1000 + inner_face.first,
498 comm,
499 &status);
500 AssertThrowMPI(ierr);
501 AssertDimension(n_faces_lower_proc, othersize);
502 ierr = MPI_Sendrecv(&n_faces_higher_proc,
503 1,
504 MPI_UNSIGNED,
505 inner_face.first,
506 1100 + my_domain,
507 &othersize,
508 1,
509 MPI_UNSIGNED,
510 inner_face.first,
511 1100 + inner_face.first,
512 comm,
513 &status);
514 AssertThrowMPI(ierr);
515 AssertDimension(n_faces_higher_proc, othersize);
516# endif
517
518 // collect the faces on both sides
519 std::vector<std::pair<CellId, CellId>> owned_faces_lower,
520 owned_faces_higher;
521 for (unsigned int i = 0; i < assignment.size(); ++i)
522 if (assignment[i] < 0)
523 owned_faces_lower.push_back(
524 inner_face.second.shared_faces[i]);
525 else if (assignment[i] > 0)
526 owned_faces_higher.push_back(
527 inner_face.second.shared_faces[i]);
528 AssertIndexRange(split_index,
529 inner_face.second.shared_faces.size() + 1 -
530 owned_faces_lower.size() -
531 owned_faces_higher.size());
532
533 unsigned int i = 0, c = 0;
534 for (; i < assignment.size() && c < split_index; ++i)
535 if (assignment[i] == 0)
536 {
537 owned_faces_lower.push_back(
538 inner_face.second.shared_faces[i]);
539 ++c;
540 }
541 for (; i < assignment.size(); ++i)
542 if (assignment[i] == 0)
543 {
544 owned_faces_higher.push_back(
545 inner_face.second.shared_faces[i]);
546 }
547
548 if constexpr (running_in_debug_mode())
549 {
550 // check consistency of faces on both sides
551 std::vector<std::pair<CellId, CellId>> check_faces;
552 check_faces.insert(check_faces.end(),
553 owned_faces_lower.begin(),
554 owned_faces_lower.end());
555 check_faces.insert(check_faces.end(),
556 owned_faces_higher.begin(),
557 owned_faces_higher.end());
558 std::sort(check_faces.begin(), check_faces.end());
559 AssertDimension(check_faces.size(),
560 inner_face.second.shared_faces.size());
561 for (unsigned int i = 0; i < check_faces.size(); ++i)
562 Assert(check_faces[i] == inner_face.second.shared_faces[i],
563 ExcInternalError());
564 }
565
566 // now only set half of the faces as the ones to keep
567 if (my_domain < inner_face.first)
568 inner_face.second.shared_faces.swap(owned_faces_lower);
569 else
570 inner_face.second.shared_faces.swap(owned_faces_higher);
571
572 std::sort(inner_face.second.shared_faces.begin(),
573 inner_face.second.shared_faces.end());
574 }
575 }
576
577 // fill in the additional cells that we need access to via ghosting to
578 // cell_levels
579 std::set<std::pair<unsigned int, unsigned int>> ghost_cells;
580 for (unsigned int i = 0; i < cell_levels.size(); ++i)
581 {
582 typename ::Triangulation<dim>::cell_iterator dcell(
583 &triangulation, cell_levels[i].first, cell_levels[i].second);
584 if (use_active_cells)
585 Assert(dcell->is_active(), ExcNotImplemented());
586 for (const auto f : dcell->face_indices())
587 {
588 if (dcell->at_boundary(f) && !dcell->has_periodic_neighbor(f))
589 face_is_owned[dcell->face(f)->index()] =
590 FaceCategory::locally_active_at_boundary;
591 else if (!build_inner_faces)
592 continue;
593
594 // treat boundaries of cells of different refinement level
595 // inside the domain in case of multigrid separately
596 else if ((dcell->at_boundary(f) == false ||
597 dcell->has_periodic_neighbor(f)) &&
598 mg_level != numbers::invalid_unsigned_int &&
599 dcell->neighbor_or_periodic_neighbor(f)->level() <
600 dcell->level())
601 {
602 face_is_owned[dcell->face(f)->index()] =
603 FaceCategory::multigrid_refinement_edge;
604 }
605 else
606 {
607 typename ::Triangulation<dim>::cell_iterator neighbor =
608 dcell->neighbor_or_periodic_neighbor(f);
609
610 // neighbor is refined -> face will be treated by neighbor
611 if (use_active_cells && neighbor->has_children() &&
612 hold_all_faces_to_owned_cells == false)
613 continue;
614
615 bool add_to_ghost = false;
616 const types::subdomain_id
617 id1 = use_active_cells ? dcell->subdomain_id() :
618 dcell->level_subdomain_id(),
619 id2 = use_active_cells ?
620 (neighbor->has_children() ?
621 dcell->neighbor_child_on_subface(f, 0)
622 ->subdomain_id() :
623 neighbor->subdomain_id()) :
624 neighbor->level_subdomain_id();
625
626 // Check whether the current face should be processed
627 // locally (instead of being processed from the other
628 // side). We process a face locally when we are more refined
629 // (in the active cell case) or when the face is listed in
630 // the `shared_faces` data structure that we built above.
631 if ((id1 == id2 &&
632 (use_active_cells == false || neighbor->is_active())) ||
633 dcell->level() > neighbor->level() ||
634 std::binary_search(
635 inner_faces_at_proc_boundary[id2].shared_faces.begin(),
636 inner_faces_at_proc_boundary[id2].shared_faces.end(),
637 std::make_pair(id1 < id2 ? dcell->id() : neighbor->id(),
638 id1 < id2 ? neighbor->id() :
639 dcell->id())))
640 {
641 face_is_owned[dcell->face(f)->index()] =
642 FaceCategory::locally_active_done_here;
643 if (dcell->level() == neighbor->level() ||
644 dcell->has_periodic_neighbor(f))
645 face_is_owned
646 [neighbor
647 ->face(dcell->has_periodic_neighbor(f) ?
648 dcell->periodic_neighbor_face_no(f) :
649 dcell->neighbor_face_no(f))
650 ->index()] =
651 FaceCategory::locally_active_done_here;
652
653 // If neighbor is a ghost element (i.e.
654 // dcell->subdomain_id !
655 // dcell->neighbor(f)->subdomain_id()), we need to add its
656 // index into cell level list.
657 if (use_active_cells)
658 add_to_ghost =
659 (dcell->subdomain_id() != neighbor->subdomain_id());
660 else
661 add_to_ghost = (dcell->level_subdomain_id() !=
662 neighbor->level_subdomain_id());
663 }
664 else if (hold_all_faces_to_owned_cells == true)
665 {
666 // add all cells to ghost layer...
667 face_is_owned[dcell->face(f)->index()] =
668 FaceCategory::ghosted;
669 if (use_active_cells)
670 {
671 if (neighbor->has_children())
672 for (unsigned int s = 0;
673 s < dcell->face(f)->n_children();
674 ++s)
675 if (dcell->at_boundary(f))
676 {
677 if (dcell
678 ->periodic_neighbor_child_on_subface(f,
679 s)
680 ->subdomain_id() !=
681 dcell->subdomain_id())
682 add_to_ghost = true;
683 }
684 else
685 {
686 if (dcell->neighbor_child_on_subface(f, s)
687 ->subdomain_id() !=
688 dcell->subdomain_id())
689 add_to_ghost = true;
690 }
691 else
692 add_to_ghost = (dcell->subdomain_id() !=
693 neighbor->subdomain_id());
694 }
695 else
696 add_to_ghost = (dcell->level_subdomain_id() !=
697 neighbor->level_subdomain_id());
698 }
699
700 if (add_to_ghost)
701 {
702 if (use_active_cells && neighbor->has_children())
703 for (unsigned int s = 0;
704 s < dcell->face(f)->n_children();
705 ++s)
706 {
707 typename ::Triangulation<dim>::cell_iterator
708 neighbor_child =
709 dcell->at_boundary(f) ?
710 dcell->periodic_neighbor_child_on_subface(f,
711 s) :
712 dcell->neighbor_child_on_subface(f, s);
713 if (neighbor_child->subdomain_id() !=
714 dcell->subdomain_id())
715 ghost_cells.insert(
716 std::pair<unsigned int, unsigned int>(
717 neighbor_child->level(),
718 neighbor_child->index()));
719 }
720 else
721 ghost_cells.insert(
722 std::pair<unsigned int, unsigned int>(
723 neighbor->level(), neighbor->index()));
724 at_processor_boundary[i] = true;
725 }
726 }
727 }
728 }
729
730 // step 2: append the ghost cells at the end of the locally owned
731 // cells
732 for (const auto &ghost_cell : ghost_cells)
733 cell_levels.push_back(ghost_cell);
734 }
735
736
737
738 template <int dim>
739 void
740 FaceSetup<dim>::generate_faces(
741 const ::Triangulation<dim> &triangulation,
742 const std::vector<std::pair<unsigned int, unsigned int>> &cell_levels,
743 TaskInfo &task_info)
744 {
745 const bool is_mixed_mesh = triangulation.is_mixed_mesh();
746
747 // step 1: create the inverse map between cell iterators and the
748 // cell_level_index field
749 std::map<std::pair<unsigned int, unsigned int>, unsigned int>
750 map_to_vectorized;
751 for (unsigned int cell = 0; cell < cell_levels.size(); ++cell)
752 if (cell == 0 || cell_levels[cell] != cell_levels[cell - 1])
753 {
754 typename ::Triangulation<dim>::cell_iterator dcell(
755 &triangulation,
756 cell_levels[cell].first,
757 cell_levels[cell].second);
758 std::pair<unsigned int, unsigned int> level_index(dcell->level(),
759 dcell->index());
760 map_to_vectorized[level_index] = cell;
761 }
762
763 // step 2: fill the information about inner faces and boundary faces
764 const unsigned int vectorization_length = task_info.vectorization_length;
765 task_info.face_partition_data.resize(
766 task_info.cell_partition_data.size() - 1, 0);
767 task_info.boundary_partition_data.resize(
768 task_info.cell_partition_data.size() - 1, 0);
769 std::vector<unsigned char> face_visited(face_is_owned.size(), 0);
770 for (unsigned int partition = 0;
771 partition < task_info.cell_partition_data.size() - 2;
772 ++partition)
773 {
774 unsigned int boundary_counter = 0;
775 unsigned int inner_counter = 0;
776 for (unsigned int cell = task_info.cell_partition_data[partition] *
777 vectorization_length;
778 cell < task_info.cell_partition_data[partition + 1] *
779 vectorization_length;
780 ++cell)
781 if (cell == 0 || cell_levels[cell] != cell_levels[cell - 1])
782 {
783 typename ::Triangulation<dim>::cell_iterator dcell(
784 &triangulation,
785 cell_levels[cell].first,
786 cell_levels[cell].second);
787 for (const auto f : dcell->face_indices())
788 {
789 // boundary face
790 if (face_is_owned[dcell->face(f)->index()] ==
791 FaceCategory::locally_active_at_boundary)
792 {
793 Assert(dcell->at_boundary(f), ExcInternalError());
794 ++boundary_counter;
795 FaceToCellTopology<1> info;
796 info.cells_interior[0] = cell;
797 info.cells_exterior[0] = numbers::invalid_unsigned_int;
798 info.interior_face_no = f;
799 info.exterior_face_no = dcell->face(f)->boundary_id();
800 info.face_type =
801 is_mixed_mesh ?
802 (dcell->face(f)->reference_cell() !=
803 ReferenceCells::get_hypercube<dim - 1>()) :
804 0;
805 info.subface_index =
806 GeometryInfo<dim>::max_children_per_cell;
807 info.face_orientation = 0;
808 boundary_faces.push_back(info);
809
810 face_visited[dcell->face(f)->index()]++;
811 }
812 // interior face, including faces over periodic boundaries
813 else
814 {
815 typename ::Triangulation<dim>::cell_iterator
816 neighbor = dcell->neighbor_or_periodic_neighbor(f);
817 if (use_active_cells && neighbor->has_children())
818 {
819 const unsigned int n_children =
820 (dim == 1) ? 1 : dcell->face(f)->n_children();
821 for (unsigned int c = 0; c < n_children; ++c)
822 {
823 typename ::Triangulation<
824 dim>::cell_iterator neighbor_c;
825 if (dim > 1)
826 neighbor_c =
827 (dcell->at_boundary(f) ?
828 dcell
829 ->periodic_neighbor_child_on_subface(
830 f, c) :
831 dcell->neighbor_child_on_subface(f, c));
832 else
833 {
834 // in 1D, adjacent cells can differ by
835 // more than 1 level
836 neighbor_c = neighbor->child(1 - f);
837 while (!neighbor_c->is_active())
838 neighbor_c = neighbor_c->child(1 - f);
839 }
840 const types::subdomain_id neigh_domain =
841 neighbor_c->subdomain_id();
842 const unsigned int neighbor_face_no =
843 dcell->has_periodic_neighbor(f) ?
844 dcell->periodic_neighbor_face_no(f) :
845 dcell->neighbor_face_no(f);
846 const unsigned int child_face_index =
847 dim > 1 ? dcell->face(f)->child(c)->index() :
848 dcell->face(f)->index();
849 if (neigh_domain != dcell->subdomain_id() ||
850 face_visited[child_face_index] == 1)
851 {
852 std::pair<unsigned int, unsigned int>
853 level_index(neighbor_c->level(),
854 neighbor_c->index());
855 if (face_is_owned[child_face_index] ==
856 FaceCategory::locally_active_done_here)
857 {
858 ++inner_counter;
859 inner_faces.push_back(create_face(
860 neighbor_face_no,
861 neighbor_c,
862 map_to_vectorized[level_index],
863 dcell,
864 cell,
865 is_mixed_mesh));
866 }
867 else if (face_is_owned[child_face_index] ==
868 FaceCategory::ghosted ||
869 face_is_owned[dcell->face(f)
870 ->index()] ==
871 FaceCategory::ghosted)
872 {
873 inner_ghost_faces.push_back(create_face(
874 neighbor_face_no,
875 neighbor_c,
876 map_to_vectorized[level_index],
877 dcell,
878 cell,
879 is_mixed_mesh));
880 }
881 else
882 Assert(face_is_owned[dcell->face(f)
883 ->index()] ==
884 FaceCategory::
885 locally_active_done_elsewhere,
886 ExcInternalError());
887 }
888 else
889 {
890 face_visited[child_face_index] = 1;
891 if (dcell->has_periodic_neighbor(f))
892 face_visited[dcell->face(f)->index()] = 1;
893 }
894 }
895 }
896 else
897 {
898 const types::subdomain_id my_domain =
899 use_active_cells ? dcell->subdomain_id() :
900 dcell->level_subdomain_id();
901 const types::subdomain_id neigh_domain =
902 use_active_cells ? neighbor->subdomain_id() :
903 neighbor->level_subdomain_id();
904 const unsigned int face_index =
905 dcell->face(f)->index();
906 if (neigh_domain != my_domain ||
907 face_visited[face_index] == 1)
908 {
909 std::pair<unsigned int, unsigned int>
910 level_index(neighbor->level(),
911 neighbor->index());
912 if (face_is_owned[dcell->face(f)->index()] ==
913 FaceCategory::locally_active_done_here)
914 {
915 Assert(use_active_cells ||
916 dcell->level() ==
917 neighbor->level(),
918 ExcInternalError());
919 ++inner_counter;
920 inner_faces.push_back(create_face(
921 f,
922 dcell,
923 cell,
924 neighbor,
925 map_to_vectorized[level_index],
926 is_mixed_mesh));
927 }
928 else if (face_is_owned[face_index] ==
929 FaceCategory::ghosted)
930 {
931 inner_ghost_faces.push_back(create_face(
932 f,
933 dcell,
934 cell,
935 neighbor,
936 map_to_vectorized[level_index],
937 is_mixed_mesh));
938 }
939 }
940 else
941 {
942 face_visited[face_index] = 1;
943 if (dcell->has_periodic_neighbor(f))
944 face_visited
945 [neighbor
946 ->face(
947 dcell->periodic_neighbor_face_no(f))
948 ->index()] = 1;
949 }
950 if (face_is_owned[face_index] ==
951 FaceCategory::multigrid_refinement_edge)
952 {
953 refinement_edge_faces.push_back(
954 create_face(f,
955 dcell,
956 cell,
957 neighbor,
958 refinement_edge_faces.size(),
959 is_mixed_mesh));
960 }
961 }
962 }
963 }
964 }
965 task_info.face_partition_data[partition + 1] =
966 task_info.face_partition_data[partition] + inner_counter;
967 task_info.boundary_partition_data[partition + 1] =
968 task_info.boundary_partition_data[partition] + boundary_counter;
969 }
970 task_info.ghost_face_partition_data.resize(2);
971 task_info.ghost_face_partition_data[0] = 0;
972 task_info.ghost_face_partition_data[1] = inner_ghost_faces.size();
973 task_info.refinement_edge_face_partition_data.resize(2);
974 task_info.refinement_edge_face_partition_data[0] = 0;
975 task_info.refinement_edge_face_partition_data[1] =
976 refinement_edge_faces.size();
977 }
978
979
980
981 template <int dim>
982 FaceToCellTopology<1>
983 FaceSetup<dim>::create_face(
984 const unsigned int face_no,
985 const typename ::Triangulation<dim>::cell_iterator &cell,
986 const unsigned int number_cell_interior,
987 const typename ::Triangulation<dim>::cell_iterator &neighbor,
988 const unsigned int number_cell_exterior,
989 const bool is_mixed_mesh)
990 {
991 FaceToCellTopology<1> info;
992 info.cells_interior[0] = number_cell_interior;
993 info.cells_exterior[0] = number_cell_exterior;
994 info.interior_face_no = face_no;
995 if (cell->has_periodic_neighbor(face_no))
996 info.exterior_face_no = cell->periodic_neighbor_face_no(face_no);
997 else
998 info.exterior_face_no = cell->neighbor_face_no(face_no);
999
1000 info.face_type = is_mixed_mesh ?
1001 (cell->face(face_no)->reference_cell() !=
1002 ReferenceCells::get_hypercube<dim - 1>()) :
1003 0;
1004
1005 info.subface_index = GeometryInfo<dim>::max_children_per_cell;
1006 Assert(neighbor->level() <= cell->level(), ExcInternalError());
1007
1008 // for dim > 1 and hanging faces we must set a subface index
1009 if (dim > 1 && cell->level() > neighbor->level())
1010 {
1011 if (cell->has_periodic_neighbor(face_no))
1012 info.subface_index =
1013 cell->periodic_neighbor_of_coarser_periodic_neighbor(face_no)
1014 .second;
1015 else
1016 info.subface_index =
1017 cell->neighbor_of_coarser_neighbor(face_no).second;
1018 }
1019
1020 // special treatment of periodic boundaries
1021 if (cell->has_periodic_neighbor(face_no))
1022 {
1023 info.face_orientation = cell->get_triangulation()
1024 .get_periodic_face_map()
1025 .at({cell, face_no})
1026 .second;
1027 }
1028 else
1029 {
1030 const auto interior_face_orientation =
1031 cell->combined_face_orientation(face_no);
1032 const auto exterior_face_orientation =
1033 neighbor->combined_face_orientation(info.exterior_face_no);
1034 if (interior_face_orientation !=
1035 numbers::default_geometric_orientation)
1036 {
1037 info.face_orientation = 8 + interior_face_orientation;
1038 Assert(exterior_face_orientation ==
1039 numbers::default_geometric_orientation,
1040 ExcMessage(
1041 "Face seems to be wrongly oriented from both sides"));
1042 }
1043 else
1044 info.face_orientation = exterior_face_orientation;
1045
1046 // make sure to select correct subface index in case of non-standard
1047 // orientation of the coarser neighbor face
1048 if (cell->level() > neighbor->level() &&
1049 exterior_face_orientation > 0)
1050 {
1051 const Table<2, unsigned int> orientation =
1052 ShapeInfo<double>::compute_orientation_table(2);
1053 const auto face_reference_cell =
1054 cell->face(face_no)->reference_cell();
1055 info.subface_index = orientation(
1056 face_reference_cell.get_inverse_combined_orientation(
1057 exterior_face_orientation),
1058 info.subface_index);
1059 }
1060 }
1061
1062 return info;
1063 }
1064
1065
1066
1073 inline bool
1074 compare_faces_for_vectorization(
1075 const FaceToCellTopology<1> &face1,
1076 const FaceToCellTopology<1> &face2,
1077 const std::vector<unsigned int> &active_fe_indices,
1078 const unsigned int length)
1079 {
1080 if (face1.interior_face_no != face2.interior_face_no)
1081 return false;
1082 if (face1.exterior_face_no != face2.exterior_face_no)
1083 return false;
1084 if (face1.subface_index != face2.subface_index)
1085 return false;
1086 if (face1.face_orientation != face2.face_orientation)
1087 return false;
1088 if (face1.face_type != face2.face_type)
1089 return false;
1090
1091 if (active_fe_indices.size() > 0)
1092 {
1093 if (active_fe_indices[face1.cells_interior[0] / length] !=
1094 active_fe_indices[face2.cells_interior[0] / length])
1095 return false;
1096
1097 if (face2.cells_exterior[0] != numbers::invalid_unsigned_int)
1098 if (active_fe_indices[face1.cells_exterior[0] / length] !=
1099 active_fe_indices[face2.cells_exterior[0] / length])
1100 return false;
1101 }
1102
1103 return true;
1104 }
1105
1106
1107
1114 template <int length>
1115 struct FaceComparator
1116 {
1117 FaceComparator(const std::vector<unsigned int> &active_fe_indices)
1118 : active_fe_indices(active_fe_indices)
1119 {}
1120
1121 bool
1122 operator()(const FaceToCellTopology<length> &face1,
1123 const FaceToCellTopology<length> &face2) const
1124 {
1125 // check if active FE indices match
1126 if (face1.face_type < face2.face_type)
1127 return true;
1128 else if (face1.face_type > face2.face_type)
1129 return false;
1130
1131 // check if active FE indices match
1132 if (active_fe_indices.size() > 0)
1133 {
1134 // ... for interior faces
1135 if (active_fe_indices[face1.cells_interior[0] / length] <
1136 active_fe_indices[face2.cells_interior[0] / length])
1137 return true;
1138 else if (active_fe_indices[face1.cells_interior[0] / length] >
1139 active_fe_indices[face2.cells_interior[0] / length])
1140 return false;
1141
1142 // ... for exterior faces
1143 if (face2.cells_exterior[0] != numbers::invalid_unsigned_int)
1144 {
1145 if (active_fe_indices[face1.cells_exterior[0] / length] <
1146 active_fe_indices[face2.cells_exterior[0] / length])
1147 return true;
1148 else if (active_fe_indices[face1.cells_exterior[0] / length] >
1149 active_fe_indices[face2.cells_exterior[0] / length])
1150 return false;
1151 }
1152 }
1153
1154 for (unsigned int i = 0; i < length; ++i)
1155 if (face1.cells_interior[i] < face2.cells_interior[i])
1156 return true;
1157 else if (face1.cells_interior[i] > face2.cells_interior[i])
1158 return false;
1159 for (unsigned int i = 0; i < length; ++i)
1160 if (face1.cells_exterior[i] < face2.cells_exterior[i])
1161 return true;
1162 else if (face1.cells_exterior[i] > face2.cells_exterior[i])
1163 return false;
1164 if (face1.interior_face_no < face2.interior_face_no)
1165 return true;
1166 else if (face1.interior_face_no > face2.interior_face_no)
1167 return false;
1168 if (face1.exterior_face_no < face2.exterior_face_no)
1169 return true;
1170 else if (face1.exterior_face_no > face2.exterior_face_no)
1171 return false;
1172
1173 // we do not need to check for subface_index and orientation because
1174 // those cannot be different if when all the other values are the
1175 // same.
1176 AssertDimension(face1.subface_index, face2.subface_index);
1177 AssertDimension(face1.face_orientation, face2.face_orientation);
1178
1179 return false;
1180 }
1181
1182 private:
1183 const std::vector<unsigned int> &active_fe_indices;
1184 };
1185
1186
1187
1188 template <int vectorization_width>
1189 void
1190 collect_faces_vectorization(
1191 const std::vector<FaceToCellTopology<1>> &faces_in,
1192 const std::vector<bool> &hard_vectorization_boundary,
1193 std::vector<unsigned int> &face_partition_data,
1194 std::vector<FaceToCellTopology<vectorization_width>> &faces_out,
1195 const std::vector<unsigned int> &active_fe_indices)
1196 {
1197 FaceToCellTopology<vectorization_width> face_batch;
1198 std::vector<std::vector<unsigned int>> faces_type;
1199
1200 unsigned int face_start = face_partition_data[0],
1201 face_end = face_partition_data[0];
1202
1203 face_partition_data[0] = faces_out.size();
1204 for (unsigned int partition = 0;
1205 partition < face_partition_data.size() - 1;
1206 ++partition)
1207 {
1208 std::vector<std::vector<unsigned int>> new_faces_type;
1209
1210 // start with the end point for the last partition
1211 face_start = face_end;
1212 face_end = face_partition_data[partition + 1];
1213
1214 // set the partitioner to the new vectorized lengths
1215 face_partition_data[partition + 1] = face_partition_data[partition];
1216
1217 // loop over the faces in the current partition and reorder according
1218 // to the face type
1219 for (unsigned int face = face_start; face < face_end; ++face)
1220 {
1221 for (auto &face_type : faces_type)
1222 {
1223 // Compare current face with first face of type type
1224 if (compare_faces_for_vectorization(faces_in[face],
1225 faces_in[face_type[0]],
1226 active_fe_indices,
1227 vectorization_width))
1228 {
1229 face_type.push_back(face);
1230 goto face_found;
1231 }
1232 }
1233 faces_type.emplace_back(1, face);
1234 face_found:
1235 {}
1236 }
1237
1238 // insert new faces in sorted list to get good data locality
1239 FaceComparator<vectorization_width> face_comparator(
1240 active_fe_indices);
1241 std::set<FaceToCellTopology<vectorization_width>,
1242 FaceComparator<vectorization_width>>
1243 new_faces(face_comparator);
1244 for (const auto &face_type : faces_type)
1245 {
1246 face_batch.face_type = faces_in[face_type[0]].face_type;
1247 face_batch.interior_face_no =
1248 faces_in[face_type[0]].interior_face_no;
1249 face_batch.exterior_face_no =
1250 faces_in[face_type[0]].exterior_face_no;
1251 face_batch.subface_index = faces_in[face_type[0]].subface_index;
1252 face_batch.face_orientation =
1253 faces_in[face_type[0]].face_orientation;
1254 unsigned int no_faces = face_type.size();
1255 std::vector<unsigned char> touched(no_faces, 0);
1256
1257 // do two passes through the data. The first is to identify
1258 // similar faces within the same index range as the cells which
1259 // will allow for vectorized read operations, the second picks up
1260 // all the rest
1261 unsigned int n_vectorized = 0;
1262 for (unsigned int f = 0; f < no_faces; ++f)
1263 if (faces_in[face_type[f]].cells_interior[0] %
1264 vectorization_width ==
1265 0)
1266 {
1267 bool is_contiguous = true;
1268 if (f + vectorization_width > no_faces)
1269 is_contiguous = false;
1270 else
1271 for (unsigned int v = 1; v < vectorization_width; ++v)
1272 if (faces_in[face_type[f + v]].cells_interior[0] !=
1273 faces_in[face_type[f]].cells_interior[0] + v)
1274 is_contiguous = false;
1275 if (is_contiguous)
1276 {
1277 AssertIndexRange(f,
1278 face_type.size() -
1279 vectorization_width + 1);
1280 for (unsigned int v = 0; v < vectorization_width; ++v)
1281 {
1282 face_batch.cells_interior[v] =
1283 faces_in[face_type[f + v]].cells_interior[0];
1284 face_batch.cells_exterior[v] =
1285 faces_in[face_type[f + v]].cells_exterior[0];
1286 touched[f + v] = 1;
1287 }
1288 new_faces.insert(face_batch);
1289 f += vectorization_width - 1;
1290 n_vectorized += vectorization_width;
1291 }
1292 }
1293
1294 std::vector<unsigned int> untouched;
1295 untouched.reserve(no_faces - n_vectorized);
1296 for (unsigned int f = 0; f < no_faces; ++f)
1297 if (touched[f] == 0)
1298 untouched.push_back(f);
1299 unsigned int v = 0;
1300 for (const auto f : untouched)
1301 {
1302 face_batch.cells_interior[v] =
1303 faces_in[face_type[f]].cells_interior[0];
1304 face_batch.cells_exterior[v] =
1305 faces_in[face_type[f]].cells_exterior[0];
1306 ++v;
1307 if (v == vectorization_width)
1308 {
1309 new_faces.insert(face_batch);
1310 v = 0;
1311 }
1312 }
1313 if (v > 0 && v < vectorization_width)
1314 {
1315 // must add non-filled face
1316 if (hard_vectorization_boundary[partition + 1] ||
1317 partition == face_partition_data.size() - 2)
1318 {
1319 for (; v < vectorization_width; ++v)
1320 {
1321 // Dummy cell, not used
1322 face_batch.cells_interior[v] =
1323 numbers::invalid_unsigned_int;
1324 face_batch.cells_exterior[v] =
1325 numbers::invalid_unsigned_int;
1326 }
1327 new_faces.insert(face_batch);
1328 }
1329 else
1330 {
1331 // postpone to the next partition
1332 std::vector<unsigned int> untreated(v);
1333 for (unsigned int f = 0; f < v; ++f)
1334 untreated[f] = face_type[*(untouched.end() - 1 - f)];
1335 new_faces_type.push_back(untreated);
1336 }
1337 }
1338 }
1339
1340 // insert sorted list to vector of faces
1341 for (auto it = new_faces.begin(); it != new_faces.end(); ++it)
1342 faces_out.push_back(*it);
1343 face_partition_data[partition + 1] += new_faces.size();
1344
1345 // set the faces that were left over to faces_type for the next round
1346 faces_type = std::move(new_faces_type);
1347 }
1348
1349 if constexpr (running_in_debug_mode())
1350 {
1351 // final safety checks
1352 for (const auto &face_type : faces_type)
1353 AssertDimension(face_type.size(), 0U);
1354
1355 AssertDimension(faces_out.size(), face_partition_data.back());
1356 unsigned int nfaces = 0;
1357 for (unsigned int i = face_partition_data[0];
1358 i < face_partition_data.back();
1359 ++i)
1360 for (unsigned int v = 0; v < vectorization_width; ++v)
1361 nfaces += (faces_out[i].cells_interior[v] !=
1362 numbers::invalid_unsigned_int);
1363 AssertDimension(nfaces, faces_in.size());
1364
1365 std::vector<std::pair<unsigned int, unsigned int>> in_faces,
1366 out_faces;
1367 for (const auto &face_in : faces_in)
1368 in_faces.emplace_back(face_in.cells_interior[0],
1369 face_in.cells_exterior[0]);
1370 for (unsigned int i = face_partition_data[0];
1371 i < face_partition_data.back();
1372 ++i)
1373 for (unsigned int v = 0;
1374 v < vectorization_width && faces_out[i].cells_interior[v] !=
1375 numbers::invalid_unsigned_int;
1376 ++v)
1377 out_faces.emplace_back(faces_out[i].cells_interior[v],
1378 faces_out[i].cells_exterior[v]);
1379 std::sort(in_faces.begin(), in_faces.end());
1380 std::sort(out_faces.begin(), out_faces.end());
1381 AssertDimension(in_faces.size(), out_faces.size());
1382 for (unsigned int i = 0; i < in_faces.size(); ++i)
1383 {
1384 AssertDimension(in_faces[i].first, out_faces[i].first);
1385 AssertDimension(in_faces[i].second, out_faces[i].second);
1386 }
1387 }
1388 }
1389
1390#endif // ifndef DOXYGEN
1391
1392 } // namespace MatrixFreeFunctions
1393} // namespace internal
1394
1395
1396DEAL_II_NAMESPACE_CLOSE
1397
1398#endif
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:40
running_in_debug_mode
constexpr bool running_in_debug_mode()
Definition config.h:78
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:41
face_info.h
StandardExceptions::ExcNotImplemented
static ::ExceptionBase & ExcNotImplemented()
Assert
#define Assert(cond, exc)
Definition exception_macros.h:559
AssertDimension
#define AssertDimension(dim1, dim2)
Definition exception_macros.h:851
AssertThrowMPI
#define AssertThrowMPI(error_code)
Definition exception_macros.h:982
AssertIndexRange
#define AssertIndexRange(index, range)
Definition exception_macros.h:923
StandardExceptions::ExcInternalError
static ::ExceptionBase & ExcInternalError()
StandardExceptions::ExcMessage
static ::ExceptionBase & ExcMessage(std::string arg1)
AssertThrow
#define AssertThrow(cond, exc)
Definition exception_macros.h:656
cell_iterator
typename ActiveSelector::cell_iterator cell_iterator
Definition dof_handler.h:476
memory_consumption.h
GridGenerator::reference_cell
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
ReferenceCells::get_hypercube
constexpr const ReferenceCell & get_hypercube()
Definition reference_cell.h:2804
SparsityTools::partition
void partition(const SparsityPattern &sparsity_pattern, const unsigned int n_partitions, std::vector< unsigned int > &partition_indices, const Partitioner partitioner=Partitioner::metis)
Definition sparsity_tools.cc:396
internal::MatrixFreeFunctions::collect_faces_vectorization
void collect_faces_vectorization(const std::vector< FaceToCellTopology< 1 > > &faces_in, const std::vector< bool > &hard_vectorization_boundary, std::vector< unsigned int > &face_partition_data, std::vector< FaceToCellTopology< vectorization_width > > &faces_out)
numbers::invalid_unsigned_int
constexpr unsigned int invalid_unsigned_int
Definition types.h:238
numbers::invalid_subdomain_id
constexpr types::subdomain_id invalid_subdomain_id
Definition types.h:381
numbers::default_geometric_orientation
constexpr types::geometric_orientation default_geometric_orientation
Definition types.h:352
types::subdomain_id
unsigned int subdomain_id
Definition types.h:52
shape_info.h
GeometryInfo::max_children_per_cell
static constexpr unsigned int max_children_per_cell
Definition geometry_info.h:1973
internal::MatrixFreeFunctions::FaceIdentifier::shared_faces
std::vector< std::pair< CellId, CellId > > shared_faces
Definition face_setup_internal.h:54
internal::MatrixFreeFunctions::FaceSetup::inner_faces
std::vector< FaceToCellTopology< 1 > > inner_faces
Definition face_setup_internal.h:132
internal::MatrixFreeFunctions::FaceSetup::boundary_faces
std::vector< FaceToCellTopology< 1 > > boundary_faces
Definition face_setup_internal.h:133
internal::MatrixFreeFunctions::FaceSetup::refinement_edge_faces
std::vector< FaceToCellTopology< 1 > > refinement_edge_faces
Definition face_setup_internal.h:135
internal::MatrixFreeFunctions::FaceSetup::create_face
FaceToCellTopology< 1 > create_face(const unsigned int face_no, const typename ::Triangulation< dim >::cell_iterator &cell, const unsigned int number_cell_interior, const typename ::Triangulation< dim >::cell_iterator &neighbor, const unsigned int number_cell_exterior, const bool is_mixed_mesh)
internal::MatrixFreeFunctions::FaceSetup::initialize
void initialize(const ::Triangulation< dim > &triangulation, const unsigned int mg_level, const bool hold_all_faces_to_owned_cells, const bool build_inner_faces, std::vector< std::pair< unsigned int, unsigned int > > &cell_levels)
internal::MatrixFreeFunctions::FaceSetup::inner_ghost_faces
std::vector< FaceToCellTopology< 1 > > inner_ghost_faces
Definition face_setup_internal.h:134
internal::MatrixFreeFunctions::FaceSetup::generate_faces
void generate_faces(const ::Triangulation< dim > &triangulation, const std::vector< std::pair< unsigned int, unsigned int > > &cell_levels, TaskInfo &task_info)
task_info.h
tria_accessor.h
tria_base.h