From: maier Date: Sat, 1 Dec 2012 21:25:44 +0000 (+0000) Subject: Reintroduce the old make_periodicity_constraints interface as well (because this... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=ef45aeff10f8703221a4b6b91749c3cfd5179522;p=dealii-svn.git Reintroduce the old make_periodicity_constraints interface as well (because this one is _very_ convenient) git-svn-id: https://svn.dealii.org/trunk@27722 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/dofs/dof_tools.h b/deal.II/include/deal.II/dofs/dof_tools.h index 5a9b032536..53ebedd75b 100644 --- a/deal.II/include/deal.II/dofs/dof_tools.h +++ b/deal.II/include/deal.II/dofs/dof_tools.h @@ -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 2*dimension and boundary + * indicator @p b_id and, similarly, a 'right' boundary consisting of all + * face with local face index 2*dimension+1 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 + 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 + 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()); + + //@} diff --git a/deal.II/include/deal.II/grid/grid_tools.h b/deal.II/include/deal.II/grid/grid_tools.h index 8be9ef13e9..6abe799679 100644 --- a/deal.II/include/deal.II/grid/grid_tools.h +++ b/deal.II/include/deal.II/grid/grid_tools.h @@ -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 2*dimension and boundary + * indicator @p b_id and, similarly, a 'right' boundary consisting of all + * face with local face index 2*dimension+1 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 + std::map + 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 diff --git a/deal.II/source/dofs/dof_tools.cc b/deal.II/source/dofs/dof_tools.cc index e6133255f6..874a78c99d 100644 --- a/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/source/dofs/dof_tools.cc @@ -3556,6 +3556,88 @@ namespace DoFTools + template + 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 + 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*> + (&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 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 diff --git a/deal.II/source/dofs/dof_tools.inst.in b/deal.II/source/dofs/dof_tools.inst.in index 7aa2cb456a..0285982028 100644 --- a/deal.II/source/dofs/dof_tools.inst.in +++ b/deal.II/source/dofs/dof_tools.inst.in @@ -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 } diff --git a/deal.II/source/grid/grid_tools.cc b/deal.II/source/grid/grid_tools.cc index a01fff9f02..5ff4f0bd56 100644 --- a/deal.II/source/grid/grid_tools.cc +++ b/deal.II/source/grid/grid_tools.cc @@ -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 + std::map + 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 faces1; + std::set 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 > > FaceMap; + FaceMap matching = match_periodic_face_pairs(faces1, faces2, direction, offset); + + std::map + 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 */ diff --git a/deal.II/source/grid/grid_tools.inst.in b/deal.II/source/grid/grid_tools.inst.in index fe5ce669bb..2343745ff9 100644 --- a/deal.II/source/grid/grid_tools.inst.in +++ b/deal.II/source/grid/grid_tools.inst.in @@ -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 + collect_periodic_face_pairs + (const X &, + const types::boundary_id, + int, + const Tensor<1,deal_II_space_dimension> &); + #endif \}