]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
introduce DoFInfoBox
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 16 Feb 2010 04:00:16 +0000 (04:00 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 16 Feb 2010 04:00:16 +0000 (04:00 +0000)
git-svn-id: https://svn.dealii.org/trunk@20645 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/mesh_worker_info.h
deal.II/deal.II/include/numerics/mesh_worker_info.templates.h
deal.II/deal.II/include/numerics/mesh_worker_loop.h
deal.II/examples/step-38/step-38.cc
deal.II/examples/step-39/step-39.cc

index 40d9445f59ee463904bc5b302cb11afa3938c28e..8be76611c1d8c254bd5220f2218ab6e24f026bde 100644 (file)
@@ -193,6 +193,7 @@ namespace MeshWorker
       std::vector<MatrixBlock<FullMatrix<number> > > M2;
   };
 
+  template <int dim, int spacedim> class DoFInfoBox;
 
 /**
  * A class containing information on geometry and degrees of freedom
@@ -277,8 +278,8 @@ namespace MeshWorker
                                        * the #aux_local_indices.
                                        */
       template <class DH>
-      DoFInfo(const DH& dof_handler);
-
+      DoFInfo (const DH& dof_handler);
+      
                                       /**
                                        * Set the current cell and
                                        * fill #indices.
@@ -314,7 +315,18 @@ namespace MeshWorker
       SmartPointer<const BlockInfo,DoFInfo<dim,spacedim> > block_info;
 
     private:
-                                      /// Fill index vector with active indices
+                                      /**
+                                       * Standard constructor, not
+                                       * setting any block
+                                       * indices. Use of this
+                                       * constructor is not
+                                       * recommended, but it is
+                                       * needed for the arrays in
+                                       * DoFInfoBox.
+                                       */
+      DoFInfo ();
+      
+                                      /// Fill index vector with active indices
       void get_indices(const typename DoFHandler<dim, spacedim>::cell_iterator& c);
 
                                       /// Fill index vector with level indices
@@ -332,8 +344,66 @@ namespace MeshWorker
                                        * degrees of freedom per cell.
                                        */
       BlockIndices aux_local_indices;
+      
+      friend class DoFInfoBox<dim, spacedim>;
+  };
+
+
+  /**
+ * A class bundling the MeshWorker::DoFInfo objects used on a cell.
+ *
+ * @todo Currently, we are storing an object for the cells and two for
+ * each face. We could gather all face data pertaining to the cell
+ * itself in one object, saving a bit of memory and a few operations,
+ * but sacrificing some cleanliness.
+ *
+ * @author Guido Kanschat, 2010
+ */
+  template < int dim, int spacedim=dim>
+  class DoFInfoBox
+  {
+    public:
+                                      /**
+                                       * Constructor copying the seed
+                                       * into all other objects.
+                                       */
+      DoFInfoBox(const MeshWorker::DoFInfo<dim, spacedim>& seed);
+
+                                      /**
+                                       * Reset all the availability flags.
+                                       */
+      void reset();
+      
+                                      /**
+                                       * The data for the cell.
+                                       */
+      MeshWorker::DoFInfo<dim> cell;
+                                      /**
+                                       * The data for the faces from inside.
+                                       */
+      MeshWorker::DoFInfo<dim> interior[GeometryInfo<dim>::faces_per_cell];
+                                      /**
+                                       * The data for the faces from outside.
+                                       */
+      MeshWorker::DoFInfo<dim> exterior[GeometryInfo<dim>::faces_per_cell];
+
+                                      /**
+                                       * A set of flags, indicating
+                                       * whether data on an interior
+                                       * face is available.
+                                       */
+      bool interior_face_available[GeometryInfo<dim>::faces_per_cell];
+                                      /**
+                                       * A set of flags, indicating
+                                       * whether data on an exterior
+                                       * face is available.
+                                       */
+      bool exterior_face_available[GeometryInfo<dim>::faces_per_cell];
   };
 
