Reference documentation for deal.II version 9.4.0
\(\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\}}\)
loop.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2006 - 2020 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_loop_h
18#define dealii_mesh_worker_loop_h
19
20#include <deal.II/base/config.h>
21
24
26#include <deal.II/grid/tria.h>
27
31
32#include <functional>
33
35
36// Forward declaration
37#ifndef DOXYGEN
38template <typename>
40#endif
41
42namespace internal
43{
47 template <class DI>
48 inline bool
50 {
51 return false;
52 }
53
54 template <class ACCESSOR>
55 inline bool
57 {
58 return true;
59 }
60
61 template <class ACCESSOR>
62 inline bool
64 const ::FilteredIterator<TriaActiveIterator<ACCESSOR>> &)
65 {
66 return true;
67 }
68
69 template <int dim, class DOFINFO, class A>
70 void
72 {
73 dinfo.assemble(*assembler);
74 }
75} // namespace internal
76
77
78
79namespace MeshWorker
80{
85 {
86 public:
91 : own_cells(true)
92 , ghost_cells(false)
95 , cells_first(true)
96 {}
97
102
108
115 {
127 both
128 };
129
143
154
160 };
161
162
163
189 template <class INFOBOX, class DOFINFO, int dim, int spacedim, class ITERATOR>
190 void
192 ITERATOR cell,
193 DoFInfoBox<dim, DOFINFO> &dof_info,
194 INFOBOX & info,
195 const std::function<void(DOFINFO &, typename INFOBOX::CellInfo &)>
196 &cell_worker,
197 const std::function<void(DOFINFO &, typename INFOBOX::CellInfo &)>
198 & boundary_worker,
199 const std::function<void(DOFINFO &,
200 DOFINFO &,
201 typename INFOBOX::CellInfo &,
202 typename INFOBOX::CellInfo &)> &face_worker,
203 const LoopControl & loop_control)
204 {
205 const bool ignore_subdomain =
206 (cell->get_triangulation().locally_owned_subdomain() ==
208
209 types::subdomain_id csid = (cell->is_level_cell()) ?
210 cell->level_subdomain_id() :
211 cell->subdomain_id();
212
213 const bool own_cell =
214 ignore_subdomain ||
215 (csid == cell->get_triangulation().locally_owned_subdomain());
216
217 dof_info.reset();
218
219 if ((!ignore_subdomain) && (csid == numbers::artificial_subdomain_id))
220 return;
221
222 dof_info.cell.reinit(cell);
223 dof_info.cell_valid = true;
224
225 const bool integrate_cell = (cell_worker != nullptr);
226 const bool integrate_boundary = (boundary_worker != nullptr);
227 const bool integrate_interior_face = (face_worker != nullptr);
228
229 if (integrate_cell)
230 info.cell.reinit(dof_info.cell);
231 // Execute this, if cells
232 // have to be dealt with
233 // before faces
234 if (integrate_cell && loop_control.cells_first &&
235 ((loop_control.own_cells && own_cell) ||
236 (loop_control.ghost_cells && !own_cell)))
237 cell_worker(dof_info.cell, info.cell);
238
239 // Call the callback function in
240 // the info box to do
241 // computations between cell and
242 // face action.
243 info.post_cell(dof_info);
244
245 if (integrate_interior_face || integrate_boundary)
246 for (const unsigned int face_no : cell->face_indices())
247 {
248 typename ITERATOR::AccessorType::Container::face_iterator face =
249 cell->face(face_no);
250 if (cell->at_boundary(face_no) &&
251 !cell->has_periodic_neighbor(face_no))
252 {
253 // only integrate boundary faces of own cells
254 if (integrate_boundary && own_cell)
255 {
256 dof_info.interior_face_available[face_no] = true;
257 dof_info.interior[face_no].reinit(cell, face, face_no);
258 info.boundary.reinit(dof_info.interior[face_no]);
259 boundary_worker(dof_info.interior[face_no], info.boundary);
260 }
261 }
262 else if (integrate_interior_face)
263 {
264 // Interior face
266 cell->neighbor_or_periodic_neighbor(face_no);
267
269 if (neighbor->is_level_cell())
270 neighbid = neighbor->level_subdomain_id();
271 // subdomain id is only valid for active cells
272 else if (neighbor->is_active())
273 neighbid = neighbor->subdomain_id();
274
275 const bool own_neighbor =
276 ignore_subdomain ||
277 (neighbid ==
278 cell->get_triangulation().locally_owned_subdomain());
279
280 // skip all faces between two ghost cells
281 if (!own_cell && !own_neighbor)
282 continue;
283
284 // skip if the user doesn't want faces between own cells
285 if (own_cell && own_neighbor &&
286 loop_control.own_faces == LoopControl::never)
287 continue;
288
289 // skip face to ghost
290 if (own_cell != own_neighbor &&
291 loop_control.faces_to_ghost == LoopControl::never)
292 continue;
293
294 // Deal with refinement edges from the refined side. Assuming
295 // one-irregular meshes, this situation should only occur if both
296 // cells are active.
297 const bool periodic_neighbor =
298 cell->has_periodic_neighbor(face_no);
299
300 if ((!periodic_neighbor && cell->neighbor_is_coarser(face_no)) ||
301 (periodic_neighbor &&
302 cell->periodic_neighbor_is_coarser(face_no)))
303 {
304 Assert(cell->is_active(), ExcInternalError());
305 Assert(neighbor->is_active(), ExcInternalError());
306
307 // skip if only one processor needs to assemble the face
308 // to a ghost cell and the fine cell is not ours.
309 if (!own_cell &&
310 loop_control.faces_to_ghost == LoopControl::one)
311 continue;
312
313 const std::pair<unsigned int, unsigned int> neighbor_face_no =
314 periodic_neighbor ?
315 cell->periodic_neighbor_of_coarser_periodic_neighbor(
316 face_no) :
317 cell->neighbor_of_coarser_neighbor(face_no);
318 const typename ITERATOR::AccessorType::Container::
319 face_iterator nface =
320 neighbor->face(neighbor_face_no.first);
321
322 dof_info.interior_face_available[face_no] = true;
323 dof_info.exterior_face_available[face_no] = true;
324 dof_info.interior[face_no].reinit(cell, face, face_no);
325 info.face.reinit(dof_info.interior[face_no]);
326 dof_info.exterior[face_no].reinit(neighbor,
327 nface,
328 neighbor_face_no.first,
329 neighbor_face_no.second);
330 info.subface.reinit(dof_info.exterior[face_no]);
331
332 face_worker(dof_info.interior[face_no],
333 dof_info.exterior[face_no],
334 info.face,
335 info.subface);
336 }
337 else
338 {
339 // If iterator is active and neighbor is refined, skip
340 // internal face.
342 neighbor->has_children())
343 {
344 Assert(
345 loop_control.own_faces != LoopControl::both,
347 "Assembling from both sides for own_faces is not "
348 "supported with hanging nodes!"));
349 continue;
350 }
351
352 // Now neighbor is on same level, double-check this:
353 Assert(cell->level() == neighbor->level(),
355
356 // If we own both cells only do faces from one side (unless
357 // LoopControl says otherwise). Here, we rely on cell
358 // comparison that will look at cell->index().
359 if (own_cell && own_neighbor &&
360 loop_control.own_faces == LoopControl::one &&
361 (neighbor < cell))
362 continue;
363
364 // independent of loop_control.faces_to_ghost,
365 // we only look at faces to ghost on the same level once
366 // (only where own_cell=true and own_neighbor=false)
367 if (!own_cell)
368 continue;
369
370 // now only one processor assembles faces_to_ghost. We let the
371 // processor with the smaller (level-)subdomain id assemble
372 // the face.
373 if (own_cell && !own_neighbor &&
374 loop_control.faces_to_ghost == LoopControl::one &&
375 (neighbid < csid))
376 continue;
377
378 const unsigned int neighbor_face_no =
379 periodic_neighbor ?
380 cell->periodic_neighbor_face_no(face_no) :
381 cell->neighbor_face_no(face_no);
382 Assert(periodic_neighbor ||
383 neighbor->face(neighbor_face_no) == face,
385 // Regular interior face
386 dof_info.interior_face_available[face_no] = true;
387 dof_info.exterior_face_available[face_no] = true;
388 dof_info.interior[face_no].reinit(cell, face, face_no);
389 info.face.reinit(dof_info.interior[face_no]);
390 dof_info.exterior[face_no].reinit(neighbor,
391 neighbor->face(
392 neighbor_face_no),
393 neighbor_face_no);
394 info.neighbor.reinit(dof_info.exterior[face_no]);
395
396 face_worker(dof_info.interior[face_no],
397 dof_info.exterior[face_no],
398 info.face,
399 info.neighbor);
400 }
401 }
402 } // faces
403 // Call the callback function in
404 // the info box to do
405 // computations between face and
406 // cell action.
407 info.post_faces(dof_info);
408
409 // Execute this, if faces
410 // have to be handled first
411 if (integrate_cell && !loop_control.cells_first &&
412 ((loop_control.own_cells && own_cell) ||
413 (loop_control.ghost_cells && !own_cell)))
414 cell_worker(dof_info.cell, info.cell);
415 }
416
417
432 template <int dim,
433 int spacedim,
434 class DOFINFO,
435 class INFOBOX,
436 class ASSEMBLER,
437 class ITERATOR>
438 void
439 loop(ITERATOR begin,
441 DOFINFO & dinfo,
442 INFOBOX & info,
443 const std::function<void(DOFINFO &, typename INFOBOX::CellInfo &)>
444 &cell_worker,
445 const std::function<void(DOFINFO &, typename INFOBOX::CellInfo &)>
446 & boundary_worker,
447 const std::function<void(DOFINFO &,
448 DOFINFO &,
449 typename INFOBOX::CellInfo &,
450 typename INFOBOX::CellInfo &)> &face_worker,
451 ASSEMBLER & assembler,
452 const LoopControl &lctrl = LoopControl())
453 {
454 DoFInfoBox<dim, DOFINFO> dof_info(dinfo);
455
456 assembler.initialize_info(dof_info.cell, false);
457 for (unsigned int i : GeometryInfo<dim>::face_indices())
458 {
459 assembler.initialize_info(dof_info.interior[i], true);
460 assembler.initialize_info(dof_info.exterior[i], true);
461 }
462
463 // Loop over all cells
465 begin,
466 end,
467 [&cell_worker, &boundary_worker, &face_worker, &lctrl](
468 ITERATOR cell, INFOBOX &info, DoFInfoBox<dim, DOFINFO> &dof_info) {
469 cell_action<INFOBOX, DOFINFO, dim, spacedim, ITERATOR>(cell,
470 dof_info,
471 info,
472 cell_worker,
473 boundary_worker,
474 face_worker,
475 lctrl);
476 },
477 [&assembler](const MeshWorker::DoFInfoBox<dim, DOFINFO> &dinfo) {
478 ::internal::assemble<dim, DOFINFO, ASSEMBLER>(dinfo, &assembler);
479 },
480 info,
481 dof_info);
482 }
483
484
491 template <int dim, int spacedim, class ITERATOR, class ASSEMBLER>
492 void
495 DoFInfo<dim, spacedim> & dof_info,
497 const LocalIntegrator<dim, spacedim> &integrator,
498 ASSEMBLER & assembler,
499 const LoopControl & lctrl = LoopControl())
500 {
501 std::function<void(DoFInfo<dim, spacedim> &,
503 cell_worker;
504 std::function<void(DoFInfo<dim, spacedim> &,
506 boundary_worker;
507 std::function<void(DoFInfo<dim, spacedim> &,
511 face_worker;
512 if (integrator.use_cell)
513 cell_worker =
514 [&integrator](DoFInfo<dim, spacedim> & dof_info,
515 IntegrationInfo<dim, spacedim> &integration_info) {
516 integrator.cell(dof_info, integration_info);
517 };
518 if (integrator.use_boundary)
519 boundary_worker =
520 [&integrator](DoFInfo<dim, spacedim> & dof_info,
521 IntegrationInfo<dim, spacedim> &integration_info) {
522 integrator.boundary(dof_info, integration_info);
523 };
524 if (integrator.use_face)
525 face_worker =
526 [&integrator](DoFInfo<dim, spacedim> & dof_info_1,
527 DoFInfo<dim, spacedim> & dof_info_2,
528 IntegrationInfo<dim, spacedim> &integration_info_1,
529 IntegrationInfo<dim, spacedim> &integration_info_2) {
530 integrator.face(dof_info_1,
531 dof_info_2,
532 integration_info_1,
533 integration_info_2);
534 };
535
536 loop<dim, spacedim>(begin,
537 end,
538 dof_info,
539 box,
540 cell_worker,
541 boundary_worker,
542 face_worker,
543 assembler,
544 lctrl);
545 }
546
547} // namespace MeshWorker
548
550
551#endif
DOFINFO interior[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:264
void assemble(ASSEMBLER &ass) const
Definition: dof_info.h:493
bool exterior_face_available[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:280
bool interior_face_available[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:274
DOFINFO exterior[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:268
virtual void cell(DoFInfo< dim, spacedim, number > &dinfo, IntegrationInfo< dim, spacedim > &info) const
Definition: mesh_worker.cc:81
virtual void face(DoFInfo< dim, spacedim, number > &dinfo1, DoFInfo< dim, spacedim, number > &dinfo2, IntegrationInfo< dim, spacedim > &info1, IntegrationInfo< dim, spacedim > &info2) const
Definition: mesh_worker.cc:101
virtual void boundary(DoFInfo< dim, spacedim, number > &dinfo, IntegrationInfo< dim, spacedim > &info) const
Definition: mesh_worker.cc:91
FaceOption own_faces
Definition: loop.h:153
FaceOption faces_to_ghost
Definition: loop.h:142
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
#define Assert(cond, exc)
Definition: exceptions.h:1473
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
void integration_loop(ITERATOR begin, typename identity< ITERATOR >::type end, DoFInfo< dim, spacedim > &dof_info, IntegrationInfoBox< dim, spacedim > &box, const LocalIntegrator< dim, spacedim > &integrator, ASSEMBLER &assembler, const LoopControl &lctrl=LoopControl())
Definition: loop.h:493
void loop(ITERATOR begin, typename identity< ITERATOR >::type end, DOFINFO &dinfo, 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, ASSEMBLER &assembler, const LoopControl &lctrl=LoopControl())
Definition: loop.h:439
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:191
static const char A
VectorType::value_type * end(VectorType &V)
VectorType::value_type * begin(VectorType &V)
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:1275
void assemble(const MeshWorker::DoFInfoBox< dim, DOFINFO > &dinfo, A *assembler)
Definition: loop.h:71
bool is_active_iterator(const DI &)
Definition: loop.h:49
bool is_active_iterator(const ::FilteredIterator< TriaActiveIterator< ACCESSOR > > &)
Definition: loop.h:63
const types::subdomain_id artificial_subdomain_id
Definition: types.h:298
const types::subdomain_id invalid_subdomain_id
Definition: types.h:281