]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Reintroduce the old make_periodicity_constraints interface as well (because this...
authormaier <maier@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 1 Dec 2012 21:25:44 +0000 (21:25 +0000)
committermaier <maier@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 1 Dec 2012 21:25:44 +0000 (21:25 +0000)
git-svn-id: https://svn.dealii.org/trunk@27722 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/dofs/dof_tools.h
deal.II/include/deal.II/grid/grid_tools.h
deal.II/source/dofs/dof_tools.cc
deal.II/source/dofs/dof_tools.inst.in
deal.II/source/grid/grid_tools.cc
deal.II/source/grid/grid_tools.inst.in

index 5a9b032536ba87bbc4263ecee2007b1695c6c4ab..53ebedd75b2543d972e34bfde39784355416c6ad 100644 (file)
@@ -1163,6 +1163,61 @@ namespace DoFTools
                                 dealii::Tensor<1,DH::space_dimension> &offset,
                                 dealii::ConstraintMatrix &constraint_matrix,
                                 const ComponentMask      &component_mask = ComponentMask());
+
+
+  /**
+   * This compatibility version of make_periodicity_constraints only works
+   * on grids with cells in @ref GlossFaceOrientation "standard orientation".
+   *
+   * Instead of defining a 'first' and 'second' boundary with the help of
+   * two boundary_indicators this function defines a 'left' boundary as all
+   * faces with local face index <code>2*dimension</code> and boundary
+   * indicator @p b_id and, similarly, a 'right' boundary consisting of all
+   * face with local face index <code>2*dimension+1</code> and boundary
+   * indicator @p b_id.
+   *
+   * @note This version of make_periodicity_constraints  will not work on
+   * meshes with cells not in @ref GlossFaceOrientation
+   * "standard orientation".
+   *
+   * @note This function will not work for DoFHandler objects that are
+   * built on a parallel::distributed::Triangulation object.
+   */
+  template<typename DH>
+  void
+  make_periodicity_constraints (const DH                 &dof_handler,
+                                const types::boundary_id b_id,
+                                const int                direction,
+                                dealii::ConstraintMatrix &constraint_matrix,
+                                const ComponentMask      &component_mask = ComponentMask());
+
+
+  /**
+   * Same as above but with an optional argument @p offset.
+   *
+   * The @p offset is a vector tangential to the faces that is added to
+   * the location of vertices of the 'first' boundary when attempting to
+   * match them to the corresponding vertices of the 'second' boundary via
+   * @p orthogonal_equality. This can be used to implement conditions such
+   * as $u(0,y)=u(1,y+1)$.
+   *
+   * @note This version of make_periodicity_constraints  will not work on
+   * meshes with cells not in @ref GlossFaceOrientation
+   * "standard orientation".
+   *
+   * @note This function will not work for DoFHandler objects that are
+   * built on a parallel::distributed::Triangulation object.
+   */
+  template<typename DH>
+  void
+  make_periodicity_constraints (const DH                 &dof_handler,
+                                const types::boundary_id b_id,
+                                const int                direction,
+                                dealii::Tensor<1,DH::space_dimension> &offset,
+                                dealii::ConstraintMatrix &constraint_matrix,
+                                const ComponentMask      &component_mask = ComponentMask());
+
+
   //@}
 
 
index 8be9ef13e97e7516e37d9e1bba4aca211ecfd9f6..6abe79967948528289eec8a62288b9680eff34ac 100644 (file)
@@ -1046,6 +1046,34 @@ namespace GridTools
                                const dealii::Tensor<1,DH::space_dimension> &offset);
 
 
