]> https://gitweb.dealii.org/ - dealii.git/commitdiff
bundle cell and face data in loop
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Mon, 15 Feb 2010 16:24:04 +0000 (16:24 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Mon, 15 Feb 2010 16:24:04 +0000 (16:24 +0000)
git-svn-id: https://svn.dealii.org/trunk@20633 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 8ab1652271950e2d8766f79bcf9076106df6f29f..40d9445f59ee463904bc5b302cb11afa3938c28e 100644 (file)
@@ -568,7 +568,14 @@ namespace MeshWorker
 
 /**
  * A simple container collecting the five info objects required by the
- * integration loops.
+ * integration loops. In addition to providing these objects, this
+ * class offers two functions to initialize them with reasonable
+ * contents all at once.
+ *
+ * MeshWorker::loop() requires five info objects collected into a
+ * single class. These are the data members #cell, #boundary, #face,
+ * #subface, and #neighbor. The names of those are expected by
+ * MeshWorker::loop().
  *
  * @ingroup MeshWorker
  * @author Guido Kanschat, 2009
@@ -577,16 +584,13 @@ namespace MeshWorker
   class IntegrationInfoBox
   {
     public:
-      typedef IntegrationInfo<dim, FEValuesBase<dim, spacedim>, spacedim> CellInfo;
-      typedef IntegrationInfo<dim, FEFaceValuesBase<dim, spacedim>, spacedim> FaceInfo;
 
-                                      /**
-                                       * Initialize all members
-                                       * with the same argument.
-                                       */
-      template <typename T>
-      IntegrationInfoBox(const T&);
+/// The type of the info object for cells
+      typedef IntegrationInfo<dim, FEValuesBase<dim, spacedim>, spacedim> CellInfo;
 
+/// The type of the info objects for faces.
+      typedef IntegrationInfo<dim, FEFaceValuesBase<dim, spacedim>, spacedim> FaceInfo;
+      
       template <class WORKER>
       void initialize(const WORKER&,
                      const FiniteElement<dim, spacedim>& el,
@@ -598,18 +602,20 @@ namespace MeshWorker
                      const Mapping<dim, spacedim>& mapping,
                      const NamedData<VECTOR*>& data);
 
-//      private:
-
       boost::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > cell_data;
       boost::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > boundary_data;
       boost::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > face_data;
 
-      DoFInfo<dim, spacedim> dof_info;
-      CellInfo cell_info;
-      FaceInfo boundary_info;
-      FaceInfo face_info;
-      FaceInfo subface_info;
-      FaceInfo neighbor_info;
+/// The info object for a cell
+      CellInfo cell;
+/// The info object for a boundary face
+      FaceInfo boundary;
+/// The info object for a regular interior face, seen from the first cell
+      FaceInfo face;
+/// The info object for the refined side of an interior face seen from the first cell
+      FaceInfo subface;
+/// The info object for an interior face, seen from the other cell
+      FaceInfo neighbor;
   };
 
 
@@ -925,14 +931,6 @@ namespace MeshWorker
 
 //----------------------------------------------------------------------//
 
-  template <int dim, int sdim>
-  template <typename T>
-  IntegrationInfoBox<dim,sdim>::IntegrationInfoBox(const T& t)
-                 :
-                 dof_info(t)
-  {}
-
-
   template <int dim, int sdim>
   template <class WORKER>
   void
