Reference documentation for deal.II version GIT b6bf1e606d 2022-08-11 15:25:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
face_setup_internal.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2018 - 2022 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 
17 #ifndef dealii_face_setup_internal_h
18 #define dealii_face_setup_internal_h
19 
20 #include <deal.II/base/config.h>
21 
23 #include <deal.II/base/utilities.h>
24 
26 
27 #include <deal.II/grid/tria.h>
29 
33 
34 #include <fstream>
35 
36 
38 
39 
40 namespace internal
41 {
42  namespace MatrixFreeFunctions
43  {
49  {
53  {}
54 
55  std::vector<std::pair<CellId, CellId>> shared_faces;
58  };
59 
60 
61 
70  template <int dim>
71  struct FaceSetup
72  {
74 
81  void
83  const ::Triangulation<dim> &triangulation,
84  const unsigned int mg_level,
85  const bool hold_all_faces_to_owned_cells,
86  const bool build_inner_faces,
87  std::vector<std::pair<unsigned int, unsigned int>> &cell_levels);
88 
95  void
97  const ::Triangulation<dim> & triangulation,
98  const std::vector<std::pair<unsigned int, unsigned int>> &cell_levels,
99  TaskInfo & task_info);
100 
109  const unsigned int face_no,
110  const typename ::Triangulation<dim>::cell_iterator &cell,
111  const unsigned int number_cell_interior,
112  const typename ::Triangulation<dim>::cell_iterator &neighbor,
113  const unsigned int number_cell_exterior,
114  const bool is_mixed_mesh);
115 
117 
122  enum class FaceCategory : char
123  {
127  ghosted,
129  };
130 
131  std::vector<FaceCategory> face_is_owned;
132  std::vector<bool> at_processor_boundary;
133  std::vector<FaceToCellTopology<1>> inner_faces;
134  std::vector<FaceToCellTopology<1>> boundary_faces;
135  std::vector<FaceToCellTopology<1>> inner_ghost_faces;
136  std::vector<FaceToCellTopology<1>> refinement_edge_faces;
137  };
138 
139 
140 
144  template <int vectorization_width>
145  void
147  const std::vector<FaceToCellTopology<1>> &faces_in,
148  const std::vector<bool> & hard_vectorization_boundary,
149  std::vector<unsigned int> & face_partition_data,
150  std::vector<FaceToCellTopology<vectorization_width>> &faces_out);
151 
152 
153 
154  /* -------------------------------------------------------------------- */
155 
156 #ifndef DOXYGEN
157 
158  template <int dim>
160  : use_active_cells(true)
161  {}
162 
163 
164 
165  template <int dim>
166  void
167  FaceSetup<dim>::initialize(
168  const ::Triangulation<dim> &triangulation,
169  const unsigned int mg_level,
170  const bool hold_all_faces_to_owned_cells,
171  const bool build_inner_faces,
172  std::vector<std::pair<unsigned int, unsigned int>> &cell_levels)
173  {
174  use_active_cells = mg_level == numbers::invalid_unsigned_int;
175 
176 # ifdef DEBUG
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 # endif
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() !=
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_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<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 # ifdef DEBUG
549  // check consistency of faces on both sides
550  std::vector<std::pair<CellId, CellId>> check_faces;
551  check_faces.insert(check_faces.end(),
552  owned_faces_lower.begin(),
553  owned_faces_lower.end());
554  check_faces.insert(check_faces.end(),
555  owned_faces_higher.begin(),
556  owned_faces_higher.end());
557  std::sort(check_faces.begin(), check_faces.end());
558  AssertDimension(check_faces.size(),
559  inner_face.second.shared_faces.size());
560  for (unsigned int i = 0; i < check_faces.size(); ++i)
561  Assert(check_faces[i] == inner_face.second.shared_faces[i],
562  ExcInternalError());
563 # endif
564 
565  // now only set half of the faces as the ones to keep
566  if (my_domain < inner_face.first)
567  inner_face.second.shared_faces.swap(owned_faces_lower);
568  else
569  inner_face.second.shared_faces.swap(owned_faces_higher);
570 
571  std::sort(inner_face.second.shared_faces.begin(),
572  inner_face.second.shared_faces.end());
573  }
574  }
575 
576  // fill in the additional cells that we need access to via ghosting to
577  // cell_levels
578  std::set<std::pair<unsigned int, unsigned int>> ghost_cells;
579  for (unsigned int i = 0; i < cell_levels.size(); ++i)
580  {
581  typename ::Triangulation<dim>::cell_iterator dcell(
582  &triangulation, cell_levels[i].first, cell_levels[i].second);
583  if (use_active_cells)
584  Assert(dcell->is_active(), ExcNotImplemented());
585  for (const auto f : dcell->face_indices())
586  {
587  if (dcell->at_boundary(f) && !dcell->has_periodic_neighbor(f))
588  face_is_owned[dcell->face(f)->index()] =
589  FaceCategory::locally_active_at_boundary;
590  else if (!build_inner_faces)
591  continue;
592 
593  // treat boundaries of cells of different refinement level
594  // inside the domain in case of multigrid separately
595  else if ((dcell->at_boundary(f) == false ||
596  dcell->has_periodic_neighbor(f)) &&
597  mg_level != numbers::invalid_unsigned_int &&
598  dcell->neighbor_or_periodic_neighbor(f)->level() <
599  dcell->level())
600  {
601  face_is_owned[dcell->face(f)->index()] =
602  FaceCategory::multigrid_refinement_edge;
603  }
604  else
605  {
606  typename ::Triangulation<dim>::cell_iterator neighbor =
607  dcell->neighbor_or_periodic_neighbor(f);
608 
609  // neighbor is refined -> face will be treated by neighbor
610  if (use_active_cells && neighbor->has_children() &&
611  hold_all_faces_to_owned_cells == false)
612  continue;
613 
614  bool add_to_ghost = false;
615  const types::subdomain_id
616  id1 = use_active_cells ? dcell->subdomain_id() :
617  dcell->level_subdomain_id(),
618  id2 = use_active_cells ?
619  (neighbor->has_children() ?
620  dcell->neighbor_child_on_subface(f, 0)
621  ->subdomain_id() :
622  neighbor->subdomain_id()) :
623  neighbor->level_subdomain_id();
624 
625  // Check whether the current face should be processed
626  // locally (instead of being processed from the other
627  // side). We process a face locally when we are more refined
628  // (in the active cell case) or when the face is listed in
629  // the `shared_faces` data structure that we built above.
630  if ((id1 == id2 &&
631  (use_active_cells == false || neighbor->is_active())) ||
632  dcell->level() > neighbor->level() ||
633  std::binary_search(
634  inner_faces_at_proc_boundary[id2].shared_faces.begin(),
635  inner_faces_at_proc_boundary[id2].shared_faces.end(),
636  std::make_pair(id1 < id2 ? dcell->id() : neighbor->id(),
637  id1 < id2 ? neighbor->id() :
638  dcell->id())))
639  {
640  face_is_owned[dcell->face(f)->index()] =
641  FaceCategory::locally_active_done_here;
642  if (dcell->level() == neighbor->level())
643  face_is_owned
644  [neighbor
645  ->face(dcell->has_periodic_neighbor(f) ?
646  dcell->periodic_neighbor_face_no(f) :
647  dcell->neighbor_face_no(f))
648  ->index()] =
649  FaceCategory::locally_active_done_here;
650 
651  // If neighbor is a ghost element (i.e.
652  // dcell->subdomain_id !
653  // dcell->neighbor(f)->subdomain_id()), we need to add its
654  // index into cell level list.
655  if (use_active_cells)
656  add_to_ghost =
657  (dcell->subdomain_id() != neighbor->subdomain_id());
658  else
659  add_to_ghost = (dcell->level_subdomain_id() !=
660  neighbor->level_subdomain_id());
661  }
662  else if (hold_all_faces_to_owned_cells == true)
663  {
664  // add all cells to ghost layer...
665  face_is_owned[dcell->face(f)->index()] =
666  FaceCategory::ghosted;
667  if (use_active_cells)
668  {
669  if (neighbor->has_children())
670  for (unsigned int s = 0;
671  s < dcell->face(f)->n_children();
672  ++s)
673  if (dcell->at_boundary(f))
674  {
675  if (dcell
676  ->periodic_neighbor_child_on_subface(f,
677  s)
678  ->subdomain_id() !=
679  dcell->subdomain_id())
680  add_to_ghost = true;
681  }
682  else
683  {
684  if (dcell->neighbor_child_on_subface(f, s)
685  ->subdomain_id() !=
686  dcell->subdomain_id())
687  add_to_ghost = true;
688  }
689  else
690  add_to_ghost = (dcell->subdomain_id() !=
691  neighbor->subdomain_id());
692  }
693  else
694  add_to_ghost = (dcell->level_subdomain_id() !=
695  neighbor->level_subdomain_id());
696  }
697 
698  if (add_to_ghost)
699  {
700  if (use_active_cells && neighbor->has_children())
701  for (unsigned int s = 0;
702  s < dcell->face(f)->n_children();
703  ++s)
704  {
705  typename ::Triangulation<dim>::cell_iterator
706  neighbor_child =
707  dcell->at_boundary(f) ?
708  dcell->periodic_neighbor_child_on_subface(f,
709  s) :
710  dcell->neighbor_child_on_subface(f, s);
711  if (neighbor_child->subdomain_id() !=
712  dcell->subdomain_id())
713  ghost_cells.insert(
714  std::pair<unsigned int, unsigned int>(
715  neighbor_child->level(),
716  neighbor_child->index()));
717  }
718  else
719  ghost_cells.insert(
720  std::pair<unsigned int, unsigned int>(
721  neighbor->level(), neighbor->index()));
722  at_processor_boundary[i] = true;
723  }
724  }
725  }
726  }
727 
728  // step 2: append the ghost cells at the end of the locally owned
729  // cells
730  for (const auto &ghost_cell : ghost_cells)
731  cell_levels.push_back(ghost_cell);
732  }
733 
734 
735 
736  template <int dim>
737  void
738  FaceSetup<dim>::generate_faces(
739  const ::Triangulation<dim> & triangulation,
740  const std::vector<std::pair<unsigned int, unsigned int>> &cell_levels,
741  TaskInfo & task_info)
742  {
743  const bool is_mixed_mesh = triangulation.is_mixed_mesh();
744 
745  // step 1: create the inverse map between cell iterators and the
746  // cell_level_index field
747  std::map<std::pair<unsigned int, unsigned int>, unsigned int>
748  map_to_vectorized;
749  for (unsigned int cell = 0; cell < cell_levels.size(); ++cell)
750  if (cell == 0 || cell_levels[cell] != cell_levels[cell - 1])
751  {
752  typename ::Triangulation<dim>::cell_iterator dcell(
753  &triangulation,
754  cell_levels[cell].first,
755  cell_levels[cell].second);
756  std::pair<unsigned int, unsigned int> level_index(dcell->level(),
757  dcell->index());
758  map_to_vectorized[level_index] = cell;
759  }
760 
761  // step 2: fill the information about inner faces and boundary faces
762  const unsigned int vectorization_length = task_info.vectorization_length;
763  task_info.face_partition_data.resize(
764  task_info.cell_partition_data.size() - 1, 0);
765  task_info.boundary_partition_data.resize(
766  task_info.cell_partition_data.size() - 1, 0);
767  std::vector<unsigned char> face_visited(face_is_owned.size(), 0);
768  for (unsigned int partition = 0;
769  partition < task_info.cell_partition_data.size() - 2;
770  ++partition)
771  {
772  unsigned int boundary_counter = 0;
773  unsigned int inner_counter = 0;
774  for (unsigned int cell = task_info.cell_partition_data[partition] *
775  vectorization_length;
776  cell < task_info.cell_partition_data[partition + 1] *
777  vectorization_length;
778  ++cell)
779  if (cell == 0 || cell_levels[cell] != cell_levels[cell - 1])
780  {
781  typename ::Triangulation<dim>::cell_iterator dcell(
782  &triangulation,
783  cell_levels[cell].first,
784  cell_levels[cell].second);
785  for (const auto f : dcell->face_indices())
786  {
787  // boundary face
788  if (face_is_owned[dcell->face(f)->index()] ==
789  FaceCategory::locally_active_at_boundary)
790  {
791  Assert(dcell->at_boundary(f), ExcInternalError());
792  ++boundary_counter;
793  FaceToCellTopology<1> info;
794  info.cells_interior[0] = cell;
795  info.cells_exterior[0] = numbers::invalid_unsigned_int;
796  info.interior_face_no = f;
797  info.exterior_face_no = dcell->face(f)->boundary_id();
798  info.face_type =
799  is_mixed_mesh ?
800  (dcell->face(f)->reference_cell() !=
801  ::ReferenceCells::get_hypercube<dim - 1>()) :
802  0;
803  info.subface_index =
805  info.face_orientation = 0;
806  boundary_faces.push_back(info);
807 
808  face_visited[dcell->face(f)->index()]++;
809  }
810  // interior face, including faces over periodic boundaries
811  else
812  {
813  typename ::Triangulation<dim>::cell_iterator
814  neighbor = dcell->neighbor_or_periodic_neighbor(f);
815  if (use_active_cells && neighbor->has_children())
816  {
817  for (unsigned int c = 0;
818  c < dcell->face(f)->n_children();
819  ++c)
820  {
821  typename ::Triangulation<
822  dim>::cell_iterator neighbor_c =
823  dcell->at_boundary(f) ?
824  dcell->periodic_neighbor_child_on_subface(
825  f, c) :
826  dcell->neighbor_child_on_subface(f, c);
827  const types::subdomain_id neigh_domain =
828  neighbor_c->subdomain_id();
829  const unsigned int neighbor_face_no =
830  dcell->has_periodic_neighbor(f) ?
831  dcell->periodic_neighbor_face_no(f) :
832  dcell->neighbor_face_no(f);
833  if (neigh_domain != dcell->subdomain_id() ||
834  face_visited
835  [dcell->face(f)->child(c)->index()] ==
836  1)
837  {
838  std::pair<unsigned int, unsigned int>
839  level_index(neighbor_c->level(),
840  neighbor_c->index());
841  if (face_is_owned
842  [dcell->face(f)->child(c)->index()] ==
843  FaceCategory::locally_active_done_here)
844  {
845  ++inner_counter;
846  inner_faces.push_back(create_face(
847  neighbor_face_no,
848  neighbor_c,
849  map_to_vectorized[level_index],
850  dcell,
851  cell,
852  is_mixed_mesh));
853  }
854  else if (face_is_owned[dcell->face(f)
855  ->child(c)
856  ->index()] ==
857  FaceCategory::ghosted)
858  {
859  inner_ghost_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
868  Assert(
869  face_is_owned[dcell->face(f)
870  ->index()] ==
871  FaceCategory::
872  locally_active_done_elsewhere ||
873  face_is_owned[dcell->face(f)
874  ->index()] ==
875  FaceCategory::ghosted,
876  ExcInternalError());
877  }
878  else
879  {
880  face_visited
881  [dcell->face(f)->child(c)->index()] = 1;
882  }
883  }
884  }
885  else
886  {
887  const types::subdomain_id my_domain =
888  use_active_cells ? dcell->subdomain_id() :
889  dcell->level_subdomain_id();
890  const types::subdomain_id neigh_domain =
891  use_active_cells ? neighbor->subdomain_id() :
892  neighbor->level_subdomain_id();
893  if (neigh_domain != my_domain ||
894  face_visited[dcell->face(f)->index()] == 1)
895  {
896  std::pair<unsigned int, unsigned int>
897  level_index(neighbor->level(),
898  neighbor->index());
899  if (face_is_owned[dcell->face(f)->index()] ==
900  FaceCategory::locally_active_done_here)
901  {
902  Assert(use_active_cells ||
903  dcell->level() ==
904  neighbor->level(),
905  ExcInternalError());
906  ++inner_counter;
907  inner_faces.push_back(create_face(
908  f,
909  dcell,
910  cell,
911  neighbor,
912  map_to_vectorized[level_index],
913  is_mixed_mesh));
914  }
915  else if (face_is_owned[dcell->face(f)
916  ->index()] ==
917  FaceCategory::ghosted)
918  {
919  inner_ghost_faces.push_back(create_face(
920  f,
921  dcell,
922  cell,
923  neighbor,
924  map_to_vectorized[level_index],
925  is_mixed_mesh));
926  }
927  }
928  else
929  {
930  face_visited[dcell->face(f)->index()] = 1;
931  if (dcell->has_periodic_neighbor(f))
932  face_visited
933  [neighbor
934  ->face(
935  dcell->periodic_neighbor_face_no(f))
936  ->index()] = 1;
937  }
938  if (face_is_owned[dcell->face(f)->index()] ==
939  FaceCategory::multigrid_refinement_edge)
940  {
941  refinement_edge_faces.push_back(
942  create_face(f,
943  dcell,
944  cell,
945  neighbor,
946  refinement_edge_faces.size(),
947  is_mixed_mesh));
948  }
949  }
950  }
951  }
952  }
953  task_info.face_partition_data[partition + 1] =
954  task_info.face_partition_data[partition] + inner_counter;
955  task_info.boundary_partition_data[partition + 1] =
956  task_info.boundary_partition_data[partition] + boundary_counter;
957  }
958  task_info.ghost_face_partition_data.resize(2);
959  task_info.ghost_face_partition_data[0] = 0;
960  task_info.ghost_face_partition_data[1] = inner_ghost_faces.size();
961  task_info.refinement_edge_face_partition_data.resize(2);
962  task_info.refinement_edge_face_partition_data[0] = 0;
963  task_info.refinement_edge_face_partition_data[1] =
964  refinement_edge_faces.size();
965  }
966 
967 
968 
969  template <int dim>
970  FaceToCellTopology<1>
971  FaceSetup<dim>::create_face(
972  const unsigned int face_no,
973  const typename ::Triangulation<dim>::cell_iterator &cell,
974  const unsigned int number_cell_interior,
975  const typename ::Triangulation<dim>::cell_iterator &neighbor,
976  const unsigned int number_cell_exterior,
977  const bool is_mixed_mesh)
978  {
979  FaceToCellTopology<1> info;
980  info.cells_interior[0] = number_cell_interior;
981  info.cells_exterior[0] = number_cell_exterior;
982  info.interior_face_no = face_no;
983  if (cell->has_periodic_neighbor(face_no))
984  info.exterior_face_no = cell->periodic_neighbor_face_no(face_no);
985  else
986  info.exterior_face_no = cell->neighbor_face_no(face_no);
987 
988  info.face_type = is_mixed_mesh ?
989  (cell->face(face_no)->reference_cell() !=
990  ::ReferenceCells::get_hypercube<dim - 1>()) :
991  0;
992 
993  info.subface_index = GeometryInfo<dim>::max_children_per_cell;
994  Assert(neighbor->level() <= cell->level(), ExcInternalError());
995  if (cell->level() > neighbor->level())
996  {
997  if (cell->has_periodic_neighbor(face_no))
998  info.subface_index =
999  cell->periodic_neighbor_of_coarser_periodic_neighbor(face_no)
1000  .second;
1001  else
1002  info.subface_index =
1003  cell->neighbor_of_coarser_neighbor(face_no).second;
1004  }
1005 
1006  // special treatment of periodic boundaries
1007  if (dim == 3 && cell->has_periodic_neighbor(face_no))
1008  {
1009  const unsigned int exterior_face_orientation =
1010  !cell->get_triangulation()
1011  .get_periodic_face_map()
1012  .at({cell, face_no})
1013  .second[0] +
1014  2 * cell->get_triangulation()
1015  .get_periodic_face_map()
1016  .at({cell, face_no})
1017  .second[1] +
1018  4 * cell->get_triangulation()
1019  .get_periodic_face_map()
1020  .at({cell, face_no})
1021  .second[2];
1022 
1023  info.face_orientation = exterior_face_orientation;
1024 
1025  return info;
1026  }
1027 
1028  info.face_orientation = 0;
1029  const unsigned int interior_face_orientation =
1030  !cell->face_orientation(face_no) + 2 * cell->face_flip(face_no) +
1031  4 * cell->face_rotation(face_no);
1032  const unsigned int exterior_face_orientation =
1033  !neighbor->face_orientation(info.exterior_face_no) +
1034  2 * neighbor->face_flip(info.exterior_face_no) +
1035  4 * neighbor->face_rotation(info.exterior_face_no);
1036  if (interior_face_orientation != 0)
1037  {
1038  info.face_orientation = 8 + interior_face_orientation;
1039  Assert(exterior_face_orientation == 0,
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() && exterior_face_orientation > 0)
1049  {
1050  const Table<2, unsigned int> orientation =
1051  ShapeInfo<double>::compute_orientation_table(2);
1052  const std::array<unsigned int, 8> inverted_orientations{
1053  {0, 1, 2, 3, 6, 5, 4, 7}};
1054  info.subface_index =
1055  orientation[inverted_orientations[exterior_face_orientation]]
1056  [info.subface_index];
1057  }
1058 
1059  return info;
1060  }
1061 
1062 
1063 
1070  inline bool
1071  compare_faces_for_vectorization(
1072  const FaceToCellTopology<1> & face1,
1073  const FaceToCellTopology<1> & face2,
1074  const std::vector<unsigned int> &active_fe_indices,
1075  const unsigned int length)
1076  {
1077  if (face1.interior_face_no != face2.interior_face_no)
1078  return false;
1079  if (face1.exterior_face_no != face2.exterior_face_no)
1080  return false;
1081  if (face1.subface_index != face2.subface_index)
1082  return false;
1083  if (face1.face_orientation != face2.face_orientation)
1084  return false;
1085  if (face1.face_type != face2.face_type)
1086  return false;
1087 
1088  if (active_fe_indices.size() > 0)
1089  {
1090  if (active_fe_indices[face1.cells_interior[0] / length] !=
1091  active_fe_indices[face2.cells_interior[0] / length])
1092  return false;
1093 
1094  if (face2.cells_exterior[0] != numbers::invalid_unsigned_int)
1095  if (active_fe_indices[face1.cells_exterior[0] / length] !=
1096  active_fe_indices[face2.cells_exterior[0] / length])
1097  return false;
1098  }
1099 
1100  return true;
1101  }
1102 
1103 
1104 
1111  template <int length>
1112  struct FaceComparator
1113  {
1114  FaceComparator(const std::vector<unsigned int> &active_fe_indices)
1115  : active_fe_indices(active_fe_indices)
1116  {}
1117 
1118  bool
1119  operator()(const FaceToCellTopology<length> &face1,
1120  const FaceToCellTopology<length> &face2) const
1121  {
1122  // check if active FE indices match
1123  if (face1.face_type < face2.face_type)
1124  return true;
1125  else if (face1.face_type > face2.face_type)
1126  return false;
1127 
1128  // check if active FE indices match
1129  if (active_fe_indices.size() > 0)
1130  {
1131  // ... for interior faces
1132  if (active_fe_indices[face1.cells_interior[0] / length] <
1133  active_fe_indices[face2.cells_interior[0] / length])
1134  return true;
1135  else if (active_fe_indices[face1.cells_interior[0] / length] >
1136  active_fe_indices[face2.cells_interior[0] / length])
1137  return false;
1138 
1139  // ... for exterior faces
1140  if (face2.cells_exterior[0] != numbers::invalid_unsigned_int)
1141  {
1142  if (active_fe_indices[face1.cells_exterior[0] / length] <
1143  active_fe_indices[face2.cells_exterior[0] / length])
1144  return true;
1145  else if (active_fe_indices[face1.cells_exterior[0] / length] >
1146  active_fe_indices[face2.cells_exterior[0] / length])
1147  return false;
1148  }
1149  }
1150 
1151  for (unsigned int i = 0; i < length; ++i)
1152  if (face1.cells_interior[i] < face2.cells_interior[i])
1153  return true;
1154  else if (face1.cells_interior[i] > face2.cells_interior[i])
1155  return false;
1156  for (unsigned int i = 0; i < length; ++i)
1157  if (face1.cells_exterior[i] < face2.cells_exterior[i])
1158  return true;
1159  else if (face1.cells_exterior[i] > face2.cells_exterior[i])
1160  return false;
1161  if (face1.interior_face_no < face2.interior_face_no)
1162  return true;
1163  else if (face1.interior_face_no > face2.interior_face_no)
1164  return false;
1165  if (face1.exterior_face_no < face2.exterior_face_no)
1166  return true;
1167  else if (face1.exterior_face_no > face2.exterior_face_no)
1168  return false;
1169 
1170  // we do not need to check for subface_index and orientation because
1171  // those cannot be different if when all the other values are the
1172  // same.
1173  AssertDimension(face1.subface_index, face2.subface_index);
1174  AssertDimension(face1.face_orientation, face2.face_orientation);
1175 
1176  return false;
1177  }
1178 
1179  private:
1180  const std::vector<unsigned int> &active_fe_indices;
1181  };
1182 
1183 
1184 
1185  template <int vectorization_width>
1186  void
1188  const std::vector<FaceToCellTopology<1>> &faces_in,
1189  const std::vector<bool> & hard_vectorization_boundary,
1190  std::vector<unsigned int> & face_partition_data,
1191  std::vector<FaceToCellTopology<vectorization_width>> &faces_out,
1192  const std::vector<unsigned int> & active_fe_indices)
1193  {
1194  FaceToCellTopology<vectorization_width> face_batch;
1195  std::vector<std::vector<unsigned int>> faces_type;
1196 
1197  unsigned int face_start = face_partition_data[0],
1198  face_end = face_partition_data[0];
1199 
1200  face_partition_data[0] = faces_out.size();
1201  for (unsigned int partition = 0;
1202  partition < face_partition_data.size() - 1;
1203  ++partition)
1204  {
1205  std::vector<std::vector<unsigned int>> new_faces_type;
1206 
1207  // start with the end point for the last partition
1208  face_start = face_end;
1209  face_end = face_partition_data[partition + 1];
1210 
1211  // set the partitioner to the new vectorized lengths
1212  face_partition_data[partition + 1] = face_partition_data[partition];
1213 
1214  // loop over the faces in the current partition and reorder according
1215  // to the face type
1216  for (unsigned int face = face_start; face < face_end; ++face)
1217  {
1218  for (auto &face_type : faces_type)
1219  {
1220  // Compare current face with first face of type type
1221  if (compare_faces_for_vectorization(faces_in[face],
1222  faces_in[face_type[0]],
1223  active_fe_indices,
1224  vectorization_width))
1225  {
1226  face_type.push_back(face);
1227  goto face_found;
1228  }
1229  }
1230  faces_type.emplace_back(1, face);
1231  face_found:
1232  {}
1233  }
1234 
1235  // insert new faces in sorted list to get good data locality
1236  FaceComparator<vectorization_width> face_comparator(
1237  active_fe_indices);
1238  std::set<FaceToCellTopology<vectorization_width>,
1239  FaceComparator<vectorization_width>>
1240  new_faces(face_comparator);
1241  for (const auto &face_type : faces_type)
1242  {
1243  face_batch.face_type = faces_in[face_type[0]].face_type;
1244  face_batch.interior_face_no =
1245  faces_in[face_type[0]].interior_face_no;
1246  face_batch.exterior_face_no =
1247  faces_in[face_type[0]].exterior_face_no;
1248  face_batch.subface_index = faces_in[face_type[0]].subface_index;
1249  face_batch.face_orientation =
1250  faces_in[face_type[0]].face_orientation;
1251  unsigned int no_faces = face_type.size();
1252  std::vector<unsigned char> touched(no_faces, 0);
1253 
1254  // do two passes through the data. The first is to identify
1255  // similar faces within the same index range as the cells which
1256  // will allow for vectorized read operations, the second picks up
1257  // all the rest
1258  unsigned int n_vectorized = 0;
1259  for (unsigned int f = 0; f < no_faces; ++f)
1260  if (faces_in[face_type[f]].cells_interior[0] %
1261  vectorization_width ==
1262  0)
1263  {
1264  bool is_contiguous = true;
1265  if (f + vectorization_width > no_faces)
1266  is_contiguous = false;
1267  else
1268  for (unsigned int v = 1; v < vectorization_width; ++v)
1269  if (faces_in[face_type[f + v]].cells_interior[0] !=
1270  faces_in[face_type[f]].cells_interior[0] + v)
1271  is_contiguous = false;
1272  if (is_contiguous)
1273  {
1274  AssertIndexRange(f,
1275  face_type.size() -
1276  vectorization_width + 1);
1277  for (unsigned int v = 0; v < vectorization_width; ++v)
1278  {
1279  face_batch.cells_interior[v] =
1280  faces_in[face_type[f + v]].cells_interior[0];
1281  face_batch.cells_exterior[v] =
1282  faces_in[face_type[f + v]].cells_exterior[0];
1283  touched[f + v] = 1;
1284  }
1285  new_faces.insert(face_batch);
1286  f += vectorization_width - 1;
1287  n_vectorized += vectorization_width;
1288  }
1289  }
1290 
1291  std::vector<unsigned int> untouched;
1292  untouched.reserve(no_faces - n_vectorized);
1293  for (unsigned int f = 0; f < no_faces; ++f)
1294  if (touched[f] == 0)
1295  untouched.push_back(f);
1296  unsigned int v = 0;
1297  for (const auto f : untouched)
1298  {
1299  face_batch.cells_interior[v] =
1300  faces_in[face_type[f]].cells_interior[0];
1301  face_batch.cells_exterior[v] =
1302  faces_in[face_type[f]].cells_exterior[0];
1303  ++v;
1304  if (v == vectorization_width)
1305  {
1306  new_faces.insert(face_batch);
1307  v = 0;
1308  }
1309  }
1310  if (v > 0 && v < vectorization_width)
1311  {
1312  // must add non-filled face
1313  if (hard_vectorization_boundary[partition + 1] ||
1314  partition == face_partition_data.size() - 2)
1315  {
1316  for (; v < vectorization_width; ++v)
1317  {
1318  // Dummy cell, not used
1319  face_batch.cells_interior[v] =
1321  face_batch.cells_exterior[v] =
1323  }
1324  new_faces.insert(face_batch);
1325  }
1326  else
1327  {
1328  // postpone to the next partition
1329  std::vector<unsigned int> untreated(v);
1330  for (unsigned int f = 0; f < v; ++f)
1331  untreated[f] = face_type[*(untouched.end() - 1 - f)];
1332  new_faces_type.push_back(untreated);
1333  }
1334  }
1335  }
1336 
1337  // insert sorted list to vector of faces
1338  for (auto it = new_faces.begin(); it != new_faces.end(); ++it)
1339  faces_out.push_back(*it);
1340  face_partition_data[partition + 1] += new_faces.size();
1341 
1342  // set the faces that were left over to faces_type for the next round
1343  faces_type = std::move(new_faces_type);
1344  }
1345 
1346 # ifdef DEBUG
1347  // final safety checks
1348  for (const auto &face_type : faces_type)
1349  AssertDimension(face_type.size(), 0U);
1350 
1351  AssertDimension(faces_out.size(), face_partition_data.back());
1352  unsigned int nfaces = 0;
1353  for (unsigned int i = face_partition_data[0];
1354  i < face_partition_data.back();
1355  ++i)
1356  for (unsigned int v = 0; v < vectorization_width; ++v)
1357  nfaces +=
1358  (faces_out[i].cells_interior[v] != numbers::invalid_unsigned_int);
1359  AssertDimension(nfaces, faces_in.size());
1360 
1361  std::vector<std::pair<unsigned int, unsigned int>> in_faces, out_faces;
1362  for (const auto &face_in : faces_in)
1363  in_faces.emplace_back(face_in.cells_interior[0],
1364  face_in.cells_exterior[0]);
1365  for (unsigned int i = face_partition_data[0];
1366  i < face_partition_data.back();
1367  ++i)
1368  for (unsigned int v = 0;
1369  v < vectorization_width &&
1370  faces_out[i].cells_interior[v] != numbers::invalid_unsigned_int;
1371  ++v)
1372  out_faces.emplace_back(faces_out[i].cells_interior[v],
1373  faces_out[i].cells_exterior[v]);
1374  std::sort(in_faces.begin(), in_faces.end());
1375  std::sort(out_faces.begin(), out_faces.end());
1376  AssertDimension(in_faces.size(), out_faces.size());
1377  for (unsigned int i = 0; i < in_faces.size(); ++i)
1378  {
1379  AssertDimension(in_faces[i].first, out_faces[i].first);
1380  AssertDimension(in_faces[i].second, out_faces[i].second);
1381  }
1382 # endif
1383  }
1384 
1385 #endif // ifndef DOXYGEN
1386 
1387  } // namespace MatrixFreeFunctions
1388 } // namespace internal
1389 
1390 
1392 
1393 #endif
Definition: cell_id.h:71
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
Point< 2 > second
Definition: grid_out.cc:4605
Point< 2 > first
Definition: grid_out.cc:4604
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
Definition: exceptions.h:1473
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1667
#define AssertThrowMPI(error_code)
Definition: exceptions.h:1790
#define AssertIndexRange(index, range)
Definition: exceptions.h:1732
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1583
static const char U
constexpr const ReferenceCell & get_hypercube()
void partition(const SparsityPattern &sparsity_pattern, const unsigned int n_partitions, std::vector< unsigned int > &partition_indices, const Partitioner partitioner=Partitioner::metis)
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)
const types::subdomain_id invalid_subdomain_id
Definition: types.h:281
static const unsigned int invalid_unsigned_int
Definition: types.h:201
unsigned int subdomain_id
Definition: types.h:43
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
std::vector< std::pair< CellId, CellId > > shared_faces
std::vector< FaceToCellTopology< 1 > > inner_faces
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)
std::vector< FaceToCellTopology< 1 > > boundary_faces
std::vector< FaceToCellTopology< 1 > > refinement_edge_faces
void generate_faces(const ::Triangulation< dim > &triangulation, const std::vector< std::pair< unsigned int, unsigned int >> &cell_levels, TaskInfo &task_info)
std::vector< FaceToCellTopology< 1 > > inner_ghost_faces
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)
const MPI_Comm & comm