+  /**
+   * This compatibility version of collect_periodic_face_pairs only works
+   * on grids with cells in @ref GlossFaceOrientation "standard orientation".
+   *
+   * Instead of defining a 'first' and 'second' boundary with the help of
+   * two boundary_indicators this function defines a 'left' boundary as all
+   * faces with local face index <code>2*dimension</code> and boundary
+   * indicator @p b_id and, similarly, a 'right' boundary consisting of all
+   * 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.
+   *
+   * @note This version of collect_periodic_face_pairs  will not work on
+   * meshes with cells not in @ref GlossFaceOrientation
+   * "standard orientation".
+   *
+   * @author Matthias Maier, 2012
+   */
+  template<typename DH>
+  std::map<typename DH::face_iterator, typename DH::face_iterator>
+  collect_periodic_face_pairs (const DH                 &dof_handler, /*TODO: Name*/
+                               const types::boundary_id b_id,
+                               int                      direction,
+                               const dealii::Tensor<1,DH::space_dimension> &offset);
+
+
 
   /**
    * Exception
index e6133255f69c9090d510df9810bfcc394924e054..874a78c99d7f704f2fe15617e213dc17e0ccb432 100644 (file)
@@ -3556,6 +3556,88 @@ namespace DoFTools
 
 
 
+  template<typename DH>
+  void
+  make_periodicity_constraints (const DH                       &dof_handler,
+                                const types::boundary_id       b_id,
+                                const int                      direction,
+                                dealii::ConstraintMatrix       &constraint_matrix,
+                                const ComponentMask            &component_mask)
+  {
+    Tensor<1,DH::space_dimension> dummy;
+    make_periodicity_constraints (dof_handler,
+                                  b_id,
+                                  direction,
+                                  dummy,
+                                  constraint_matrix,
+                                  component_mask);
+  }
+
+
+
+  template<typename DH>
+  void
+  make_periodicity_constraints (const DH                  &dof_handler,
+                                const types::boundary_id  b_id,
+                                const int                 direction,
+                                dealii::Tensor<1,DH::space_dimension> &offset,
+                                dealii::ConstraintMatrix  &constraint_matrix,
+                                const ComponentMask       &component_mask)
+  {
+    static const int dim = DH::dimension;
+    static const int space_dim = DH::space_dimension;
+
+    Assert (0<=direction && direction<space_dim,
+            ExcIndexRange (direction, 0, space_dim));
+
+    Assert(dim == space_dim,
+           ExcNotImplemented());
+
+    Assert ((dynamic_cast<const parallel::distributed::Triangulation<DH::dimension,DH::space_dimension>*>
+             (&dof_handler.get_tria())
+             ==
+             0),
+            ExcMessage ("This function can not be used with distributed triangulations."
+                        "See the documentation for more information."));
+
+    typedef typename DH::face_iterator FaceIterator;
+    typedef std::map<FaceIterator, FaceIterator> FaceMap;
+
+    // Collect matching periodic cells on the coarsest level:
+    FaceMap matched_cells =
+      GridTools::collect_periodic_face_pairs(dof_handler,
+                                             b_id,
+                                             direction, offset);
+
+    // And apply the low level make_periodicity_constraints function to
+    // every matching pair:
+    for (typename FaceMap::iterator it = matched_cells.begin();
+         it != matched_cells.end(); ++it)
+      {
+        typedef typename DH::face_iterator FaceIterator;
+        const FaceIterator &face_1 = it->first;
+        const FaceIterator &face_2 = it->second;
+
+        Assert(face_1->at_boundary() && face_2->at_boundary(),
+               ExcInternalError());
+
+        Assert (face_1->boundary_indicator() == b_id &&
+                face_2->boundary_indicator() == b_id,
+                ExcInternalError());
+
+        Assert (face_1 != face_2,
+                ExcInternalError());
+
+        make_periodicity_constraints(face_1,
+                                     face_2,
+                                     constraint_matrix,
+                                     component_mask
+                                     /* standard orientation */);
+      }
+  }
+
+
+
   namespace internal
   {
     // return an array that for each dof on the reference cell
index 7aa2cb456a8d38c237236d2b71f21b24ebbfacc8..02859820284c7e57bcea717b9c01a7766e359ef1 100644 (file)
@@ -347,6 +347,23 @@ for (DH : DOFHANDLERS; deal_II_dimension : DIMENSIONS)
                                          dealii::Tensor<1,DH::space_dimension> &,
                                          dealii::ConstraintMatrix &,
                                          const ComponentMask &);