@@ -941,15 +939,15 @@ namespace MeshWorker
             const FiniteElement<dim,sdim>& el,
             const Mapping<dim,sdim>& mapping)
   {
-    cell_info.initialize<FEValues<dim,sdim> >(el, mapping, integrator.cell_quadrature,
+    cell.initialize<FEValues<dim,sdim> >(el, mapping, integrator.cell_quadrature,
                         integrator.cell_flags);
-    boundary_info.initialize<FEFaceValues<dim,sdim> >(el, mapping, integrator.boundary_quadrature,
+    boundary.initialize<FEFaceValues<dim,sdim> >(el, mapping, integrator.boundary_quadrature,
                         integrator.boundary_flags);
-    face_info.initialize<FEFaceValues<dim,sdim> >(el, mapping, integrator.face_quadrature,
+    face.initialize<FEFaceValues<dim,sdim> >(el, mapping, integrator.face_quadrature,
                         integrator.face_flags);
-    subface_info.initialize<FESubfaceValues<dim,sdim> >(el, mapping, integrator.face_quadrature,
+    subface.initialize<FESubfaceValues<dim,sdim> >(el, mapping, integrator.face_quadrature,
                            integrator.face_flags);
-    neighbor_info.initialize<FEFaceValues<dim,sdim> >(el, mapping, integrator.face_quadrature,
+    neighbor.initialize<FEFaceValues<dim,sdim> >(el, mapping, integrator.face_quadrature,
                                                      integrator.neighbor_flags);
   }
 
@@ -969,19 +967,19 @@ namespace MeshWorker
     p = boost::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (integrator.cell_selector));
     p->initialize(data);
     cell_data = p;
-    cell_info.initialize_data(p);
+    cell.initialize_data(p);
 
     p = boost::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (integrator.boundary_selector));
     p->initialize(data);
     boundary_data = p;
-    boundary_info.initialize_data(p);
+    boundary.initialize_data(p);
 
     p = boost::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (integrator.face_selector));
     p->initialize(data);
     face_data = p;
-    face_info.initialize_data(p);
-    subface_info.initialize_data(p);
-    neighbor_info.initialize_data(p);
+    face.initialize_data(p);
+    subface.initialize_data(p);
+    neighbor.initialize_data(p);
   }
 }
 
index b48f92cc374b069244a057ac12479f8c1cb404ce..02978105f00356146ce3d0f411ce43e67ce3323a 100644 (file)
@@ -24,27 +24,29 @@ DEAL_II_NAMESPACE_OPEN
 
 template <typename> class TriaActiveIterator;
 
