]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
improve comments for periodicity stuff
authorheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 4 Dec 2013 05:30:53 +0000 (05:30 +0000)
committerheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 4 Dec 2013 05:30:53 +0000 (05:30 +0000)
git-svn-id: https://svn.dealii.org/trunk@31868 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/distributed/tria.h
deal.II/include/deal.II/grid/grid_tools.h
deal.II/source/grid/grid_tools.cc

index a4c05265295cf948b108302aa03c4b09063cca8c..1bfe32cd5af6efaff259a6aad13704f55a515723 100644 (file)
@@ -706,17 +706,20 @@ namespace parallel
 
 
       /**
-       * Join faces in the p4est forest due to periodic boundary conditions.
+       * Join faces in the p4est forest for periodic boundary conditions. As a
+       * result, each pair of faces will differ by at most one refinement level
+       * and ghost neighbors will be available across these faces.
        *
        * The vector can be filled by the function
        * GridTools::collect_periodic_faces.
        *
        * @todo At the moment just default orientation is implemented.
        *
-       * @note Before this function can be used the triangulation has to be
+       * @note Before this function can be used the Triangulation has to be
        * initialized and must not be refined.
        * Calling this function more than once is possible, but not recommended:
-       * The function destroys and rebuilds the p4est forest each time it is called.
+       * The function destroys and rebuilds the p4est forest each time it is
+       * called.
        */
       void
       add_periodicity
index 61ce67c53d81592a452e52c520d5480b75fa57bf..29cefdfe526f1312476c4f12e8215d044ef3d12b 100644 (file)
@@ -1065,8 +1065,8 @@ namespace GridTools
 
 
   /**
-   * This function will collect periodic face pairs on the highest (i.e.
-   * coarsest) mesh level.
+   * This function will collect periodic face pairs on the
+   * coarsest mesh level of the given @p container (a Triangulation or DoFHandler).
    *
    * Define a 'first' boundary as all boundary faces having boundary_id
    * @p b_id1 and a 'second' boundary consisting of all faces belonging
@@ -1095,14 +1095,14 @@ namespace GridTools
    *
    * @author Daniel Arndt, Matthias Maier, 2013
    */
-  template<typename DH>
-  std::vector<PeriodicFacePair<typename DH::cell_iterator> >
+  template<typename CONTAINER>
+  std::vector<PeriodicFacePair<typename CONTAINER::cell_iterator> >
   collect_periodic_faces
-  (const DH                 &dof_handler,
+  (const CONTAINER &container,
    const types::boundary_id b_id1,
    const types::boundary_id b_id2,
    const int                direction,
-   const dealii::Tensor<1,DH::space_dimension> &offset = dealii::Tensor<1,DH::space_dimension>());
+   const dealii::Tensor<1,CONTAINER::space_dimension> &offset = dealii::Tensor<1,CONTAINER::space_dimension>());
 
 
   /**
@@ -1116,8 +1116,7 @@ namespace GridTools
    * face with local face index <code>2*dimension+1</code> and boundary
    * indicator @p b_id.
    *
-   * This function will collect periodic face pairs on the highest (i.e.
-   * coarsest) mesh level.
+   * This function will collect periodic face pairs on the coarsest mesh level.
    *
    * @note This version of collect_periodic_face_pairs  will not work on
    * meshes with cells not in @ref GlossFaceOrientation
@@ -1125,13 +1124,13 @@ namespace GridTools
    *
    * @author Daniel Arndt, Matthias Maier, 2013
    */
-  template<typename DH>
-  std::vector<PeriodicFacePair<typename DH::cell_iterator> >
+  template<typename CONTAINER>
+  std::vector<PeriodicFacePair<typename CONTAINER::cell_iterator> >
   collect_periodic_faces