+
+  template
+  void
+  DoFTools::make_periodicity_constraints(const DH &,
+                                         const types::boundary_id,
+                                         const int,
+                                         dealii::ConstraintMatrix &,
+                                         const ComponentMask &);
+
+  template
+  void
+  DoFTools::make_periodicity_constraints(const DH &,
+                                         const types::boundary_id,
+                                         const int,
+                                         dealii::Tensor<1,DH::space_dimension> &,
+                                         dealii::ConstraintMatrix &,
+                                         const ComponentMask &);
 #endif
 }
 
index a01fff9f0292d9f9b6f54c934b8158cbb4368b77..5ff4f0bd56088beb2a2739bc1713d62d180c7a56 100644 (file)
@@ -2402,7 +2402,6 @@ next_cell:
             ExcMessage ("Unmatched faces on periodic boundaries"));
 
     // and call match_periodic_face_pairs that does the actual matching:
-
     return match_periodic_face_pairs(faces1, faces2, direction, offset);
   }
 
@@ -2446,11 +2445,72 @@ next_cell:
             ExcMessage ("Unmatched faces on periodic boundaries"));
 
     // and call match_periodic_face_pairs that does the actual matching:
-
     return match_periodic_face_pairs(faces1, faces2, direction, offset);
   }
 
 
+
+  template<typename DH>
+  std::map<typename DH::face_iterator, typename DH::face_iterator>
+  collect_periodic_face_pairs (const DH                 &dof_handler,
+                               const types::boundary_id b_id,
+                               int                      direction,
+                               const dealii::Tensor<1,DH::space_dimension> &offset)
+  {
+    static const int dim = DH::dimension;
+    static const int space_dim = DH::space_dimension;
+    Assert (0<=direction && direction<space_dim,
+            ExcIndexRange (direction, 0, space_dim));
+
+    Assert(dim == space_dim,
+           ExcNotImplemented());
+
+    // Loop over all cells on the highest level and collect all boundary
+    // faces 2*direction and 2*direction*1:
+
+    std::set<typename DH::face_iterator> faces1;
+    std::set<typename DH::face_iterator> faces2;
+
+    for (typename DH::cell_iterator cell = dof_handler.begin(0);
+         cell != dof_handler.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);
+
+        if (face_1->at_boundary() && face_1->boundary_indicator() == b_id)
+          faces1.insert(face_1);
+
+        if (face_2->at_boundary() && face_2->boundary_indicator() == b_id)
+          faces2.insert(face_2);
+      }
+
+    Assert (faces1.size() == faces2.size(),
+            ExcMessage ("Unmatched faces on periodic boundaries"));
+
+    // and call match_periodic_face_pairs that does the actual matching:
+
+    typedef std::map<typename DH::face_iterator, std::pair<typename DH::face_iterator, std::bitset<3> > > FaceMap;
+    FaceMap matching = match_periodic_face_pairs(faces1, faces2, direction, offset);
+
+    std::map<typename DH::face_iterator, typename DH::face_iterator>
+      return_value;
+
+    for(typename FaceMap::iterator pairing = matching.begin();
+        pairing != matching.end(); ++pairing)
+      {
+        Assert(pairing->second.second == 1,
+               ExcMessage("Found a face match with non standard orientation. "
+                          "This function is only suitable for meshes with cells "
+                          "in default orientation"));
+
+        return_value[pairing->first] = pairing->second.first;
+      }
+
+    return return_value;
+  }
+
+
+
 } /* namespace GridTools */
 
 
index fe5ce669bb94fafa249fd85fc9e624b3a4374f6b..2343745ff905c6215cc8cb282b55081244349d44 100644 (file)
@@ -240,6 +240,14 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS ; deal_II
                                   int,
                                   const Tensor<1,deal_II_space_dimension> &);
 
+      template
+      std::map<X::face_iterator, X::face_iterator>
+      collect_periodic_face_pairs
+                                 (const X &,
+                                  const types::boundary_id,
+                                  int,
+                                  const Tensor<1,deal_II_space_dimension> &);
+
     #endif
 
   \}

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.