From: maier <maier@0785d39b-7218-0410-832d-ea1e28bc413d> Date: Mon, 30 Sep 2013 23:49:41 +0000 (+0000) Subject: Take over revisions 31045 - 31051 from branch_port_the_testsuite X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=eb5a74a80f25edf88f0761949ccbfae8a1757b67;p=dealii-svn.git Take over revisions 31045 - 31051 from branch_port_the_testsuite Cleanup GridTools::collect_periodic_*: - Remove GridTools::collect_periodic_face_pairs - Unify interface for remaining GridTools::collect_periodic_faces - Remove deprecated functions that aren't even released - Bugfix: Instantiate collect_periodic_faces, this resolves a compilation issue with icc git-svn-id: https://svn.dealii.org/trunk@31052 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 16257eeffb..ab041d8c68 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -24,6 +24,14 @@ inconvenience this causes. </p> <ol> + <li> + Removed: GridTools::collect_periodic_face_pairs. This function is superseded + by GridTools::collect_periodic_faces which exports an + std::vector<PeriodicFacepair<...>> instead. + <br> + (Matthias Maier, 2013/09/30) + </li> + <li> Removed: The member function face_to_equivalent_cell_index() in FiniteElementData has been removed. It had been deprecated a while diff --git a/deal.II/include/deal.II/distributed/tria.h b/deal.II/include/deal.II/distributed/tria.h index 5d647258f1..575116d455 100644 --- a/deal.II/include/deal.II/distributed/tria.h +++ b/deal.II/include/deal.II/distributed/tria.h @@ -719,28 +719,6 @@ namespace parallel add_periodicity (const std::vector<GridTools::PeriodicFacePair<cell_iterator> >&); - /** - * Same as the function above, but takes a different argument. - * - * The entries in the std::vector should have the form - * std_cxx1x::tuple<cell1, face_no1, cell2, face_no2>. - * - * The vector can be filled by the function - * GridTools::identify_periodic_face_pairs. - * - * @note This function can only be used if the faces are in - * default orientation. - * - * @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. - */ - void - add_periodicity - (const std::vector<std_cxx1x::tuple<cell_iterator, unsigned int, - cell_iterator, unsigned int> >&); - private: diff --git a/deal.II/include/deal.II/grid/grid_tools.h b/deal.II/include/deal.II/grid/grid_tools.h index 29900ec738..ac6080e565 100644 --- a/deal.II/include/deal.II/grid/grid_tools.h +++ b/deal.II/include/deal.II/grid/grid_tools.h @@ -1061,7 +1061,7 @@ namespace GridTools bool orthogonal_equality (const FaceIterator &face1, const FaceIterator &face2, - const int direction, + const int direction, const dealii::Tensor<1,FaceIterator::AccessorType::space_dimension> &offset); @@ -1094,7 +1094,7 @@ namespace GridTools * parallel::distributed::Triangulation::add_periodicity to enforce * periodicity algebraically. * - * @author Daniel Arndt, 2013 + * @author Daniel Arndt, Matthias Maier, 2013 */ template<typename DH> std::vector<PeriodicFacePair<typename DH::cell_iterator> > @@ -1102,9 +1102,10 @@ namespace GridTools (const DH &dof_handler, const types::boundary_id b_id1, const types::boundary_id b_id2, - const unsigned int direction, + const int direction, 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". @@ -1123,102 +1124,16 @@ namespace GridTools * meshes with cells not in @ref GlossFaceOrientation * "standard orientation". * - * @author Daniel Arndt, 2013 + * @author Daniel Arndt, Matthias Maier, 2013 */ template<typename DH> std::vector<PeriodicFacePair<typename DH::cell_iterator> > collect_periodic_faces (const DH &dof_handler, const types::boundary_id b_id, - const unsigned int direction, + const int direction, const dealii::Tensor<1,DH::space_dimension> &offset); - /** - * This function does the same as collect_periodic_faces but returns a - * different data type. - * - * @author Matthias Maier, 2012 - * - * @note The returned data type is not compatible with - * DoFTools::make_periodicity_constraints - * - * @deprecated - */ - template<typename DH> - std::map<typename DH::face_iterator, std::pair<typename DH::face_iterator, std::bitset<3> > > - collect_periodic_face_pairs - (const DH &container, - const types::boundary_id b_id1, - const types::boundary_id b_id2, - int direction, - const dealii::Tensor<1,DH::space_dimension> &offset) DEAL_II_DEPRECATED; - - /** - * This compatibility version of collect_periodic_face_pairs only works - * on grids with cells in @ref GlossFaceOrientation "standard orientation". - * - * @author Matthias Maier, 2012 - * - * @note The returned data type is not compatible with - * DoFTools::make_periodicity_constraints - * - * @deprecated - */ - 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) DEAL_II_DEPRECATED; - - /** - * This version loops over all faces from @p begin to @p end - * instead of accepting a DoFHandler or a Triangulation. - * - * @author Matthias Maier, 2012 - * - * @note This function cannot produce the return as the other - * collect_periodic_faces functions. - * - * @deprecated - */ - template<typename FaceIterator> - std::map<FaceIterator, std::pair<FaceIterator, std::bitset<3> > > - collect_periodic_face_pairs - (const FaceIterator &begin, - const typename identity<FaceIterator>::type &end, - const types::boundary_id b_id1, - const types::boundary_id b_id2, - const int direction, - const dealii::Tensor<1,FaceIterator::AccessorType::space_dimension> &offset) - DEAL_II_DEPRECATED; - - /** - * This function does the same as collect_periodic_faces but returns a - * different data type. - * - * @author Daniel Arndt, 2013 - * - * @note The returned data type is not compatible with - * DoFTools::make_periodicity_constraints, but with a version of - * parallel::distributed::Triangulation::add_periodicity - * - * @note Use collect_periodic_faces instead. - * - * @deprecated - */ - template<typename DH> - void - identify_periodic_face_pairs - (const DH &dof_handler, - const types::boundary_id b_id1, - const types::boundary_id b_id2, - const unsigned int direction, - std::vector<std_cxx1x::tuple<typename DH::cell_iterator, unsigned int, - typename DH::cell_iterator, unsigned int> > - &periodicity_vector) DEAL_II_DEPRECATED; - /** * Exception diff --git a/deal.II/source/distributed/tria.cc b/deal.II/source/distributed/tria.cc index 5b0caa979f..fae79cacfa 100644 --- a/deal.II/source/distributed/tria.cc +++ b/deal.II/source/distributed/tria.cc @@ -3529,41 +3529,6 @@ namespace parallel #endif } - - template <int dim, int spacedim> - void - Triangulation<dim,spacedim>::add_periodicity - (const std::vector<std_cxx1x::tuple<cell_iterator, unsigned int, - cell_iterator, unsigned int> >& - periodicity_vector) - { -#if DEAL_II_P4EST_VERSION_GTE(0,3,4,1) - typedef std::vector<std_cxx1x::tuple <cell_iterator, unsigned int, - cell_iterator, unsigned int> > - FaceVector; - - typename FaceVector::const_iterator it, end_periodic; - it = periodicity_vector.begin(); - end_periodic = periodicity_vector.end(); - - std::vector<GridTools::PeriodicFacePair<cell_iterator> > periodic_faces; - - for(; it!=end_periodic; ++it) - { - const cell_iterator cell1 = std_cxx1x::get<0> (*it); - const cell_iterator cell2 = std_cxx1x::get<2> (*it); - const unsigned int face_idx1 = std_cxx1x::get<1> (*it); - const unsigned int face_idx2 = std_cxx1x::get<3> (*it); - const GridTools::PeriodicFacePair<cell_iterator> matched_face - = {{cell1, cell2},{face_idx1, face_idx2}, 1}; - periodic_faces.push_back(matched_face); - } - - add_periodicity(periodic_faces); -#else - Assert(false, ExcMessage ("Need p4est version >= 0.3.4.1!")); -#endif - } template <int dim, int spacedim> diff --git a/deal.II/source/grid/grid_tools.cc b/deal.II/source/grid/grid_tools.cc index 4f6981479d..df146a150b 100644 --- a/deal.II/source/grid/grid_tools.cc +++ b/deal.II/source/grid/grid_tools.cc @@ -2199,7 +2199,7 @@ next_cell: template<int spacedim> inline bool orthogonal_equality (const dealii::Point<spacedim> &point1, const dealii::Point<spacedim> &point2, - const int direction, + const int direction, const dealii::Tensor<1,spacedim> &offset) { Assert (0<=direction && direction<spacedim, @@ -2343,7 +2343,7 @@ next_cell: inline bool orthogonal_equality (const FaceIterator &face1, const FaceIterator &face2, - const int direction, + const int direction, const dealii::Tensor<1,FaceIterator::AccessorType::space_dimension> &offset) { // Call the function above with a dummy orientation array @@ -2359,13 +2359,10 @@ next_cell: template<typename CellIterator> std::vector<PeriodicFacePair<CellIterator> > match_periodic_face_pairs - (std::set<std::pair<CellIterator, unsigned int> > - &pairs1, - std::set<std::pair<typename identity<CellIterator>::type, unsigned int> > - &pairs2, - int direction, - const dealii::Tensor<1,CellIterator::AccessorType::space_dimension> - &offset) + (std::set<std::pair<CellIterator, unsigned int> > &pairs1, + std::set<std::pair<typename identity<CellIterator>::type, unsigned int> > &pairs2, + const int direction, + const dealii::Tensor<1,CellIterator::AccessorType::space_dimension> &offset) { static const int space_dim = CellIterator::AccessorType::space_dimension; Assert (0<=direction && direction<space_dim, @@ -2411,61 +2408,6 @@ next_cell: return matched_faces; } - /* Deprecated version of the function above with different return value. - * It is used the deprecated collect_periodic_face_pairs. - */ - template<typename FaceIterator> - std::map<FaceIterator, std::pair<FaceIterator, std::bitset<3> > > - match_periodic_face_pairs - (std::set<FaceIterator> &faces1, /* not const! */ - std::set<typename identity<FaceIterator>::type> &faces2, /* not const! */ - int direction, - const dealii::Tensor<1,FaceIterator::AccessorType::space_dimension> &offset) - DEAL_II_DEPRECATED; - - template<typename FaceIterator> - std::map<FaceIterator, std::pair<FaceIterator, std::bitset<3> > > - match_periodic_face_pairs - (std::set<FaceIterator> &faces1, /* not const! */ - std::set<typename identity<FaceIterator>::type> &faces2, /* not const! */ - int direction, - const dealii::Tensor<1,FaceIterator::AccessorType::space_dimension> &offset) - { - static const int space_dim = FaceIterator::AccessorType::space_dimension; - Assert (0<=direction && direction<space_dim, - ExcIndexRange (direction, 0, space_dim)); - - Assert (faces1.size() == faces2.size(), - ExcMessage ("Unmatched faces on periodic boundaries")); - - typedef std::pair<FaceIterator, std::bitset<3> > ResultPair; - std::map<FaceIterator, ResultPair> matched_faces; - - // Match with a complexity of O(n^2). This could be improved... - std::bitset<3> orientation; - typedef typename std::set<FaceIterator>::const_iterator SetIterator; - for (SetIterator it1 = faces1.begin(); it1 != faces1.end(); ++it1) - { - for (SetIterator it2 = faces2.begin(); it2 != faces2.end(); ++it2) - { - if (GridTools::orthogonal_equality(orientation, *it1, *it2, - direction, offset)) - { - // We have a match, so insert the matching pairs and - // remove the matched cell in faces2 to speed up the - // matching: - matched_faces[*it1] = std::make_pair(*it2, orientation); - faces2.erase(it2); - break; - } - } - } - - AssertThrow (matched_faces.size() == faces1.size() && faces2.size() == 0, - ExcMessage ("Unmatched faces on periodic boundaries")); - - return matched_faces; - } template<typename DH> @@ -2474,11 +2416,11 @@ next_cell: (const DH &dof_handler, const types::boundary_id b_id1, const types::boundary_id b_id2, - const unsigned int direction, + const int direction, const dealii::Tensor<1,DH::space_dimension> &offset) { - static const unsigned int dim = DH::dimension; - static const unsigned int space_dim = DH::space_dimension; + 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)); @@ -2517,40 +2459,42 @@ next_cell: return match_periodic_face_pairs(pairs1, pairs2, direction, offset); } + + template<typename DH> std::vector<PeriodicFacePair<typename DH::cell_iterator> > collect_periodic_faces (const DH &dof_handler, const types::boundary_id b_id, - const unsigned int direction, + const int direction, const dealii::Tensor<1,DH::space_dimension> &offset) { - static const unsigned int dim = DH::dimension; - static const unsigned int space_dim = DH::space_dimension; + 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<std::pair<typename DH::cell_iterator, unsigned int> > pairs1; std::set<std::pair<typename DH::cell_iterator, unsigned int> > pairs2; - + 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) { const std::pair<typename DH::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 @@ -2558,10 +2502,10 @@ next_cell: pairs2.insert(pair2); } } - + Assert (pairs1.size() == pairs2.size(), ExcMessage ("Unmatched faces on periodic boundaries")); - + // and call match_periodic_face_pairs that does the actual matching: typedef std::vector<PeriodicFacePair<typename DH::cell_iterator> > @@ -2570,137 +2514,20 @@ next_cell: FaceVector matching = match_periodic_face_pairs(pairs1, pairs2, direction, offset); +#ifdef DEBUG for (typename FaceVector::iterator pairing = matching.begin(); pairing != matching.end(); ++pairing) { Assert(pairing->orientation == 1, - ExcMessage("Found a face match with non standard orientation. " - "This function is only suitable for meshes with cells " - "in default orientation")); + ExcMessage("Found a face match with non standard orientation. " + "This function is only suitable for meshes with cells " + "in default orientation")); } - - return matching; - } - - template<typename FaceIterator> - std::map<FaceIterator, std::pair<FaceIterator, std::bitset<3> > > - collect_periodic_face_pairs (const FaceIterator &begin, - const typename identity<FaceIterator>::type &end, - const types::boundary_id b_id1, - const types::boundary_id b_id2, - int direction, - const dealii::Tensor<1,FaceIterator::AccessorType::space_dimension> &offset) - { - static const int space_dim = FaceIterator::AccessorType::space_dimension; - Assert (0<=direction && direction<space_dim, - ExcIndexRange (direction, 0, space_dim)); - - // Collect all faces belonging to b_id1 and b_id2: - - std::set<FaceIterator> faces1; - std::set<FaceIterator> faces2; - - for (FaceIterator face = begin; face!= end; ++face) - { - if (face->at_boundary() && face->boundary_indicator() == b_id1) - faces1.insert(face); - - if (face->at_boundary() && face->boundary_indicator() == b_id2) - faces2.insert(face); - } - - Assert (faces1.size() == faces2.size(), - 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, std::pair<typename DH::face_iterator, std::bitset<3> > > - collect_periodic_face_pairs (const DH &dof_handler, - const types::boundary_id b_id1, - const types::boundary_id b_id2, - int direction, - const dealii::Tensor<1,DH::space_dimension> &offset) - { - typedef std::vector<PeriodicFacePair<typename DH::cell_iterator> > FaceVector; - - const FaceVector face_vector - = collect_periodic_faces (dof_handler, b_id1, b_id2, direction, offset); - - std::map<typename DH::face_iterator, - std::pair<typename DH::face_iterator, std::bitset<3> > > - return_value; - for(typename FaceVector::const_iterator it = face_vector.begin(); - it != face_vector.end(); ++it) - { - const typename DH::face_iterator face1 = it->cell[0]->face(it->face_idx[0]); - const typename DH::face_iterator face2 = it->cell[1]->face(it->face_idx[1]); - return_value[face1] = std::make_pair(face2, it->orientation); - } - - return return_value; - } +#endif - - 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) - { - typedef std::vector<PeriodicFacePair<typename DH::cell_iterator> > FaceVector; - - const FaceVector face_vector - = collect_periodic_faces (dof_handler, b_id, direction, offset); - - std::map<typename DH::face_iterator, typename DH::face_iterator> return_value; - for(typename FaceVector::const_iterator it = face_vector.begin(); - it != face_vector.end(); ++it) - { - const typename DH::face_iterator face1 = it->cell[0]->face(it->face_idx[0]); - const typename DH::face_iterator face2 = it->cell[1]->face(it->face_idx[1]); - return_value[face1] = face2; - } - - return return_value; + return matching; } - - template<typename DH> - void - identify_periodic_face_pairs - (const DH &dof_handler, - const types::boundary_id b_id1, - const types::boundary_id b_id2, - const unsigned int direction, - std::vector<std_cxx1x::tuple<typename DH::cell_iterator, unsigned int, - typename DH::cell_iterator, unsigned int> > - &periodicity_vector) - { - typedef std::vector<PeriodicFacePair<typename DH::cell_iterator> > - FaceVector; - const FaceVector periodic_faces - = collect_periodic_faces(dof_handler, b_id1, b_id2, direction, - dealii::Tensor<1,DH::space_dimension> ()); - - typename FaceVector::const_iterator it, end_faces; - it = periodic_faces.begin(); - end_faces = periodic_faces.end(); - - for(; it!=end_faces; ++it) - { - - const std_cxx1x::tuple<typename DH::cell_iterator, unsigned int, - typename DH::cell_iterator, unsigned int> - periodic_tuple (it->cell[0], it->face_idx[0], - it->cell[1], it->face_idx[1]); - - periodicity_vector.push_back(periodic_tuple); - } - } } /* namespace GridTools */ diff --git a/deal.II/source/grid/grid_tools.inst.in b/deal.II/source/grid/grid_tools.inst.in index 4f16d66392..6f4d0f0d6e 100644 --- a/deal.II/source/grid/grid_tools.inst.in +++ b/deal.II/source/grid/grid_tools.inst.in @@ -253,92 +253,24 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS ; deal_II template std::vector<PeriodicFacePair<X::cell_iterator> > - collect_periodic_faces + collect_periodic_faces<X> (const X &, const types::boundary_id, const types::boundary_id, - unsigned int, + const int, const Tensor<1,X::space_dimension> &); - - template<typename DH> - std::vector<PeriodicFacePair<X::cell_iterator> > - collect_periodic_faces - (const X &, - const types::boundary_id, - unsigned int, - const Tensor<1,X::space_dimension> &); - - template - std::map<X::face_iterator, std::pair<X::face_iterator, std::bitset<3> > > - collect_periodic_face_pairs - (const X::face_iterator &, - const X::face_iterator &, - const types::boundary_id, - const types::boundary_id, - int, - const Tensor<1,X::space_dimension> &); - template - std::map<X::active_face_iterator, std::pair<X::active_face_iterator, std::bitset<3> > > - collect_periodic_face_pairs - (const X::active_face_iterator &, - const X::active_face_iterator &, - const types::boundary_id, - const types::boundary_id, - int, - const Tensor<1,X::active_face_iterator::AccessorType::space_dimension> &); template - std::map<X::face_iterator, std::pair<X::face_iterator, std::bitset<3> > > - collect_periodic_face_pairs + std::vector<PeriodicFacePair<X::cell_iterator> > + collect_periodic_faces<X> (const X &, const types::boundary_id, - const types::boundary_id, - int, + const int, const Tensor<1,X::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> &); - - template - void - identify_periodic_face_pairs - (const X &, - const types::boundary_id, - const types::boundary_id, - const unsigned int, - std::vector<std_cxx1x::tuple<X::cell_iterator, unsigned int, - X::cell_iterator, unsigned int> > &); - #endif \} #endif } -for (deal_II_dimension : DIMENSIONS ; deal_II_space_dimension : SPACE_DIMENSIONS) -{ -#if deal_II_dimension <= deal_II_space_dimension - #if deal_II_dimension >= 2 - namespace GridTools \{ - template - void - identify_periodic_face_pairs - (const parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension> &, - const types::boundary_id, - const types::boundary_id, - const unsigned int, - std::vector<std_cxx1x::tuple<parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension>::cell_iterator, - unsigned int, - parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension>::cell_iterator, - unsigned int> > &); - \} - #endif -#endif -} - -