Reference documentation for deal.II version 8.1.0
loop.h
1 // ---------------------------------------------------------------------
2 // @f$Id: loop.h 31932 2013-12-08 02:15:54Z heister @f$
3 //
4 // Copyright (C) 2006 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 
18 #ifndef __deal2__mesh_worker_loop_h
19 #define __deal2__mesh_worker_loop_h
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/std_cxx1x/function.h>
23 #include <deal.II/base/work_stream.h>
24 #include <deal.II/base/template_constraints.h>
25 #include <deal.II/grid/tria.h>
26 #include <deal.II/meshworker/local_integrator.h>
27 #include <deal.II/meshworker/dof_info.h>
28 #include <deal.II/meshworker/integration_info.h>
29 
30 
31 #define DEAL_II_MESHWORKER_PARALLEL 1
32 
34 
35 template <typename> class TriaActiveIterator;
36 template <typename> class FilteredIterator;
37 
38 namespace internal
39 {
43  template <class DI>
44  inline bool is_active_iterator(const DI &)
45  {
46  return false;
47  }
48 
49  template <class ACCESSOR>
51  {
52  return true;
53  }
54 
55  template <class ACCESSOR>
57  {
58  return true;
59  }
60 
61  template<int dim, class DOFINFO, class A>
62  void assemble(const MeshWorker::DoFInfoBox<dim, DOFINFO> &dinfo, A *assembler)
63  {
64  dinfo.assemble(*assembler);
65  }
66 }
67 
68 
69 
70 namespace MeshWorker
71 {
108  template<class INFOBOX, class DOFINFO, int dim, int spacedim, class ITERATOR>
110  ITERATOR cell,
111  DoFInfoBox<dim, DOFINFO> &dof_info,
112  INFOBOX &info,
113  const std_cxx1x::function<void (DOFINFO &, typename INFOBOX::CellInfo &)> &cell_worker,
114  const std_cxx1x::function<void (DOFINFO &, typename INFOBOX::CellInfo &)> &boundary_worker,
115  const std_cxx1x::function<void (DOFINFO &, DOFINFO &,
116  typename INFOBOX::CellInfo &,
117  typename INFOBOX::CellInfo &)> &face_worker,
118  const bool cells_first,
119  const bool unique_faces_only)
120  {
121  const bool integrate_cell = (cell_worker != 0);
122  const bool integrate_boundary = (boundary_worker != 0);
123  const bool integrate_interior_face = (face_worker != 0);
124 
125  dof_info.reset();
126 
127  dof_info.cell.reinit(cell);
128  if (integrate_cell)
129  info.cell.reinit(dof_info.cell);
130  // Execute this, if cells
131  // have to be dealt with
132  // before faces
133  if (integrate_cell && cells_first)
134  cell_worker(dof_info.cell, info.cell);
135 
136  // Call the callback function in
137  // the info box to do
138  // computations between cell and
139  // face action.
140  info.post_cell(dof_info);
141 
142  if (integrate_interior_face || integrate_boundary)
143  for (unsigned int face_no=0; face_no < GeometryInfo<ITERATOR::AccessorType::Container::dimension>::faces_per_cell; ++face_no)
144  {
145  typename ITERATOR::AccessorType::Container::face_iterator face = cell->face(face_no);
146  if (cell->at_boundary(face_no))
147  {
148  if (integrate_boundary)
149  {
150  dof_info.interior_face_available[face_no] = true;
151  dof_info.interior[face_no].reinit(cell, face, face_no);
152  info.boundary.reinit(dof_info.interior[face_no]);
153  boundary_worker(dof_info.interior[face_no], info.boundary);
154  }
155  }
156  else if (integrate_interior_face)
157  {
158  // Interior face
159  TriaIterator<typename ITERATOR::AccessorType> neighbor = cell->neighbor(face_no);
160 
161  // Deal with
162  // refinement edges
163  // from the refined
164  // side. Assuming
165  // one-irregular
166  // meshes, this
167  // situation should
168  // only occur if
169  // both cells are
170  // active.
171  if (cell->neighbor_is_coarser(face_no))
172  {
173  Assert(!cell->has_children(), ExcInternalError());
174  Assert(!neighbor->has_children(), ExcInternalError());
175 
176  const std::pair<unsigned int, unsigned int> neighbor_face_no
177  = cell->neighbor_of_coarser_neighbor(face_no);
178  const typename ITERATOR::AccessorType::Container::face_iterator nface
179  = neighbor->face(neighbor_face_no.first);
180 
181  dof_info.interior_face_available[face_no] = true;
182  dof_info.exterior_face_available[face_no] = true;
183  dof_info.interior[face_no].reinit(cell, face, face_no);
184  info.face.reinit(dof_info.interior[face_no]);
185  dof_info.exterior[face_no].reinit(
186  neighbor, nface, neighbor_face_no.first, neighbor_face_no.second);
187  info.subface.reinit(dof_info.exterior[face_no]);
188 
189  face_worker(dof_info.interior[face_no], dof_info.exterior[face_no],
190  info.face, info.subface);
191  }
192  else
193  {
194  // Neighbor is
195  // on same
196  // level, but
197  // only do this
198  // from one side.
199  if (unique_faces_only && (neighbor < cell)) continue;
200 
201  // If iterator
202  // is active
203  // and neighbor
204  // is refined,
205  // skip
206  // internal face.
207  if (internal::is_active_iterator(cell) && neighbor->has_children())
208  continue;
209 
210  const unsigned int neighbor_face_no = cell->neighbor_face_no(face_no);
211  Assert (neighbor->face(neighbor_face_no) == face, ExcInternalError());
212  // Regular interior face
213  dof_info.interior_face_available[face_no] = true;
214  dof_info.exterior_face_available[face_no] = true;
215  dof_info.interior[face_no].reinit(cell, face, face_no);
216  info.face.reinit(dof_info.interior[face_no]);
217  dof_info.exterior[face_no].reinit(
218  neighbor, neighbor->face(neighbor_face_no), neighbor_face_no);
219  info.neighbor.reinit(dof_info.exterior[face_no]);
220 
221  face_worker(dof_info.interior[face_no], dof_info.exterior[face_no],
222  info.face, info.neighbor);
223  }
224  }
225  } // faces
226  // Call the callback function in
227  // the info box to do
228  // computations between face and
229  // cell action.
230  info.post_faces(dof_info);
231 
232  // Execute this, if faces
233  // have to be handled first
234  if (integrate_cell && !cells_first)
235  cell_worker(dof_info.cell, info.cell);
236  }
237 
254  template<int dim, int spacedim, class DOFINFO, class INFOBOX, class ASSEMBLER, class ITERATOR>
255  void loop(ITERATOR begin,
256  typename identity<ITERATOR>::type end,
257  DOFINFO &dinfo,
258  INFOBOX &info,
259  const std_cxx1x::function<void (DOFINFO &, typename INFOBOX::CellInfo &)> &cell_worker,
260  const std_cxx1x::function<void (DOFINFO &, typename INFOBOX::CellInfo &)> &boundary_worker,
261  const std_cxx1x::function<void (DOFINFO &, DOFINFO &,
262  typename INFOBOX::CellInfo &,
263  typename INFOBOX::CellInfo &)> &face_worker,
264  ASSEMBLER &assembler,
265  bool cells_first = true,
266  bool unique_faces_only = true)
267  {
268  DoFInfoBox<dim, DOFINFO> dof_info(dinfo);
269 
270  assembler.initialize_info(dof_info.cell, false);
271  for (unsigned int i=0; i<GeometryInfo<dim>::faces_per_cell; ++i)
272  {
273  assembler.initialize_info(dof_info.interior[i], true);
274  assembler.initialize_info(dof_info.exterior[i], true);
275  }
276 
277  // Loop over all cells
278 #ifdef DEAL_II_MESHWORKER_PARALLEL
279  WorkStream::run(begin, end,
280  std_cxx1x::bind(&cell_action<INFOBOX, DOFINFO, dim, spacedim, ITERATOR>,
281  std_cxx1x::_1, std_cxx1x::_3, std_cxx1x::_2,
282  cell_worker, boundary_worker, face_worker, cells_first, unique_faces_only),
283  std_cxx1x::bind(&internal::assemble<dim,DOFINFO,ASSEMBLER>, std_cxx1x::_1, &assembler),
284  info, dof_info);
285 #else
286  for (ITERATOR cell = begin; cell != end; ++cell)
287  {
288  cell_action<INFOBOX,DOFINFO,dim,spacedim>(cell, dof_info,
289  info, cell_worker,
290  boundary_worker, face_worker,
291  cells_first,
292  unique_faces_only);
293  dof_info.assemble(assembler);
294  }
295 #endif
296  }
297 
307  template<int dim, int spacedim, class ITERATOR, class ASSEMBLER>
308  void integration_loop(ITERATOR begin,
309  typename identity<ITERATOR>::type end,
310  DoFInfo<dim, spacedim> &dof_info,
312  const std_cxx1x::function<void (DoFInfo<dim>&, IntegrationInfo<dim, spacedim>&)> &cell_worker,
313  const std_cxx1x::function<void (DoFInfo<dim>&, IntegrationInfo<dim, spacedim>&)> &boundary_worker,
314  const std_cxx1x::function<void (DoFInfo<dim> &, DoFInfo<dim> &,
316  IntegrationInfo<dim, spacedim> &)> &face_worker,
317  ASSEMBLER &assembler,
318  bool cells_first = true) DEAL_II_DEPRECATED;
319 
320 
321  template<int dim, int spacedim, class ITERATOR, class ASSEMBLER>
322  void integration_loop(ITERATOR begin,
323  typename identity<ITERATOR>::type end,
324  DoFInfo<dim, spacedim> &dof_info,
326  const std_cxx1x::function<void (DoFInfo<dim>&, IntegrationInfo<dim, spacedim>&)> &cell_worker,
327  const std_cxx1x::function<void (DoFInfo<dim>&, IntegrationInfo<dim, spacedim>&)> &boundary_worker,
328  const std_cxx1x::function<void (DoFInfo<dim> &, DoFInfo<dim> &,
330  IntegrationInfo<dim, spacedim> &)> &face_worker,
331  ASSEMBLER &assembler,
332  bool cells_first)
333  {
334  loop<dim, spacedim>
335  (begin, end,
336  dof_info,
337  box,
338  cell_worker,
339  boundary_worker,
340  face_worker,
341  assembler,
342  cells_first);
343  }
344 
345 
353  template<int dim, int spacedim, class ITERATOR, class ASSEMBLER>
354  void integration_loop(ITERATOR begin,
355  typename identity<ITERATOR>::type end,
356  DoFInfo<dim, spacedim> &dof_info,
358  const LocalIntegrator<dim, spacedim> &integrator,
359  ASSEMBLER &assembler,
360  bool cells_first = true)
361  {
362  std_cxx1x::function<void (DoFInfo<dim>&, IntegrationInfo<dim, spacedim>&)> cell_worker;
363  std_cxx1x::function<void (DoFInfo<dim>&, IntegrationInfo<dim, spacedim>&)> boundary_worker;
364  std_cxx1x::function<void (DoFInfo<dim> &, DoFInfo<dim> &,
366  IntegrationInfo<dim, spacedim> &)> face_worker;
367  if (integrator.use_cell)
368  cell_worker = std_cxx1x::bind(&LocalIntegrator<dim, spacedim>::cell, &integrator, std_cxx1x::_1, std_cxx1x::_2);
369  if (integrator.use_boundary)
370  boundary_worker = std_cxx1x::bind(&LocalIntegrator<dim, spacedim>::boundary, &integrator, std_cxx1x::_1, std_cxx1x::_2);
371  if (integrator.use_face)
372  face_worker = std_cxx1x::bind(&LocalIntegrator<dim, spacedim>::face, &integrator, std_cxx1x::_1, std_cxx1x::_2, std_cxx1x::_3, std_cxx1x::_4);
373 
374  loop<dim, spacedim>
375  (begin, end,
376  dof_info,
377  box,
378  cell_worker,
379  boundary_worker,
380  face_worker,
381  assembler,
382  cells_first);
383  }
384 
385 
386 
387 }
388 
389 DEAL_II_NAMESPACE_CLOSE
390 
391 #endif
void loop(ITERATOR begin, typename identity< ITERATOR >::type end, DOFINFO &dinfo, INFOBOX &info, const std_cxx1x::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &cell_worker, const std_cxx1x::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &boundary_worker, const std_cxx1x::function< void(DOFINFO &, DOFINFO &, typename INFOBOX::CellInfo &, typename INFOBOX::CellInfo &)> &face_worker, ASSEMBLER &assembler, bool cells_first=true, bool unique_faces_only=true)
Definition: loop.h:255
bool exterior_face_available[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:282
void cell_action(ITERATOR cell, DoFInfoBox< dim, DOFINFO > &dof_info, INFOBOX &info, const std_cxx1x::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &cell_worker, const std_cxx1x::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &boundary_worker, const std_cxx1x::function< void(DOFINFO &, DOFINFO &, typename INFOBOX::CellInfo &, typename INFOBOX::CellInfo &)> &face_worker, const bool cells_first, const bool unique_faces_only)
Definition: loop.h:109
bool interior_face_available[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:276
#define Assert(cond, exc)
Definition: exceptions.h:299
void assemble(ASSEMBLER &ass) const
Definition: dof_info.h:463
bool is_active_iterator(const DI &)
Definition: loop.h:44
BlockCompressedSparsityPattern CompressedBlockSparsityPattern DEAL_II_DEPRECATED
void run(const Iterator &begin, const typename identity< Iterator >::type &end, Worker worker, Copier copier, const ScratchData &sample_scratch_data, const CopyData &sample_copy_data, const unsigned int queue_length=2 *multithread_info.n_threads(), const unsigned int chunk_size=8)
Definition: work_stream.h:1206
DOFINFO interior[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:265
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, bool cells_first=true)
Definition: loop.h:354
DOFINFO exterior[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:269
::ExceptionBase & ExcInternalError()