-  (const DH                 &dof_handler,
+  (const CONTAINER                 &container,
    const types::boundary_id b_id,
    const int                direction,
-   const dealii::Tensor<1,DH::space_dimension> &offset = dealii::Tensor<1,DH::space_dimension>());
+   const dealii::Tensor<1,CONTAINER::space_dimension> &offset = dealii::Tensor<1,CONTAINER::space_dimension>());
 
 
   /**
index 669dccb71ba7d2e12d19b20dae822c151ea39211..1316387b538c2a8c1f4282ff15d26c9014afe3ac 100644 (file)
@@ -2405,42 +2405,42 @@ next_cell:
 
 
 
-  template<typename DH>
-  std::vector<PeriodicFacePair<typename DH::cell_iterator> >
+  template<typename CONTAINER>
+  std::vector<PeriodicFacePair<typename CONTAINER::cell_iterator> >
   collect_periodic_faces
-    (const DH                 &dof_handler,
+    (const CONTAINER                 &container,
      const types::boundary_id b_id1,
      const types::boundary_id b_id2,
      const int                direction,
-     const dealii::Tensor<1,DH::space_dimension> &offset)
+     const dealii::Tensor<1,CONTAINER::space_dimension> &offset)
   {
-    static const int dim = DH::dimension;
-    static const int space_dim = DH::space_dimension;
+    static const int dim = CONTAINER::dimension;
+    static const int space_dim = CONTAINER::space_dimension;
     Assert (0<=direction && direction<space_dim,
             ExcIndexRange (direction, 0, space_dim));
 
     // Loop over all cells on the highest level and collect all boundary
     // faces belonging to b_id1 and b_id2:
 
-    std::set<std::pair<typename DH::cell_iterator, unsigned int> > pairs1;
-    std::set<std::pair<typename DH::cell_iterator, unsigned int> > pairs2;
+    std::set<std::pair<typename CONTAINER::cell_iterator, unsigned int> > pairs1;
+    std::set<std::pair<typename CONTAINER::cell_iterator, unsigned int> > pairs2;
 
-    for (typename DH::cell_iterator cell = dof_handler.begin(0);
-         cell != dof_handler.end(0); ++cell)
+    for (typename CONTAINER::cell_iterator cell = container.begin(0);
+         cell != container.end(0); ++cell)
       {
         for (unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell; ++i)
           {
-            const typename DH::face_iterator face = cell->face(i);
+            const typename CONTAINER::face_iterator face = cell->face(i);
             if (face->at_boundary() && face->boundary_indicator() == b_id1)
             {
-              const std::pair<typename DH::cell_iterator, unsigned int> pair1
+              const std::pair<typename CONTAINER::cell_iterator, unsigned int> pair1
                 = std::make_pair(cell, i);
               pairs1.insert(pair1);
             } 
 
             if (face->at_boundary() && face->boundary_indicator() == b_id2)
             {
-              const std::pair<typename DH::cell_iterator, unsigned int> pair2
+              const std::pair<typename CONTAINER::cell_iterator, unsigned int> pair2
                 = std::make_pair(cell, i);
               pairs2.insert(pair2);
             }
@@ -2456,15 +2456,15 @@ next_cell:
 
 
 
-  template<typename DH>
-  std::vector<PeriodicFacePair<typename DH::cell_iterator> >
-  collect_periodic_faces (const DH                 &dof_handler,
+  template<typename CONTAINER>
+  std::vector<PeriodicFacePair<typename CONTAINER::cell_iterator> >
+  collect_periodic_faces (const CONTAINER                 &container,
                           const types::boundary_id b_id,
                           const int                direction,
-                          const dealii::Tensor<1,DH::space_dimension> &offset)
+                          const dealii::Tensor<1,CONTAINER::space_dimension> &offset)
   {
-    static const int dim = DH::dimension;
-    static const int space_dim = DH::space_dimension;
+    static const int dim = CONTAINER::dimension;
+    static const int space_dim = CONTAINER::space_dimension;
     Assert (0<=direction && direction<space_dim,
             ExcIndexRange (direction, 0, space_dim));
 
@@ -2474,25 +2474,25 @@ next_cell:
     // Loop over all cells on the highest level and collect all boundary
     // faces 2*direction and 2*direction*1:
 
-    std::set<std::pair<typename DH::cell_iterator, unsigned int> > pairs1;
-    std::set<std::pair<typename DH::cell_iterator, unsigned int> > pairs2;
+    std::set<std::pair<typename CONTAINER::cell_iterator, unsigned int> > pairs1;
+    std::set<std::pair<typename CONTAINER::cell_iterator, unsigned int> > pairs2;
 
-    for (typename DH::cell_iterator cell = dof_handler.begin(0);
-         cell != dof_handler.end(0); ++cell)
+    for (typename CONTAINER::cell_iterator cell = container.begin(0);
+         cell != container.end(0); ++cell)
          {
-           const typename DH::face_iterator face_1 = cell->face(2*direction);
-           const typename DH::face_iterator face_2 = cell->face(2*direction+1);
+           const typename CONTAINER::face_iterator face_1 = cell->face(2*direction);
+           const typename CONTAINER::face_iterator face_2 = cell->face(2*direction+1);
 
            if (face_1->at_boundary() && face_1->boundary_indicator() == b_id)
            {
-             const std::pair<typename DH::cell_iterator, unsigned int> pair1
+             const std::pair<typename CONTAINER::cell_iterator, unsigned int> pair1
                = std::make_pair(cell, 2*direction);
              pairs1.insert(pair1);
            }
 
            if (face_2->at_boundary() && face_2->boundary_indicator() == b_id)
            {
-             const std::pair<typename DH::cell_iterator, unsigned int> pair2
+             const std::pair<typename CONTAINER::cell_iterator, unsigned int> pair2
                = std::make_pair(cell, 2*direction+1);
              pairs2.insert(pair2);
            }
@@ -2503,7 +2503,7 @@ next_cell:
 
     // and call match_periodic_face_pairs that does the actual matching:
 
-    typedef std::vector<PeriodicFacePair<typename DH::cell_iterator> >
+    typedef std::vector<PeriodicFacePair<typename CONTAINER::cell_iterator> >
       FaceVector;
 
     FaceVector matching = match_periodic_face_pairs(pairs1, pairs2,

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.