-namespace MeshWorker
+namespace internal
 {
-  namespace internal
-  {
 /**
  * Find out if an iterator supports inactive cells.
  */
-    template <class DI>
-    inline bool is_active_iterator(const DI&)
-    {
-      return false;
-    }
-
-    template <class ACCESSOR>
-    inline bool is_active_iterator(const TriaActiveIterator<ACCESSOR>&)
-    {
-      return true;
-    }
+  template <class DI>
+  inline bool is_active_iterator(const DI&)
+  {
+    return false;
+  }
+  
+  template <class ACCESSOR>
+  inline bool is_active_iterator(const TriaActiveIterator<ACCESSOR>&)
+  {
+    return true;
   }
 
+}
+
+
 
+namespace MeshWorker
+{
 /**
  * The main work function of this namespace. Its action consists of two
  * loops.
@@ -66,15 +68,16 @@ namespace MeshWorker
  * @ingroup MeshWorker
  * @author Guido Kanschat, 2009
  */
-  template<class DOFINFO, class CELLINFO, class FACEINFO, class ITERATOR, typename ASSEMBLER>
+  template<class DOFINFO, class INFOBOX, class ITERATOR, typename ASSEMBLER>
   void loop(ITERATOR begin,
            typename identity<ITERATOR>::type end,
            DOFINFO cell_dof_info,
-           CELLINFO& cell_info, FACEINFO& boundary_info,
-           FACEINFO& face_info, FACEINFO& subface_info, FACEINFO& neighbor_info,
-           const std_cxx1x::function<void (DOFINFO&, CELLINFO &)> &cell_worker,
-           const std_cxx1x::function<void (DOFINFO&, FACEINFO &)> &boundary_worker,
-           const std_cxx1x::function<void (DOFINFO&, DOFINFO&, FACEINFO &, FACEINFO &)> &face_worker,
+           INFOBOX& info,
+           const std_cxx1x::function<void (DOFINFO&, typename INFOBOX::CellInfo &)> &cell_worker,
+           const std_cxx1x::function<void (DOFINFO&, typename INFOBOX::FaceInfo &)> &boundary_worker,
+           const std_cxx1x::function<void (DOFINFO&, DOFINFO&,
+                                           typename INFOBOX::FaceInfo &,
+                                           typename INFOBOX::FaceInfo &)> &face_worker,
            ASSEMBLER &assembler,
            bool cells_first = true)
   {
@@ -98,8 +101,8 @@ namespace MeshWorker
        if (integrate_cell && cells_first)
          {
            cell_dof_info.reinit(cell);
-           cell_info.reinit(cell_dof_info);
-           cell_worker(cell_dof_info, cell_info);
+           info.cell.reinit(cell_dof_info);
+           cell_worker(cell_dof_info, info.cell);
            assembler.assemble(cell_dof_info);
          }
 
@@ -112,8 +115,8 @@ namespace MeshWorker
                  if (integrate_boundary)
                    {
                      cell_dof_info.reinit(cell, face, face_no);
-                     boundary_info.reinit(cell_dof_info);
-                     boundary_worker(cell_dof_info, boundary_info);
+                     info.boundary.reinit(cell_dof_info);
+                     boundary_worker(cell_dof_info, info.boundary);
                      assembler.assemble(cell_dof_info);
                    }
                }
@@ -144,15 +147,15 @@ namespace MeshWorker
                        = neighbor->face(neighbor_face_no.first);
                      
                      face_dof_info.reinit(cell, face, face_no);
-                     face_info.reinit(face_dof_info);
+                     info.face.reinit(face_dof_info);
                      neighbor_dof_info.reinit(neighbor, nface,
                                               neighbor_face_no.first, neighbor_face_no.second);
-                     subface_info.reinit(neighbor_dof_info);
+                     info.subface.reinit(neighbor_dof_info);
                                                       // Neighbor
                                                       // first to
                                                       // conform to
                                                       // old version
-                     face_worker(face_dof_info, neighbor_dof_info, face_info, subface_info);
+                     face_worker(face_dof_info, neighbor_dof_info, info.face, info.subface);
                      assembler.assemble (face_dof_info, neighbor_dof_info);
                    }
                  else
@@ -177,11 +180,11 @@ namespace MeshWorker
                      Assert (neighbor->face(neighbor_face_no) == face, ExcInternalError());
                                                       // Regular interior face
                      face_dof_info.reinit(cell, face, face_no);
-                     face_info.reinit(face_dof_info);
+                     info.face.reinit(face_dof_info);
                      neighbor_dof_info.reinit(neighbor, neighbor->face(neighbor_face_no),
                                               neighbor_face_no);
-                     neighbor_info.reinit(neighbor_dof_info);
-                     face_worker(face_dof_info, neighbor_dof_info, face_info, neighbor_info);
+                     info.neighbor.reinit(neighbor_dof_info);
+                     face_worker(face_dof_info, neighbor_dof_info, info.face, info.neighbor);
                      assembler.assemble(face_dof_info, neighbor_dof_info);
                    }
                }
@@ -191,8 +194,8 @@ namespace MeshWorker
        if (integrate_cell && !cells_first)
          {
            cell_dof_info.reinit(cell);
-           cell_info.reinit(cell_dof_info);
-           cell_worker(cell_dof_info, cell_info);
+           info.cell.reinit(cell_dof_info);
+           cell_worker(cell_dof_info, info.cell);
            assembler.assemble (cell_dof_info);
          }
       }
@@ -207,6 +210,7 @@ namespace MeshWorker
   template<class CELLINFO, class FACEINFO, int dim, class ITERATOR, typename ASSEMBLER>
   void integration_loop(ITERATOR begin,
                        typename identity<ITERATOR>::type end,
+                       DoFInfo<dim>& dof_info,
                        IntegrationInfoBox<dim, dim>& box,
                        const std_cxx1x::function<void (DoFInfo<dim>&, CELLINFO &)> &cell_worker,
                        const std_cxx1x::function<void (DoFInfo<dim>&, FACEINFO &)> &boundary_worker,
@@ -214,12 +218,10 @@ namespace MeshWorker
                        ASSEMBLER &assembler,
                        bool cells_first = true)
   {
-    loop<DoFInfo<dim>,CELLINFO,FACEINFO>
+    loop<DoFInfo<dim>,IntegrationInfoBox<dim, dim> >
       (begin, end,
-       box.dof_info,
-       box.cell_info, box.boundary_info,
-       box.face_info,
-       box.subface_info, box.neighbor_info,
+       dof_info,
+       box,
        cell_worker,
        boundary_worker,
        face_worker,

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.