Reference documentation for deal.II version Git cb0bd54b52 2019-09-21 16:31:22 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
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 
176  template <class CellIteratorType,
177  class ScratchData,
178  class CopyData,
179  class CellIteratorBaseType =
181  void
183  const CellIteratorType & begin,
184  const typename identity<CellIteratorType>::type &end,
185 
186  const typename identity<std::function<
187  void(const CellIteratorBaseType &, ScratchData &, CopyData &)>>::type
188  &cell_worker,
189  const typename identity<std::function<void(const CopyData &)>>::type
190  &copier,
191 
192  const ScratchData &sample_scratch_data,
193  const CopyData & sample_copy_data,
194 
195  const AssembleFlags flags = assemble_own_cells,
196 
197  const typename identity<std::function<void(const CellIteratorBaseType &,
198  const unsigned int,
199  ScratchData &,
200  CopyData &)>>::type
201  &boundary_worker = std::function<void(const CellIteratorBaseType &,
202  const unsigned int,
203  ScratchData &,
204  CopyData &)>(),
205 
206  const typename identity<std::function<void(const CellIteratorBaseType &,
207  const unsigned int,
208  const unsigned int,
209  const CellIteratorBaseType &,
210  const unsigned int,
211  const unsigned int,
212  ScratchData &,
213  CopyData &)>>::type
214  &face_worker = std::function<void(const CellIteratorBaseType &,
215  const unsigned int,
216  const unsigned int,
217  const CellIteratorBaseType &,
218  const unsigned int,
219  const unsigned int,
220  ScratchData &,
221  CopyData &)>(),
222 
223  const unsigned int queue_length = 2 * MultithreadInfo::n_threads(),
224  const unsigned int chunk_size = 8)
225  {
226  Assert(
227  (!cell_worker) == !(flags & work_on_cells),
228  ExcMessage(
229  "If you specify a cell_worker, you need to set assemble_own_cells or assemble_ghost_cells."));
230 
231  Assert(
232  (flags &
235  ExcMessage(
236  "You can only specify assemble_own_interior_faces_once OR assemble_own_interior_faces_both."));
237 
238  Assert(
241  ExcMessage(
242  "You can only specify assemble_ghost_faces_once OR assemble_ghost_faces_both."));
243 
244  Assert(
245  !(flags & cells_after_faces) ||
247  ExcMessage(
248  "The option cells_after_faces only makes sense if you assemble on cells."));
249 
250  Assert((!face_worker) == !(flags & work_on_faces),
251  ExcMessage(
252  "If you specify a face_worker, assemble_face_* needs to be set."));
253 
254  Assert(
255  (!boundary_worker) == !(flags & assemble_boundary_faces),
256  ExcMessage(
257  "If you specify a boundary_worker, assemble_boundary_faces needs to be set."));
258 
259  auto cell_action = [&](const CellIteratorBaseType &cell,
260  ScratchData & scratch,
261  CopyData & copy) {
262  // First reset the CopyData class to the empty copy_data given by the
263  // user.
264  copy = sample_copy_data;
265 
266  const bool ignore_subdomain =
267  (cell->get_triangulation().locally_owned_subdomain() ==
269 
270  types::subdomain_id current_subdomain_id =
271  (cell->is_level_cell() ? cell->level_subdomain_id() :
272  cell->subdomain_id());
273 
274  const bool own_cell =
275  ignore_subdomain ||
276  (current_subdomain_id ==
277  cell->get_triangulation().locally_owned_subdomain());
278 
279  if ((!ignore_subdomain) &&
280  (current_subdomain_id == numbers::artificial_subdomain_id))
281  return;
282 
283  if (!(flags & (cells_after_faces)) &&
284  (((flags & (assemble_own_cells)) && own_cell) ||
285  ((flags & assemble_ghost_cells) && !own_cell)))
286  cell_worker(cell, scratch, copy);
287 
288  if (flags & (work_on_faces | work_on_boundary))
289  for (unsigned int face_no = 0;
290  face_no < GeometryInfo<CellIteratorBaseType::AccessorType::
291  Container::dimension>::faces_per_cell;
292  ++face_no)
293  {
294  if (cell->at_boundary(face_no) &&
295  !cell->has_periodic_neighbor(face_no))
296  {
297  // only integrate boundary faces of own cells
298  if ((flags & assemble_boundary_faces) && own_cell)
299  boundary_worker(cell, face_no, scratch, copy);
300  }
301  else
302  {
303  // interior face, potentially assemble
305  neighbor = cell->neighbor_or_periodic_neighbor(face_no);
306 
307  types::subdomain_id neighbor_subdomain_id =
309  if (neighbor->is_level_cell())
310  neighbor_subdomain_id = neighbor->level_subdomain_id();
311  // subdomain id is only valid for active cells
312  else if (neighbor->active())
313  neighbor_subdomain_id = neighbor->subdomain_id();
314 
315  const bool own_neighbor =
316  ignore_subdomain ||
317  (neighbor_subdomain_id ==
318  cell->get_triangulation().locally_owned_subdomain());
319 
320  // skip all faces between two ghost cells
321  if (!own_cell && !own_neighbor)
322  continue;
323 
324  // skip if the user doesn't want faces between own cells
325  if (own_cell && own_neighbor &&
328  continue;
329 
330  // skip face to ghost
331  if (own_cell != own_neighbor &&
332  !(flags &
334  continue;
335 
336  // Deal with refinement edges from the refined side. Assuming
337  // one-irregular meshes, this situation should only occur if
338  // both cells are active.
339  const bool periodic_neighbor =
340  cell->has_periodic_neighbor(face_no);
341 
342  if ((!periodic_neighbor &&
343  cell->neighbor_is_coarser(face_no)) ||
344  (periodic_neighbor &&
345  cell->periodic_neighbor_is_coarser(face_no)))
346  {
347  Assert(!cell->has_children(), ExcInternalError());
348  Assert(!neighbor->has_children(), ExcInternalError());
349 
350  // skip if only one processor needs to assemble the face
351  // to a ghost cell and the fine cell is not ours.
352  if (!own_cell && (flags & assemble_ghost_faces_once))
353  continue;
354 
355  const std::pair<unsigned int, unsigned int>
356  neighbor_face_no =
357  periodic_neighbor ?
358  cell->periodic_neighbor_of_coarser_periodic_neighbor(
359  face_no) :
360  cell->neighbor_of_coarser_neighbor(face_no);
361 
362  face_worker(cell,
363  face_no,
365  neighbor,
366  neighbor_face_no.first,
367  neighbor_face_no.second,
368  scratch,
369  copy);
370 
372  {
373  // If own faces are to be assembled from both sides,
374  // call the faceworker again with swapped arguments.
375  // This is because we won't be looking at an adaptively
376  // refined edge coming from the other side.
377  face_worker(neighbor,
378  neighbor_face_no.first,
379  neighbor_face_no.second,
380  cell,
381  face_no,
383  scratch,
384  copy);
385  }
386  }
387  else
388  {
389  // If iterator is active and neighbor is refined, skip
390  // internal face.
391  if (::internal::is_active_iterator(cell) &&
392  neighbor->has_children())
393  continue;
394 
395  // Now neighbor is on same level, double-check this:
396  Assert(cell->level() == neighbor->level(),
397  ExcInternalError());
398 
399  // If we own both cells only do faces from one side (unless
400  // AssembleFlags says otherwise). Here, we rely on cell
401  // comparison that will look at cell->index().
402  if (own_cell && own_neighbor &&
404  (neighbor < cell))
405  continue;
406 
407  // We only look at faces to ghost on the same level once
408  // (only where own_cell=true and own_neighbor=false)
409  if (!own_cell)
410  continue;
411 
412  // now only one processor assembles faces_to_ghost. We let
413  // the processor with the smaller (level-)subdomain id
414  // assemble the face.
415  if (own_cell && !own_neighbor &&
416  (flags & assemble_ghost_faces_once) &&
417  (neighbor_subdomain_id < current_subdomain_id))
418  continue;
419 
420  const unsigned int neighbor_face_no =
421  periodic_neighbor ?
422  cell->periodic_neighbor_face_no(face_no) :
423  cell->neighbor_face_no(face_no);
424  Assert(periodic_neighbor ||
425  neighbor->face(neighbor_face_no) ==
426  cell->face(face_no),
427  ExcInternalError());
428 
429  face_worker(cell,
430  face_no,
432  neighbor,
433  neighbor_face_no,
435  scratch,
436  copy);
437  }
438  }
439  } // faces
440 
441  // Execute the cell_worker if faces are handled before cells
442  if ((flags & cells_after_faces) &&
443  (((flags & assemble_own_cells) && own_cell) ||
444  ((flags & assemble_ghost_cells) && !own_cell)))
445  cell_worker(cell, scratch, copy);
446  };
447 
448  // Submit to workstream
449  WorkStream::run(begin,
450  end,
451  cell_action,
452  copier,
453  sample_scratch_data,
454  sample_copy_data,
455  queue_length,
456  chunk_size);
457  }
458 
528  template <class CellIteratorType,
529  class ScratchData,
530  class CopyData,
531  class CellIteratorBaseType =
533  void
535  IteratorRange<CellIteratorType> iterator_range,
536  const typename identity<std::function<
537  void(const CellIteratorBaseType &, ScratchData &, CopyData &)>>::type
538  &cell_worker,
539  const typename identity<std::function<void(const CopyData &)>>::type
540  &copier,
541 
542  const ScratchData &sample_scratch_data,
543  const CopyData & sample_copy_data,
544 
545  const AssembleFlags flags = assemble_own_cells,
546 
547  const typename identity<std::function<void(const CellIteratorBaseType &,
548  const unsigned int,
549  ScratchData &,
550  CopyData &)>>::type
551  &boundary_worker = std::function<void(const CellIteratorBaseType &,
552  const unsigned int,
553  ScratchData &,
554  CopyData &)>(),
555 
556  const typename identity<std::function<void(const CellIteratorBaseType &,
557  const unsigned int,
558  const unsigned int,
559  const CellIteratorBaseType &,
560  const unsigned int,
561  const unsigned int,
562  ScratchData &,
563  CopyData &)>>::type
564  &face_worker = std::function<void(const CellIteratorBaseType &,
565  const unsigned int,
566  const unsigned int,
567  const CellIteratorBaseType &,
568  const unsigned int,
569  const unsigned int,
570  ScratchData &,
571  CopyData &)>(),
572 
573  const unsigned int queue_length = 2 * MultithreadInfo::n_threads(),
574  const unsigned int chunk_size = 8)
575  {
576  // Call the function above
577  mesh_loop<typename IteratorRange<CellIteratorType>::IteratorOverIterators,
578  ScratchData,
579  CopyData,
580  CellIteratorBaseType>(iterator_range.begin(),
581  iterator_range.end(),
582  cell_worker,
583  copier,
584  sample_scratch_data,
585  sample_copy_data,
586  flags,
587  boundary_worker,
588  face_worker,
589  queue_length,
590  chunk_size);
591  }
592 
653  template <class CellIteratorType,
654  class ScratchData,
655  class CopyData,
656  class MainClass>
657  void
658  mesh_loop(const CellIteratorType & begin,
659  const typename identity<CellIteratorType>::type &end,
660  MainClass & main_class,
661  void (MainClass::*cell_worker)(const CellIteratorType &,
662  ScratchData &,
663  CopyData &),
664  void (MainClass::*copier)(const CopyData &),
665  const ScratchData & sample_scratch_data,
666  const CopyData & sample_copy_data,
667  const AssembleFlags flags = assemble_own_cells,
668  void (MainClass::*boundary_worker)(const CellIteratorType &,
669  const unsigned int,
670  ScratchData &,
671  CopyData &) = nullptr,
672  void (MainClass::*face_worker)(const CellIteratorType &,
673  const unsigned int,
674  const unsigned int,
675  const CellIteratorType &,
676  const unsigned int,
677  const unsigned int,
678  ScratchData &,
679  CopyData &) = nullptr,
680  const unsigned int queue_length = 2 * MultithreadInfo::n_threads(),
681  const unsigned int chunk_size = 8)
682  {
683  std::function<void(const CellIteratorType &, ScratchData &, CopyData &)>
684  f_cell_worker;
685 
686  std::function<void(
687  const CellIteratorType &, const unsigned int, ScratchData &, CopyData &)>
688  f_boundary_worker;
689 
690  std::function<void(const CellIteratorType &,
691  const unsigned int,
692  const unsigned int,
693  const CellIteratorType &,
694  const unsigned int,
695  const unsigned int,
696  ScratchData &,
697  CopyData &)>
698  f_face_worker;
699 
700  if (cell_worker != nullptr)
701  f_cell_worker = std::bind(cell_worker,
702  std::ref(main_class),
703  std::placeholders::_1,
704  std::placeholders::_2,
705  std::placeholders::_3);
706 
707  if (boundary_worker != nullptr)
708  f_boundary_worker = std::bind(boundary_worker,
709  std::ref(main_class),
710  std::placeholders::_1,
711  std::placeholders::_2,
712  std::placeholders::_3,
713  std::placeholders::_4);
714 
715  if (face_worker != nullptr)
716  f_face_worker = std::bind(face_worker,
717  std::ref(main_class),
718  std::placeholders::_1,
719  std::placeholders::_2,
720  std::placeholders::_3,
721  std::placeholders::_4,
722  std::placeholders::_5,
723  std::placeholders::_6,
724  std::placeholders::_7,
725  std::placeholders::_8);
726 
727  mesh_loop(begin,
728  end,
729  f_cell_worker,
730  std::bind(copier, main_class, std::placeholders::_1),
731  sample_scratch_data,
732  sample_copy_data,
733  flags,
734  f_boundary_worker,
735  f_face_worker,
736  queue_length,
737  chunk_size);
738  }
739 
819  template <class CellIteratorType,
820  class ScratchData,
821  class CopyData,
822  class MainClass,
823  class CellIteratorBaseType =
825  void
827  MainClass & main_class,
828  void (MainClass::*cell_worker)(const CellIteratorBaseType &,
829  ScratchData &,
830  CopyData &),
831  void (MainClass::*copier)(const CopyData &),
832  const ScratchData & sample_scratch_data,
833  const CopyData & sample_copy_data,
834  const AssembleFlags flags = assemble_own_cells,
835  void (MainClass::*boundary_worker)(const CellIteratorBaseType &,
836  const unsigned int,
837  ScratchData &,
838  CopyData &) = nullptr,
839  void (MainClass::*face_worker)(const CellIteratorBaseType &,
840  const unsigned int,
841  const unsigned int,
842  const CellIteratorBaseType &,
843  const unsigned int,
844  const unsigned int,
845  ScratchData &,
846  CopyData &) = nullptr,
847  const unsigned int queue_length = 2 * MultithreadInfo::n_threads(),
848  const unsigned int chunk_size = 8)
849  {
850  // Call the function above
851  mesh_loop<typename IteratorRange<CellIteratorType>::IteratorOverIterators,
852  ScratchData,
853  CopyData,
854  MainClass,
855  CellIteratorBaseType>(iterator_range.begin(),
856  iterator_range.end(),
857  main_class,
858  cell_worker,
859  copier,
860  sample_scratch_data,
861  sample_copy_data,
862  flags,
863  boundary_worker,
864  face_worker,
865  queue_length,
866  chunk_size);
867  }
868 } // namespace MeshWorker
869 
870 DEAL_II_NAMESPACE_CLOSE
871 
872 #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:182
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:1407
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:1167
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