+
+  
+
 /**
  * Class for objects handed to local integration functions.
  *
@@ -876,7 +946,35 @@ namespace MeshWorker
     return aux_local_indices;
   }
 
+//----------------------------------------------------------------------//
+  
+  template < int dim, int spacedim>
+  inline
+  DoFInfoBox<dim, spacedim>::DoFInfoBox(const DoFInfo<dim, spacedim>& seed)
+                 :
+                 cell(seed)
+  {
+    for (unsigned int i=0;i<GeometryInfo<dim>::faces_per_cell;++i)
+      {
+       exterior[i] = seed;
+       interior[i] = seed;
+       interior_face_available[i] = false;
+       exterior_face_available[i] = false;
+      }
+  }
+
 
+  template < int dim, int spacedim>
+  inline void
+  DoFInfoBox<dim, spacedim>::reset ()
+  {
+    for (unsigned int i=0;i<GeometryInfo<dim>::faces_per_cell;++i)
+      {
+       interior_face_available[i] = false;
+       exterior_face_available[i] = false;
+      }
+  }
+  
 //----------------------------------------------------------------------//
 
   template <int dim, class FVB, int spacedim>
index 82baecbaa82c65ac45ab5f3832d002856ed529c9..3aa4d16b039ada79fe4f00b4821279e516b9f783 100644 (file)
@@ -25,6 +25,11 @@ namespace MeshWorker
   {}
 
 
+  template <int dim, int spacedim>
+  DoFInfo<dim,spacedim>::DoFInfo()
+  {}
+
+
   template <int dim, int spacedim>
   void
   DoFInfo<dim,spacedim>::get_indices(const typename DoFHandler<dim, spacedim>::cell_iterator& c)
index 02978105f00356146ce3d0f411ce43e67ce3323a..f489e073346f67c84b9527e35f1b4c41406ceda5 100644 (file)
@@ -68,14 +68,14 @@ namespace MeshWorker
  * @ingroup MeshWorker
  * @author Guido Kanschat, 2009
  */
-  template<class DOFINFO, class INFOBOX, class ITERATOR, typename ASSEMBLER>
+  template<class INFOBOX, int dim, int spacedim, typename ASSEMBLER, class ITERATOR>
   void loop(ITERATOR begin,
            typename identity<ITERATOR>::type end,
-           DOFINFO cell_dof_info,
+           DoFInfo<dim, spacedim> dinfo,
            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&,
+           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,
            ASSEMBLER &assembler,
@@ -85,12 +85,14 @@ namespace MeshWorker
     const bool integrate_boundary      = (boundary_worker != 0);
     const bool integrate_interior_face = (face_worker != 0);
 
-    ;
-    DOFINFO face_dof_info = cell_dof_info;
-    DOFINFO neighbor_dof_info = cell_dof_info;
-    assembler.initialize_info(cell_dof_info, false);
-    assembler.initialize_info(face_dof_info, true);
-    assembler.initialize_info(neighbor_dof_info, true);
+    DoFInfoBox<dim, spacedim> dof_info(dinfo);
+    
+    assembler.initialize_info(dof_info.cell, false);
+    for (unsigned int i=0;i<GeometryInfo<dim>::faces_per_cell;++i)
+      {
+       assembler.initialize_info(dof_info.interior[i], true);
+       assembler.initialize_info(dof_info.exterior[i], true);
+      }
     
                                     // Loop over all cells
     for (ITERATOR cell = begin; cell != end; ++cell)
@@ -100,10 +102,10 @@ namespace MeshWorker
                                         // before faces
        if (integrate_cell && cells_first)
          {
-           cell_dof_info.reinit(cell);
-           info.cell.reinit(cell_dof_info);
-           cell_worker(cell_dof_info, info.cell);
-           assembler.assemble(cell_dof_info);
+           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)
@@ -114,10 +116,11 @@ namespace MeshWorker
                {
                  if (integrate_boundary)
                    {
-                     cell_dof_info.reinit(cell, face, face_no);
-                     info.boundary.reinit(cell_dof_info);
-                     boundary_worker(cell_dof_info, info.boundary);
-                     assembler.assemble(cell_dof_info);
+                     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)
@@ -146,17 +149,21 @@ namespace MeshWorker
                      typename ITERATOR::AccessorType::Container::face_iterator nface
                        = neighbor->face(neighbor_face_no.first);
                      
-                     face_dof_info.reinit(cell, face, face_no);
-                     info.face.reinit(face_dof_info);
-                     neighbor_dof_info.reinit(neighbor, nface,
-                                              neighbor_face_no.first, neighbor_face_no.second);
-                     info.subface.reinit(neighbor_dof_info);
+                     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(face_dof_info, neighbor_dof_info, info.face, info.subface);
-                     assembler.assemble (face_dof_info, neighbor_dof_info);
+                     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
                    {
@@ -179,13 +186,17 @@ namespace MeshWorker
                      unsigned int neighbor_face_no = cell->neighbor_of_neighbor(face_no);
                      Assert (neighbor->face(neighbor_face_no) == face, ExcInternalError());
                                                       // Regular interior face
-                     face_dof_info.reinit(cell, face, face_no);
-                     info.face.reinit(face_dof_info);
-                     neighbor_dof_info.reinit(neighbor, neighbor->face(neighbor_face_no),
-                                              neighbor_face_no);
-                     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);
+                     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
@@ -193,10 +204,10 @@ namespace MeshWorker
                                         // have to be handled first
        if (integrate_cell && !cells_first)
          {
-           cell_dof_info.reinit(cell);
-           info.cell.reinit(cell_dof_info);
-           cell_worker(cell_dof_info, info.cell);
-           assembler.assemble (cell_dof_info);
+           dof_info.cell.reinit(cell);
+           info.cell.reinit(dof_info.cell);
+           cell_worker(dof_info.cell, info.cell);
+           assembler.assemble(dof_info.cell);
          }
       }
   }
