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