]> https://gitweb.dealii.org/ - dealii.git/commitdiff
reduce MeshWorker::loop() to a simple loop suitable for WorkStream
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Mon, 22 Feb 2010 20:55:26 +0000 (20:55 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Mon, 22 Feb 2010 20:55:26 +0000 (20:55 +0000)
git-svn-id: https://svn.dealii.org/trunk@20673 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/mesh_worker_info.h
deal.II/deal.II/include/numerics/mesh_worker_loop.h

index 8be76611c1d8c254bd5220f2218ab6e24f026bde..54c2e54564847de85b10ba966d2443cebc9da7a6 100644 (file)
@@ -373,6 +373,19 @@ namespace MeshWorker
                                        * Reset all the availability flags.
                                        */
       void reset();
+
+                                      /**
+                                       * After all info objects have
+                                       * been filled appropriately,
+                                       * use the ASSEMBLER object
+                                       * to assemble them into the
+                                       * global data. See
+                                       * MeshWorker::Assembler for
+                                       * available classes.
+                                       */
+      template <class ASSEMBLER>
+      void assemble(ASSEMBLER& ass) const;
+      
       
                                       /**
                                        * The data for the cell.
@@ -974,7 +987,31 @@ namespace MeshWorker
        exterior_face_available[i] = false;
       }
   }
+
+  
+  template < int dim, int spacedim>
+  template <class ASSEMBLER>
+  inline void
+  DoFInfoBox<dim, spacedim>::assemble (ASSEMBLER& assembler) const
+  {
+    assembler.assemble(cell);
+    for (unsigned int i=0;i<GeometryInfo<dim>::faces_per_cell;++i)
+      {
+                                        // Only do something if data available
+       if (interior_face_available[i])
+         {
+                                            // If both data
+                                            // available, it is an
+                                            // interior face
+           if (exterior_face_available[i])
+             assembler.assemble(interior[i], exterior[i]);
+           else
+             assembler.assemble(interior[i]);
+         }
+      }
+  }
   
+
 //----------------------------------------------------------------------//
 
   template <int dim, class FVB, int spacedim>
index f489e073346f67c84b9527e35f1b4c41406ceda5..05023e5e9ad776e56eebdebe5d92a91461f9df07 100644 (file)
@@ -47,6 +47,157 @@ namespace internal
 
 namespace MeshWorker
 {
+/**
+ * The function called by loop() to perform the required actions on a
+ * cell and its faces. The three functions <tt>cell_worker</tt>,
+ * <tt>boundary_worker</tt> and <tt>face_worker</tt> are the same ones
+ * handed to loop(). While there we only run the loop over all cells,
+ * here, we do a single cell and, if necessary, its faces, interior
+ * and boundary.
+ *
+ * Upon return, the DoFInfo objects in the DoFInfoBox are filled with
+ * the data computed on the cell and each of the faces. Thus, after
+ * the execution of this function, we are ready to call
+ * DoFInfoBox::assemble() to distribute the local data into global
+ * data.
+ *
+ * @param cell is the cell we work on
+ * @param dof_info is the object into which local results are
+ * entered. It is expected to have been set up for the right types of
+ * data.
+ * @param info is the object containing additional data only needed
+ * for internal processing.
+ * @param cell_worker defines the local action on each cell.
+ * @param boundary_worker defines the local action on boundary faces
+ * @param face_worker defines the local action on interior faces.
+ * @param cells_first determines, whether faces or cells are to be
+ * dealt with first.
+ *
+ * @author Guido Kanschat, 2010
+ */
+  template<class INFOBOX, int dim, int spacedim, class ITERATOR>
+  void cell_action(
+    ITERATOR cell,
+    DoFInfoBox<dim, spacedim>& dof_info,
+    INFOBOX& info,
+    const std_cxx1x::function<void (DoFInfo<dim, spacedim>&, typename INFOBOX::CellInfo &)> &cell_worker,
+    const std_cxx1x::function<void (DoFInfo<dim, spacedim>&, typename INFOBOX::FaceInfo &)> &boundary_worker,
+    const std_cxx1x::function<void (DoFInfo<dim, spacedim>&, DoFInfo<dim, spacedim>&,
+                                   typename INFOBOX::FaceInfo &,
+                                   typename INFOBOX::FaceInfo &)> &face_worker,
+    bool cells_first = true)
+  {
+    const bool integrate_cell          = (cell_worker != 0);
+    const bool integrate_boundary      = (boundary_worker != 0);
+    const bool integrate_interior_face = (face_worker != 0);
+    
+    dof_info.reset();
+                                    // Execute this, if cells
+                                    // have to be dealt with
+                                    // before faces
+    if (integrate_cell && cells_first)
+      {
+       dof_info.cell.reinit(cell);
+       info.cell.reinit(dof_info.cell);
+       cell_worker(dof_info.cell, info.cell);
+      }
+    
+    if (integrate_interior_face || integrate_boundary)
+      for (unsigned int face_no=0; face_no < GeometryInfo<ITERATOR::AccessorType::Container::dimension>::faces_per_cell; ++face_no)
+       {
+         typename ITERATOR::AccessorType::Container::face_iterator face = cell->face(face_no);
+         if (cell->at_boundary(face_no))
+           {
+             if (integrate_boundary)
+               {
+                 dof_info.interior_face_available[face_no] = true;
+                 dof_info.interior[face_no].reinit(cell, face, face_no);
+                 info.boundary.reinit(dof_info.interior[face_no]);
+                 boundary_worker(dof_info.interior[face_no], info.boundary);
+               }
+           }
+         else if (integrate_interior_face)
+           {
+                                              // Interior face
+             typename ITERATOR::AccessorType::Container::cell_iterator
+               neighbor = cell->neighbor(face_no);
+             
+                                              // Deal with
+                                              // refinement edges
+                                              // from the refined
+                                              // side. Assuming
+                                              // one-irregular
+                                              // meshes, this
+                                              // situation should
+                                              // only occur if
+                                              // both cells are
+                                              // active.
+             if (neighbor->level() < cell->level())
+               {
+                 Assert(!cell->has_children(), ExcInternalError());
+                 Assert(!neighbor->has_children(), ExcInternalError());
+                 
+                 std::pair<unsigned int, unsigned int> neighbor_face_no
+                   = cell->neighbor_of_coarser_neighbor(face_no);
+                 typename ITERATOR::AccessorType::Container::face_iterator nface
+                   = neighbor->face(neighbor_face_no.first);
+                 
+                 dof_info.interior_face_available[face_no] = true;
+                 dof_info.exterior_face_available[face_no] = true;
+                 dof_info.interior[face_no].reinit(cell, face, face_no);
+                 info.face.reinit(dof_info.interior[face_no]);
+                 dof_info.exterior[face_no].reinit(
+                   neighbor, nface, neighbor_face_no.first, neighbor_face_no.second);
+                 info.subface.reinit(dof_info.exterior[face_no]);
+                 
+                 face_worker(dof_info.interior[face_no], dof_info.exterior[face_no],
+                             info.face, info.subface);
+               }
+             else
+               {
+                                                  // Neighbor is
+                                                  // on same
+                                                  // level, but
+                                                  // only do this
+                                                  // from one side.
+                 if (neighbor < cell) continue;
+                 
+                                                  // If iterator
+                                                  // is active
+                                                  // and neighbor
+                                                  // is refined,
+                                                  // skip
+                                                  // internal face.
+                 if (internal::is_active_iterator(cell) && neighbor->has_children())
+                   continue;
+                 
+                 unsigned int neighbor_face_no = cell->neighbor_of_neighbor(face_no);
+                 Assert (neighbor->face(neighbor_face_no) == face, ExcInternalError());
+                                                  // Regular interior face
+                 dof_info.interior_face_available[face_no] = true;
+                 dof_info.exterior_face_available[face_no] = true;
+                 dof_info.interior[face_no].reinit(cell, face, face_no);
+                 info.face.reinit(dof_info.interior[face_no]);
+                 dof_info.exterior[face_no].reinit(
+                   neighbor, neighbor->face(neighbor_face_no), neighbor_face_no);
+                 info.neighbor.reinit(dof_info.exterior[face_no]);
+                 
+                 face_worker(dof_info.interior[face_no], dof_info.exterior[face_no],
+                             info.face, info.neighbor);
+               }
+           }
+       } // faces
+                                    // Execute this, if faces
+                                    // have to be handled first
+    if (integrate_cell && !cells_first)
+      {
+       dof_info.cell.reinit(cell);
+       info.cell.reinit(dof_info.cell);
+       cell_worker(dof_info.cell, info.cell);
+      }  
+  }
+  
+  
 /**
  * The main work function of this namespace. Its action consists of two
  * loops.
@@ -81,10 +232,6 @@ namespace MeshWorker
            ASSEMBLER &assembler,
            bool cells_first = true)
   {
-    const bool integrate_cell          = (cell_worker != 0);
-    const bool integrate_boundary      = (boundary_worker != 0);
-    const bool integrate_interior_face = (face_worker != 0);
-
     DoFInfoBox<dim, spacedim> dof_info(dinfo);
     
     assembler.initialize_info(dof_info.cell, false);
@@ -97,118 +244,8 @@ namespace MeshWorker
                                     // Loop over all cells
     for (ITERATOR cell = begin; cell != end; ++cell)
       {
-                                        // Execute this, if cells
-                                        // have to be dealt with
-                                        // before faces
-       if (integrate_cell && cells_first)
-         {
-           dof_info.cell.reinit(cell);
-           info.cell.reinit(dof_info.cell);
-           cell_worker(dof_info.cell, info.cell);
-           assembler.assemble(dof_info.cell);
-         }
-
-       if (integrate_interior_face || integrate_boundary)
-         for (unsigned int face_no=0; face_no < GeometryInfo<ITERATOR::AccessorType::Container::dimension>::faces_per_cell; ++face_no)
-           {
-             typename ITERATOR::AccessorType::Container::face_iterator face = cell->face(face_no);
-             if (cell->at_boundary(face_no))
-               {
-                 if (integrate_boundary)
-                   {
-                     dof_info.interior_face_available[face_no] = true;
-                     dof_info.interior[face_no].reinit(cell, face, face_no);
-                     info.boundary.reinit(dof_info.interior[face_no]);
-                     boundary_worker(dof_info.interior[face_no], info.boundary);
-                     assembler.assemble(dof_info.interior[face_no]);
-                   }
-               }
-             else if (integrate_interior_face)
-               {
-                                                  // Interior face
-                 typename ITERATOR::AccessorType::Container::cell_iterator
-                   neighbor = cell->neighbor(face_no);
-
-                                                  // Deal with
-                                                  // refinement edges
-                                                  // from the refined
-                                                  // side. Assuming
-                                                  // one-irregular
-                                                  // meshes, this
-                                                  // situation should
-                                                  // only occur if
-                                                  // both cells are
-                                                  // active.
-                 if (neighbor->level() < cell->level())
-                   {
-                     Assert(!cell->has_children(), ExcInternalError());
-                     Assert(!neighbor->has_children(), ExcInternalError());
-
-                     std::pair<unsigned int, unsigned int> neighbor_face_no
-                       = cell->neighbor_of_coarser_neighbor(face_no);
-                     typename ITERATOR::AccessorType::Container::face_iterator nface
-                       = neighbor->face(neighbor_face_no.first);
-                     
-                     dof_info.interior_face_available[face_no] = true;
-                     dof_info.exterior_face_available[face_no] = true;
-                     dof_info.interior[face_no].reinit(cell, face, face_no);
-                     info.face.reinit(dof_info.interior[face_no]);
-                     dof_info.exterior[face_no].reinit(
-                       neighbor, nface, neighbor_face_no.first, neighbor_face_no.second);
-                     info.subface.reinit(dof_info.exterior[face_no]);
-                     
-                                                      // Neighbor
-                                                      // first to
-                                                      // conform to
-                                                      // old version
-                     face_worker(dof_info.interior[face_no], dof_info.exterior[face_no],
-                                 info.face, info.subface);
-                     assembler.assemble (dof_info.interior[face_no], dof_info.exterior[face_no]);
-                   }
-                 else
-                   {
-                                                      // Neighbor is
-                                                      // on same
-                                                      // level, but
-                                                      // only do this
-                                                      // from one side.
-                     if (neighbor < cell) continue;
-                     
-                                                      // If iterator
-                                                      // is active
-                                                      // and neighbor
-                                                      // is refined,
-                                                      // skip
-                                                      // internal face.
-                     if (internal::is_active_iterator(cell) && neighbor->has_children())
-                       continue;
-
-                     unsigned int neighbor_face_no = cell->neighbor_of_neighbor(face_no);
-                     Assert (neighbor->face(neighbor_face_no) == face, ExcInternalError());
-                                                      // Regular interior face
-                     dof_info.interior_face_available[face_no] = true;
-                     dof_info.exterior_face_available[face_no] = true;
-                     dof_info.interior[face_no].reinit(cell, face, face_no);
-                     info.face.reinit(dof_info.interior[face_no]);
-                     dof_info.exterior[face_no].reinit(
-                       neighbor, neighbor->face(neighbor_face_no), neighbor_face_no);
-                     info.neighbor.reinit(dof_info.exterior[face_no]);
-                     
-                     face_worker(dof_info.interior[face_no], dof_info.exterior[face_no],
-                                 info.face, info.neighbor);
-                     assembler.assemble(dof_info.interior[face_no], dof_info.exterior[face_no]);
-                   }
-               }
-           } // faces
-                                        // Execute this, if faces
-                                        // have to be handled first
-       if (integrate_cell && !cells_first)
-         {
-           dof_info.cell.reinit(cell);
-           info.cell.reinit(dof_info.cell);
-           cell_worker(dof_info.cell, info.cell);
-           assembler.assemble(dof_info.cell);
-         }
+       cell_action(cell, dof_info, info, cell_worker, boundary_worker, face_worker, cells_first);
+       dof_info.assemble(assembler);
       }
   }
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.