include/deal.II/meshworker/loop.h

00001 //---------------------------------------------------------------------------
00002 //    @f$Id: loop.h 25345 2012-03-31 08:37:04Z bangerth @f$
00003 //
00004 //    Copyright (C) 2006, 2007, 2008, 2009, 2010, 2011, 2012 by the deal.II authors
00005 //
00006 //    This file is subject to QPL and may not be  distributed
00007 //    without copyright and license information. Please refer
00008 //    to the file deal.II/doc/license.html for the  text  and
00009 //    further information on this license.
00010 //
00011 //---------------------------------------------------------------------------
00012 
00013 #ifndef __deal2__mesh_worker_loop_h
00014 #define __deal2__mesh_worker_loop_h
00015 
00016 #include <deal.II/base/config.h>
00017 #include <deal.II/base/std_cxx1x/function.h>
00018 #include <deal.II/base/work_stream.h>
00019 #include <deal.II/base/template_constraints.h>
00020 #include <deal.II/grid/tria.h>
00021 #include <deal.II/meshworker/dof_info.h>
00022 #include <deal.II/meshworker/integration_info.h>
00023 
00024 #define DEAL_II_MESHWORKER_PARALLEL 1
00025 
00026 DEAL_II_NAMESPACE_OPEN
00027 
00028 template <typename> class TriaActiveIterator;
00029 
00030 namespace internal
00031 {
00035   template <class DI>
00036   inline bool is_active_iterator(const DI&)
00037   {
00038     return false;
00039   }
00040 
00041   template <class ACCESSOR>
00042   inline bool is_active_iterator(const TriaActiveIterator<ACCESSOR>&)
00043   {
00044     return true;
00045   }
00046 
00047   template<int dim, class DOFINFO, class A>
00048   void assemble(const MeshWorker::DoFInfoBox<dim, DOFINFO>& dinfo, A* assembler)
00049   {
00050     dinfo.assemble(*assembler);
00051   }
00052 }
00053 
00054 
00055 
00056 namespace MeshWorker
00057 {
00094   template<class INFOBOX, class DOFINFO, int dim, int spacedim, class ITERATOR>
00095   void cell_action(
00096     ITERATOR cell,
00097     DoFInfoBox<dim, DOFINFO>& dof_info,
00098     INFOBOX& info,
00099     const std_cxx1x::function<void (DOFINFO&, typename INFOBOX::CellInfo&)>& cell_worker,
00100     const std_cxx1x::function<void (DOFINFO&, typename INFOBOX::CellInfo&)>& boundary_worker,
00101     const std_cxx1x::function<void (DOFINFO&, DOFINFO&,
00102                                     typename INFOBOX::CellInfo&,
00103                                     typename INFOBOX::CellInfo&)>& face_worker,
00104     const bool cells_first,
00105     const bool unique_faces_only)
00106   {
00107     const bool integrate_cell          = (cell_worker != 0);
00108     const bool integrate_boundary      = (boundary_worker != 0);
00109     const bool integrate_interior_face = (face_worker != 0);
00110 
00111     dof_info.reset();
00112 
00113     dof_info.cell.reinit(cell);
00114     if (integrate_cell)
00115       info.cell.reinit(dof_info.cell);
00116     // Execute this, if cells
00117     // have to be dealt with
00118     // before faces
00119     if (integrate_cell && cells_first)
00120       cell_worker(dof_info.cell, info.cell);
00121 
00122     // Call the callback function in
00123     // the info box to do
00124     // computations between cell and
00125     // face action.
00126     info.post_cell(dof_info);
00127 
00128     if (integrate_interior_face || integrate_boundary)
00129       for (unsigned int face_no=0; face_no < GeometryInfo<ITERATOR::AccessorType::Container::dimension>::faces_per_cell; ++face_no)
00130       {
00131         typename ITERATOR::AccessorType::Container::face_iterator face = cell->face(face_no);
00132         if (cell->at_boundary(face_no))
00133         {
00134           if (integrate_boundary)
00135           {
00136             dof_info.interior_face_available[face_no] = true;
00137             dof_info.interior[face_no].reinit(cell, face, face_no);
00138             info.boundary.reinit(dof_info.interior[face_no]);
00139             boundary_worker(dof_info.interior[face_no], info.boundary);
00140           }
00141         }
00142         else if (integrate_interior_face)
00143         {
00144           // Interior face
00145           typename ITERATOR::AccessorType::Container::cell_iterator
00146             neighbor = cell->neighbor(face_no);
00147 
00148           // Deal with
00149           // refinement edges
00150           // from the refined
00151           // side. Assuming
00152           // one-irregular
00153           // meshes, this
00154           // situation should
00155           // only occur if
00156           // both cells are
00157           // active.
00158           if (cell->neighbor_is_coarser(face_no))
00159           {
00160             Assert(!cell->has_children(), ExcInternalError());
00161             Assert(!neighbor->has_children(), ExcInternalError());
00162 
00163             std::pair<unsigned int, unsigned int> neighbor_face_no
00164               = cell->neighbor_of_coarser_neighbor(face_no);
00165             typename ITERATOR::AccessorType::Container::face_iterator nface
00166               = neighbor->face(neighbor_face_no.first);
00167 
00168             dof_info.interior_face_available[face_no] = true;
00169             dof_info.exterior_face_available[face_no] = true;
00170             dof_info.interior[face_no].reinit(cell, face, face_no);
00171             info.face.reinit(dof_info.interior[face_no]);
00172             dof_info.exterior[face_no].reinit(
00173                 neighbor, nface, neighbor_face_no.first, neighbor_face_no.second);
00174             info.subface.reinit(dof_info.exterior[face_no]);
00175 
00176             face_worker(dof_info.interior[face_no], dof_info.exterior[face_no],
00177                 info.face, info.subface);
00178           }
00179           else
00180           {
00181             // Neighbor is
00182             // on same
00183             // level, but
00184             // only do this
00185             // from one side.
00186             if (unique_faces_only && (neighbor < cell)) continue;
00187 
00188             // If iterator
00189             // is active
00190             // and neighbor
00191             // is refined,
00192             // skip
00193             // internal face.
00194             if (internal::is_active_iterator(cell) && neighbor->has_children())
00195               continue;
00196 
00197             unsigned int neighbor_face_no = cell->neighbor_face_no(face_no);
00198             Assert (neighbor->face(neighbor_face_no) == face, ExcInternalError());
00199             // Regular interior face
00200             dof_info.interior_face_available[face_no] = true;
00201             dof_info.exterior_face_available[face_no] = true;
00202             dof_info.interior[face_no].reinit(cell, face, face_no);
00203             info.face.reinit(dof_info.interior[face_no]);
00204             dof_info.exterior[face_no].reinit(
00205                 neighbor, neighbor->face(neighbor_face_no), neighbor_face_no);
00206             info.neighbor.reinit(dof_info.exterior[face_no]);
00207 
00208             face_worker(dof_info.interior[face_no], dof_info.exterior[face_no],
00209                 info.face, info.neighbor);
00210           }
00211         }
00212       } // faces
00213     // Call the callback function in
00214     // the info box to do
00215     // computations between face and
00216     // cell action.
00217     info.post_faces(dof_info);
00218 
00219     // Execute this, if faces
00220     // have to be handled first
00221     if (integrate_cell && !cells_first)
00222       cell_worker(dof_info.cell, info.cell);
00223   }
00224 
00241   template<int dim, int spacedim, class DOFINFO, class INFOBOX, class ASSEMBLER, class ITERATOR>
00242   void loop(ITERATOR begin,
00243             typename identity<ITERATOR>::type end,
00244             DOFINFO& dinfo,
00245             INFOBOX& info,
00246             const std_cxx1x::function<void (DOFINFO&, typename INFOBOX::CellInfo&)>& cell_worker,
00247             const std_cxx1x::function<void (DOFINFO&, typename INFOBOX::CellInfo&)>& boundary_worker,
00248             const std_cxx1x::function<void (DOFINFO&, DOFINFO&,
00249                                             typename INFOBOX::CellInfo&,
00250                                             typename INFOBOX::CellInfo&)>& face_worker,
00251             ASSEMBLER& assembler,
00252             bool cells_first = true)
00253   {
00254     DoFInfoBox<dim, DOFINFO> dof_info(dinfo);
00255 
00256     assembler.initialize_info(dof_info.cell, false);
00257     for (unsigned int i=0;i<GeometryInfo<dim>::faces_per_cell;++i)
00258       {
00259         assembler.initialize_info(dof_info.interior[i], true);
00260         assembler.initialize_info(dof_info.exterior[i], true);
00261       }
00262 
00263                                      // Loop over all cells
00264 #ifdef DEAL_II_MESHWORKER_PARALLEL
00265     WorkStream::run(begin, end,
00266                     std_cxx1x::bind(&cell_action<INFOBOX, DOFINFO, dim, spacedim, ITERATOR>,
00267                                     std_cxx1x::_1, std_cxx1x::_3, std_cxx1x::_2,
00268                                     cell_worker, boundary_worker, face_worker, cells_first, true),
00269                     std_cxx1x::bind(&internal::assemble<dim,DOFINFO,ASSEMBLER>, std_cxx1x::_1, &assembler),
00270                     info, dof_info);
00271 #else
00272     for (ITERATOR cell = begin; cell != end; ++cell)
00273       {
00274         cell_action<INFOBOX,DOFINFO,dim,spacedim>(cell, dof_info,
00275                                                   info, cell_worker,
00276                                                   boundary_worker, face_worker,
00277                                                   cells_first,
00278                                                   true);
00279         dof_info.assemble(assembler);
00280       }
00281 #endif
00282   }
00283 
00290   template<int dim, int spacedim, class ITERATOR, class ASSEMBLER>
00291   void integration_loop(ITERATOR begin,
00292                         typename identity<ITERATOR>::type end,
00293                         DoFInfo<dim, spacedim>& dof_info,
00294                         IntegrationInfoBox<dim, spacedim>& box,
00295                         const std_cxx1x::function<void (DoFInfo<dim>&, IntegrationInfo<dim, spacedim>&)> &cell_worker,
00296                         const std_cxx1x::function<void (DoFInfo<dim>&, IntegrationInfo<dim, spacedim>&)> &boundary_worker,
00297                         const std_cxx1x::function<void (DoFInfo<dim>&, DoFInfo<dim>&,
00298                                                         IntegrationInfo<dim, spacedim>&,
00299                                                         IntegrationInfo<dim, spacedim>&)> &face_worker,
00300                         ASSEMBLER &assembler,
00301                         bool cells_first = true)
00302   {
00303     loop<dim, spacedim>
00304       (begin, end,
00305        dof_info,
00306        box,
00307        cell_worker,
00308        boundary_worker,
00309        face_worker,
00310        assembler,
00311        cells_first);
00312   }
00313 }
00314 
00315 DEAL_II_NAMESPACE_CLOSE
00316 
00317 #endif
00318 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

deal.II documentation generated on Tue May 22 2012 12:06:16 by doxygen 1.7.3