Reference documentation for deal.II version Git 2276b94619 2020-02-28 10:11:55 +0100
\(\newcommand{\vcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\vcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
mesh_loop.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2017 - 2019 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_mesh_worker_mesh_loop_h
18 #define dealii_mesh_worker_mesh_loop_h
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/template_constraints.h>
23 #include <deal.II/base/work_stream.h>
24 
25 #include <deal.II/grid/filtered_iterator.h>
26 #include <deal.II/grid/tria.h>
27 
28 #include <deal.II/meshworker/assemble_flags.h>
29 #include <deal.II/meshworker/dof_info.h>
30 #include <deal.II/meshworker/integration_info.h>
31 #include <deal.II/meshworker/local_integrator.h>
32 #include <deal.II/meshworker/loop.h>
33 
34 #include <functional>
35 #include <type_traits>
36 
37 DEAL_II_NAMESPACE_OPEN
38 
39 // Forward declaration
40 #ifndef DOXYGEN
41 template <typename>
42 class TriaActiveIterator;
43 #endif
44 
45 namespace MeshWorker
46 {
47  namespace internal
48  {
53  template <class CellIteratorType>
55  {
59  using type = CellIteratorType;
60  };
61 
69  template <class CellIteratorType>
70  struct CellIteratorBaseType<IteratorOverIterators<CellIteratorType>>
71  {
75  // Since we can have filtered iterators and the like as template
76  // arguments, we recursivelyremove the template layers to retrieve the
77  // underlying iterator type.
79  };
80 
89  template <class CellIteratorType>
90  struct CellIteratorBaseType<FilteredIterator<CellIteratorType>>
91  {
95  // Since we can have nested filtered iterators, we recursively
96  // remove the template layers to retrieve the underlying iterator type.
98  };
99  } // namespace internal
100 
183  template <class CellIteratorType,
184  class ScratchData,
185  class CopyData,
186  class CellIteratorBaseType =
188  void
190  const CellIteratorType & begin,
191  const typename identity<CellIteratorType>::type &end,
192 
193  const typename identity<std::function<
194  void(const CellIteratorBaseType &, ScratchData &, CopyData &)>>::type
195  &cell_worker,
196  const typename identity<std::function<void(const CopyData &)>>::type
197  &copier,
198 
199  const ScratchData &sample_scratch_data,
200  const CopyData & sample_copy_data,
201 
202  const AssembleFlags flags = assemble_own_cells,
203 
204  const typename identity<std::function<void(const CellIteratorBaseType &,
205  const unsigned int,
206  ScratchData &,
207  CopyData &)>>::type
208  &boundary_worker = std::function<void(const CellIteratorBaseType &,
209  const unsigned int,
210  ScratchData &,
211  CopyData &)>(),
212 
213  const typename identity<std::function<void(const CellIteratorBaseType &,
214  const unsigned int,
215  const unsigned int,
216  const CellIteratorBaseType &,
217  const unsigned int,
218  const unsigned int,
219  ScratchData &,
220  CopyData &)>>::type
221  &face_worker = std::function<void(const CellIteratorBaseType &,
222  const unsigned int,
223  const unsigned int,
224  const CellIteratorBaseType &,
225  const unsigned int,
226  const unsigned int,
227  ScratchData &,
228  CopyData &)>(),
229 
230  const unsigned int queue_length = 2 * MultithreadInfo::n_threads(),
231  const unsigned int chunk_size = 8)
232  {
233  Assert(
234  (!cell_worker) == !(flags & work_on_cells),
235  ExcMessage(
236  "If you specify a cell_worker, you need to set assemble_own_cells or assemble_ghost_cells."));
237 
238  Assert(
239  (flags &
242  ExcMessage(
243  "You can only specify assemble_own_interior_faces_once OR assemble_own_interior_faces_both."));
244 
245  Assert(
248  ExcMessage(
249  "You can only specify assemble_ghost_faces_once OR assemble_ghost_faces_both."));
250 
251  Assert(
252  !(flags & cells_after_faces) ||
254  ExcMessage(
255  "The option cells_after_faces only makes sense if you assemble on cells."));
256 
257  Assert((!face_worker) == !(flags & work_on_faces),
258  ExcMessage(
259  "If you specify a face_worker, assemble_face_* needs to be set."));
260 
261  Assert(
262  (!boundary_worker) == !(flags & assemble_boundary_faces),
263  ExcMessage(
264  "If you specify a boundary_worker, assemble_boundary_faces needs to be set."));
265 
266  auto cell_action = [&](const CellIteratorBaseType &cell,
267  ScratchData & scratch,
268  CopyData & copy) {
269  // First reset the CopyData class to the empty copy_data given by the
270  // user.
271  copy = sample_copy_data;
272 
273  const bool ignore_subdomain =
274  (cell->get_triangulation().locally_owned_subdomain() ==
276 
277  types::subdomain_id current_subdomain_id =
278  (cell->is_level_cell() ? cell->level_subdomain_id() :
279  cell->subdomain_id());
280 
281  const bool own_cell =
282  ignore_subdomain ||
283  (current_subdomain_id ==
284  cell->get_triangulation().locally_owned_subdomain());
285 
286  if ((!ignore_subdomain) &&
287  (current_subdomain_id == numbers::artificial_subdomain_id))
288  return;
289 
290  if (!(flags & (cells_after_faces)) &&
291  (((flags & (assemble_own_cells)) && own_cell) ||
292  ((flags & assemble_ghost_cells) && !own_cell)))
293  cell_worker(cell, scratch, copy);
294 
295  if (flags & (work_on_faces | work_on_boundary))
296  for (unsigned int face_no = 0;
297  face_no < GeometryInfo<CellIteratorBaseType::AccessorType::
298  Container::dimension>::faces_per_cell;
299  ++face_no)
300  {
301  if (cell->at_boundary(face_no) &&
302  !cell->has_periodic_neighbor(face_no))
303  {
304  // only integrate boundary faces of own cells
305  if ((flags & assemble_boundary_faces) && own_cell)
306  boundary_worker(cell, face_no, scratch, copy);
307  }
308  else
309  {
310  // interior face, potentially assemble
312  neighbor = cell->neighbor_or_periodic_neighbor(face_no);
313 
314  types::subdomain_id neighbor_subdomain_id =
316  if (neighbor->is_level_cell())
317  neighbor_subdomain_id = neighbor->level_subdomain_id();
318  // subdomain id is only valid for active cells
319  else if (neighbor->is_active())
320  neighbor_subdomain_id = neighbor->subdomain_id();
321 
322  const bool own_neighbor =
323  ignore_subdomain ||
324  (neighbor_subdomain_id ==
325  cell->get_triangulation().locally_owned_subdomain());
326 
327  // skip all faces between two ghost cells
328  if (!own_cell && !own_neighbor)
329  continue;
330 
331  // skip if the user doesn't want faces between own cells
332  if (own_cell && own_neighbor &&
335  continue;
336 
337  // skip face to ghost
338  if (own_cell != own_neighbor &&
339  !(flags &
341  continue;
342 
343  // Deal with refinement edges from the refined side. Assuming
344  // one-irregular meshes, this situation should only occur if
345  // both cells are active.
346  const bool periodic_neighbor =
347  cell->has_periodic_neighbor(face_no);
348 
349  if ((!periodic_neighbor &&
350  cell->neighbor_is_coarser(face_no)) ||
351  (periodic_neighbor &&
352  cell->periodic_neighbor_is_coarser(face_no)))
353  {
354  Assert(cell->is_active(), ExcInternalError());
355  Assert(neighbor->is_active(), ExcInternalError());
356 
357  // skip if only one processor needs to assemble the face
358  // to a ghost cell and the fine cell is not ours.
359  if (!own_cell && (flags & assemble_ghost_faces_once))
360  continue;
361 
362  const std::pair<unsigned int, unsigned int>
363  neighbor_face_no =
364  periodic_neighbor ?
365  cell->periodic_neighbor_of_coarser_periodic_neighbor(
366  face_no) :
367  cell->neighbor_of_coarser_neighbor(face_no);
368 
369  face_worker(cell,
370  face_no,
372  neighbor,
373  neighbor_face_no.first,
374  neighbor_face_no.second,
375  scratch,
376  copy);
377 
379  {
380  // If own faces are to be assembled from both sides,
381  // call the faceworker again with swapped arguments.
382  // This is because we won't be looking at an adaptively
383  // refined edge coming from the other side.
384  face_worker(neighbor,
385  neighbor_face_no.first,
386  neighbor_face_no.second,
387  cell,
388  face_no,
390  scratch,
391  copy);
392  }
393  }
394  else
395  {
396  // If iterator is active and neighbor is refined, skip
397  // internal face.
398  if (::internal::is_active_iterator(cell) &&
399  neighbor->has_children())
400  continue;
401 
402  // Now neighbor is on same level, double-check this:
403  Assert(cell->level() == neighbor->level(),
404  ExcInternalError());
405 
406  // If we own both cells only do faces from one side (unless
407  // AssembleFlags says otherwise). Here, we rely on cell
408  // comparison that will look at cell->index().
409  if (own_cell && own_neighbor &&
411  (neighbor < cell))
412  continue;
413 
414  // We only look at faces to ghost on the same level once
415  // (only where own_cell=true and own_neighbor=false)
416  if (!own_cell)
417  continue;
418 
419  // now only one processor assembles faces_to_ghost. We let
420  // the processor with the smaller (level-)subdomain id
421  // assemble the face.
422  if (own_cell && !own_neighbor &&
423  (flags & assemble_ghost_faces_once) &&
424  (neighbor_subdomain_id < current_subdomain_id))
425  continue;
426 
427  const unsigned int neighbor_face_no =
428  periodic_neighbor ?
429  cell->periodic_neighbor_face_no(face_no) :
430  cell->neighbor_face_no(face_no);
431  Assert(periodic_neighbor ||
432  neighbor->face(neighbor_face_no) ==
433  cell->face(face_no),
434  ExcInternalError());
435 
436  face_worker(cell,
437  face_no,
439  neighbor,
440  neighbor_face_no,
442  scratch,
443  copy);
444  }
445  }
446  } // faces
447 
448  // Execute the cell_worker if faces are handled before cells
449  if ((flags & cells_after_faces) &&
450  (((flags & assemble_own_cells) && own_cell) ||
451  ((flags & assemble_ghost_cells) && !own_cell)))
452  cell_worker(cell, scratch, copy);
453  };
454 
455  // Submit to workstream
456  WorkStream::run(begin,
457  end,
458  cell_action,
459  copier,
460  sample_scratch_data,
461  sample_copy_data,
462  queue_length,
463  chunk_size);
464  }
465 
535  template <class CellIteratorType,
536  class ScratchData,
537  class CopyData,
538  class CellIteratorBaseType =
540  void
542  IteratorRange<CellIteratorType> iterator_range,
543  const typename identity<std::function<
544  void(const CellIteratorBaseType &, ScratchData &, CopyData &)>>::type
545  &cell_worker,
546  const typename identity<std::function<void(const CopyData &)>>::type
547  &copier,
548 
549  const ScratchData &sample_scratch_data,
550  const CopyData & sample_copy_data,
551 
552  const AssembleFlags flags = assemble_own_cells,
553 
554  const typename identity<std::function<void(const CellIteratorBaseType &,
555  const unsigned int,
556  ScratchData &,
557  CopyData &)>>::type
558  &boundary_worker = std::function<void(const CellIteratorBaseType &,
559  const unsigned int,
560  ScratchData &,
561  CopyData &)>(),
562 
563  const typename identity<std::function<void(const CellIteratorBaseType &,
564  const unsigned int,
565  const unsigned int,
566  const CellIteratorBaseType &,
567  const unsigned int,
568  const unsigned int,
569  ScratchData &,
570  CopyData &)>>::type
571  &face_worker = std::function<void(const CellIteratorBaseType &,
572  const unsigned int,
573  const unsigned int,
574  const CellIteratorBaseType &,
575  const unsigned int,
576  const unsigned int,
577  ScratchData &,
578  CopyData &)>(),
579 
580  const unsigned int queue_length = 2 * MultithreadInfo::n_threads(),
581  const unsigned int chunk_size = 8)
582  {
583  // Call the function above
584  mesh_loop<typename IteratorRange<CellIteratorType>::IteratorOverIterators,
585  ScratchData,
586  CopyData,
587  CellIteratorBaseType>(iterator_range.begin(),
588  iterator_range.end(),
589  cell_worker,
590  copier,
591  sample_scratch_data,
592  sample_copy_data,
593  flags,
594  boundary_worker,
595  face_worker,
596  queue_length,
597  chunk_size);
598  }
599 
660  template <class CellIteratorType,
661  class ScratchData,
662  class CopyData,
663  class MainClass>
664  void
665  mesh_loop(const CellIteratorType & begin,
666  const typename identity<CellIteratorType>::type &end,
667  MainClass & main_class,
668  void (MainClass::*cell_worker)(const CellIteratorType &,
669  ScratchData &,
670  CopyData &),
671  void (MainClass::*copier)(const CopyData &),
672  const ScratchData & sample_scratch_data,
673  const CopyData & sample_copy_data,
674  const AssembleFlags flags = assemble_own_cells,
675  void (MainClass::*boundary_worker)(const CellIteratorType &,
676  const unsigned int,
677  ScratchData &,
678  CopyData &) = nullptr,
679  void (MainClass::*face_worker)(const CellIteratorType &,
680  const unsigned int,
681  const unsigned int,
682  const CellIteratorType &,
683  const unsigned int,
684  const unsigned int,
685  ScratchData &,
686  CopyData &) = nullptr,
687  const unsigned int queue_length = 2 * MultithreadInfo::n_threads(),
688  const unsigned int chunk_size = 8)
689  {
690  std::function<void(const CellIteratorType &, ScratchData &, CopyData &)>
691  f_cell_worker;
692 
693  std::function<void(
694  const CellIteratorType &, const unsigned int, ScratchData &, CopyData &)>
695  f_boundary_worker;
696 
697  std::function<void(const CellIteratorType &,
698  const unsigned int,
699  const unsigned int,
700  const CellIteratorType &,
701  const unsigned int,
702  const unsigned int,
703  ScratchData &,
704  CopyData &)>
705  f_face_worker;
706 
707  if (cell_worker != nullptr)
708  f_cell_worker = [&main_class,
709  cell_worker](const CellIteratorType &cell_iterator,
710  ScratchData & scratch_data,
711  CopyData & copy_data) {
712  (main_class.*cell_worker)(cell_iterator, scratch_data, copy_data);
713  };
714 
715  if (boundary_worker != nullptr)
716  f_boundary_worker =
717  [&main_class, boundary_worker](const CellIteratorType &cell_iterator,
718  const unsigned int face_no,
719  ScratchData & scratch_data,
720  CopyData & copy_data) {
721  (main_class.*
722  boundary_worker)(cell_iterator, face_no, scratch_data, copy_data);
723  };
724 
725  if (face_worker != nullptr)
726  f_face_worker = [&main_class,
727  face_worker](const CellIteratorType &cell_iterator_1,
728  const unsigned int face_index_1,
729  const unsigned int subface_index_1,
730  const CellIteratorType &cell_iterator_2,
731  const unsigned int face_index_2,
732  const unsigned int subface_index_2,
733  ScratchData & scratch_data,
734  CopyData & copy_data) {
735  (main_class.*face_worker)(cell_iterator_1,
736  face_index_1,
737  subface_index_1,
738  cell_iterator_2,
739  face_index_2,
740  subface_index_2,
741  scratch_data,
742  copy_data);
743  };
744 
745  mesh_loop(begin,
746  end,
747  f_cell_worker,
748  [&main_class, copier](const CopyData &copy_data) {
749  (main_class.*copier)(copy_data);
750  },
751  sample_scratch_data,
752  sample_copy_data,
753  flags,
754  f_boundary_worker,
755  f_face_worker,
756  queue_length,
757  chunk_size);
758  }
759 
839  template <class CellIteratorType,
840  class ScratchData,
841  class CopyData,
842  class MainClass,
843  class CellIteratorBaseType =
845  void
847  MainClass & main_class,
848  void (MainClass::*cell_worker)(const CellIteratorBaseType &,
849  ScratchData &,
850  CopyData &),
851  void (MainClass::*copier)(const CopyData &),
852  const ScratchData & sample_scratch_data,
853  const CopyData & sample_copy_data,
854  const AssembleFlags flags = assemble_own_cells,
855  void (MainClass::*boundary_worker)(const CellIteratorBaseType &,
856  const unsigned int,
857  ScratchData &,
858  CopyData &) = nullptr,
859  void (MainClass::*face_worker)(const CellIteratorBaseType &,
860  const unsigned int,
861  const unsigned int,
862  const CellIteratorBaseType &,
863  const unsigned int,
864  const unsigned int,
865  ScratchData &,
866  CopyData &) = nullptr,
867  const unsigned int queue_length = 2 * MultithreadInfo::n_threads(),
868  const unsigned int chunk_size = 8)
869  {
870  // Call the function above
871  mesh_loop<typename IteratorRange<CellIteratorType>::IteratorOverIterators,
872  ScratchData,
873  CopyData,
874  MainClass,
875  CellIteratorBaseType>(iterator_range.begin(),
876  iterator_range.end(),
877  main_class,
878  cell_worker,
879  copier,
880  sample_scratch_data,
881  sample_copy_data,
882  flags,
883  boundary_worker,
884  face_worker,
885  queue_length,
886  chunk_size);
887  }
888 } // namespace MeshWorker
889 
890 DEAL_II_NAMESPACE_CLOSE
891 
892 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:187
IteratorOverIterators begin()
const types::subdomain_id invalid_subdomain_id
Definition: types.h:279
typename CellIteratorBaseType< CellIteratorType >::type type
Definition: mesh_loop.h:97
void mesh_loop(const CellIteratorType &begin, const typename identity< CellIteratorType >::type &end, const typename identity< std::function< void(const CellIteratorBaseType &, ScratchData &, CopyData &)>>::type &cell_worker, const typename identity< std::function< void(const CopyData &)>>::type &copier, const ScratchData &sample_scratch_data, const CopyData &sample_copy_data, const AssembleFlags flags=assemble_own_cells, const typename identity< std::function< void(const CellIteratorBaseType &, const unsigned int, ScratchData &, CopyData &)>>::type &boundary_worker=std::function< void(const CellIteratorBaseType &, const unsigned int, ScratchData &, CopyData &)>(), const typename identity< std::function< void(const CellIteratorBaseType &, const unsigned int, const unsigned int, const CellIteratorBaseType &, const unsigned int, const unsigned int, ScratchData &, CopyData &)>>::type &face_worker=std::function< void(const CellIteratorBaseType &, const unsigned int, const unsigned int, const CellIteratorBaseType &, const unsigned int, const unsigned int, ScratchData &, CopyData &)>(), const unsigned int queue_length=2 *MultithreadInfo::n_threads(), const unsigned int chunk_size=8)
Definition: mesh_loop.h:189
typename CellIteratorBaseType< CellIteratorType >::type type
Definition: mesh_loop.h:78
static ::ExceptionBase & ExcMessage(std::string arg1)
unsigned int subdomain_id
Definition: types.h:43
#define Assert(cond, exc)
Definition: exceptions.h:1411
bool is_active_iterator(const DI &)
Definition: loop.h:49
const types::subdomain_id artificial_subdomain_id
Definition: types.h:296
void run(const std::vector< std::vector< Iterator >> &colored_iterators, Worker worker, Copier copier, const ScratchData &sample_scratch_data, const CopyData &sample_copy_data, const unsigned int queue_length=2 *MultithreadInfo::n_threads(), const unsigned int chunk_size=8)
Definition: work_stream.h:1187
static unsigned int n_threads()
void cell_action(ITERATOR cell, DoFInfoBox< dim, DOFINFO > &dof_info, INFOBOX &info, const std::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &cell_worker, const std::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &boundary_worker, const std::function< void(DOFINFO &, DOFINFO &, typename INFOBOX::CellInfo &, typename INFOBOX::CellInfo &)> &face_worker, const LoopControl &loop_control)
Definition: loop.h:193
static ::ExceptionBase & ExcInternalError()
IteratorOverIterators end() const