index 43411b0cc60befe90564bc9f466e6f9c08d68b82..8eff89bddf586ce3be512f2ebccee5d2a13fef5c 100644 (file)
@@ -387,7 +387,7 @@ void Step12<dim>::assemble_system ()
                                   // result of std::bind if the local
                                   // integrators were, for example,
                                   // non-static member functions.
-  MeshWorker::loop<MeshWorker::DoFInfo<dim>, MeshWorker::IntegrationInfoBox<dim> >
+  MeshWorker::loop<MeshWorker::IntegrationInfoBox<dim>, dim, dim>
     (dof_handler.begin_active(), dof_handler.end(),
      dof_info, info_box,
      &Step12<dim>::integrate_cell_term,
index 12f257628d4a1d1a1c38c2688b4fa4d707ba63aa..c2edd8872066bed42dec81baf0a5cb8a17306b94 100644 (file)
@@ -446,7 +446,7 @@ Step39<dim>::assemble_matrix()
   MeshWorker::IntegrationInfoBox<dim> info_box;
   MeshWorker::DoFInfo<dim> dof_info(dof_handler);
   info_box.initialize(integration_worker, fe, mapping);
-  MeshWorker::loop<MeshWorker::DoFInfo<dim>, MeshWorker::IntegrationInfoBox<dim> >(
+  MeshWorker::loop<MeshWorker::IntegrationInfoBox<dim>, dim, dim>(
     dof_handler.begin_active(), dof_handler.end(),
     dof_info, info_box,
     &MatrixIntegrator<dim>::cell,
@@ -473,7 +473,7 @@ Step39<dim>::assemble_mg_matrix()
   MeshWorker::IntegrationInfoBox<dim> info_box;
   MeshWorker::DoFInfo<dim> dof_info(mg_dof_handler);
   info_box.initialize(integration_worker, fe, mapping);
-  MeshWorker::loop<MeshWorker::DoFInfo<dim>, MeshWorker::IntegrationInfoBox<dim> > (
+  MeshWorker::loop<MeshWorker::IntegrationInfoBox<dim>, dim, dim> (
     mg_dof_handler.begin(), mg_dof_handler.end(),
     dof_info, info_box,
     &MatrixIntegrator<dim>::cell,
@@ -503,7 +503,7 @@ Step39<dim>::assemble_right_hand_side()
   MeshWorker::DoFInfo<dim> dof_info(dof_handler);
   info_box.initialize(integration_worker, fe, mapping);
   
-  MeshWorker::loop<MeshWorker::DoFInfo<dim>, MeshWorker::IntegrationInfoBox<dim> >(
+  MeshWorker::loop<MeshWorker::IntegrationInfoBox<dim>, dim, dim>(
     dof_handler.begin_active(), dof_handler.end(),
     dof_info, info_box,
     &RHSIntegrator<dim>::cell,
@@ -529,7 +529,7 @@ Step39<dim>::solve()
   coarse_matrix.copy_from (mg_matrix[0]);
   MGCoarseGridHouseholder<double, Vector<double> > mg_coarse;
   mg_coarse.initialize(coarse_matrix);
-  typedef PreconditionSOR<SparseMatrix<double> > RELAXATION;
+  typedef PreconditionSSOR<SparseMatrix<double> > RELAXATION;
   MGSmootherRelaxation<SparseMatrix<double>, RELAXATION, Vector<double> >
     mg_smoother(mem);
   RELAXATION::AdditionalData smoother_data(1.);
@@ -624,7 +624,7 @@ Step39<dim>::estimate()
   MeshWorker::IntegrationInfoBox<dim> info_box;
   MeshWorker::DoFInfo<dim> dof_info(dof_handler);
   info_box.initialize(integration_worker, fe, mapping, solution_data);
-  MeshWorker::loop<MeshWorker::DoFInfo<dim>, MeshWorker::IntegrationInfoBox<dim> > (
+  MeshWorker::loop<MeshWorker::IntegrationInfoBox<dim>, dim, dim> (
     dof_handler.begin_active(), dof_handler.end(),
     dof_info, info_box,
     &Estimator<dim>::cell,

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.