From ab39c727c50e067b5e73d806b0d3eee41e2a59e0 Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Tue, 29 Mar 2016 19:26:54 +0200 Subject: [PATCH] periodic_neighbor in DoFTools::make_sparsity_pattern --- include/deal.II/dofs/dof_accessor.h | 26 ++++++++++++++ source/dofs/dof_accessor.cc | 37 ++++++++++++++++++++ source/dofs/dof_tools_sparsity.cc | 54 +++++++++++++++++++---------- 3 files changed, 98 insertions(+), 19 deletions(-) diff --git a/include/deal.II/dofs/dof_accessor.h b/include/deal.II/dofs/dof_accessor.h index e25e74412f..8b257ed11e 100644 --- a/include/deal.II/dofs/dof_accessor.h +++ b/include/deal.II/dofs/dof_accessor.h @@ -1264,6 +1264,22 @@ public: TriaIterator > neighbor (const unsigned int) const; + /** + * Return the @p ith periodic neighbor as a DoF cell iterator. This function + * is needed since the neighbor function of the base class returns a cell + * accessor without access to the DoF data. + */ + TriaIterator > + periodic_neighbor (const unsigned int) const; + + /** + * Return the @p ith neighbor or periodic neighbor as a DoF cell iterator. + * This function is needed since the neighbor function of the base class + * returns a cell accessor without access to the DoF data. + */ + TriaIterator > + neighbor_or_periodic_neighbor (const unsigned int) const; + /** * Return the @p ith child as a DoF cell iterator. This function is needed * since the child function of the base class returns a cell accessor @@ -1291,6 +1307,16 @@ public: neighbor_child_on_subface (const unsigned int face_no, const unsigned int subface_no) const; + /** + * Return the result of the @p periodic_neighbor_child_on_subface function + * of the base class, but convert it so that one can also access the DoF + * data (the function in the base class only returns an iterator with access + * to the triangulation data). + */ + TriaIterator > + periodic_neighbor_child_on_subface (const unsigned int face_no, + const unsigned int subface_no) const; + /** * @} */ diff --git a/source/dofs/dof_accessor.cc b/source/dofs/dof_accessor.cc index 588654d0d9..b2751a9663 100644 --- a/source/dofs/dof_accessor.cc +++ b/source/dofs/dof_accessor.cc @@ -84,6 +84,43 @@ DoFCellAccessor::neighbor_child_on_subface (const unsigned i +template +TriaIterator > +DoFCellAccessor::periodic_neighbor_child_on_subface (const unsigned int face, + const unsigned int subface) const +{ + const TriaIterator > q + = CellAccessor::periodic_neighbor_child_on_subface (face, subface); + return TriaIterator > (*q, this->dof_handler); +} + + + +template +inline +TriaIterator > +DoFCellAccessor::periodic_neighbor (const unsigned int face) const +{ + const TriaIterator > q + = CellAccessor::periodic_neighbor (face); + return TriaIterator > (*q, this->dof_handler); +} + + + +template +inline +TriaIterator > +DoFCellAccessor::neighbor_or_periodic_neighbor (const unsigned int face) const +{ + const TriaIterator > q + = CellAccessor::neighbor_or_periodic_neighbor (face); + return TriaIterator > (*q, this->dof_handler); +} + + + + // -------------------------------------------------------------------------- // explicit instantiations #include "dof_accessor.inst" diff --git a/source/dofs/dof_tools_sparsity.cc b/source/dofs/dof_tools_sparsity.cc index 5a98b4e698..7d9db73589 100644 --- a/source/dofs/dof_tools_sparsity.cc +++ b/source/dofs/dof_tools_sparsity.cc @@ -553,13 +553,15 @@ namespace DoFTools ++face) { typename DoFHandlerType::face_iterator cell_face = cell->face(face); - if (! cell->at_boundary(face) ) + const bool periodic_neighbor = cell->has_periodic_neighbor(face); + if (! cell->at_boundary(face) || periodic_neighbor) { - typename DoFHandlerType::level_cell_iterator neighbor = cell->neighbor(face); + typename DoFHandlerType::level_cell_iterator neighbor + = cell->neighbor_or_periodic_neighbor(face); // in 1d, we do not need to worry whether the neighbor // might have children and then loop over those children. - // rather, we may as well go straight to to cell behind + // rather, we may as well go straight to the cell behind // this particular cell's most terminal child if (DoFHandlerType::dimension==1) while (neighbor->has_children()) @@ -571,9 +573,10 @@ namespace DoFTools sub_nr != cell_face->number_of_children(); ++sub_nr) { - const typename DoFHandlerType::level_cell_iterator - sub_neighbor - = cell->neighbor_child_on_subface (face, sub_nr); + const typename DoFHandlerType::level_cell_iterator sub_neighbor + = periodic_neighbor? + cell->periodic_neighbor_child_on_subface (face, sub_nr): + cell->neighbor_child_on_subface (face, sub_nr); const unsigned int n_dofs_on_neighbor = sub_neighbor->get_fe().dofs_per_cell; @@ -598,9 +601,10 @@ namespace DoFTools { // Refinement edges are taken care of by coarser // cells - if (cell->neighbor_is_coarser(face) && - neighbor->subdomain_id() == cell->subdomain_id()) - continue; + if ((!periodic_neighbor && cell->neighbor_is_coarser(face)) || + (periodic_neighbor && cell->periodic_neighbor_is_coarser(face))) + if (neighbor->subdomain_id() == cell->subdomain_id()) + continue; const unsigned int n_dofs_on_neighbor = neighbor->get_fe().dofs_per_cell; @@ -616,7 +620,7 @@ namespace DoFTools // is not locally owned - otherwise, we touch each // face twice and hence put the indices the other way // around - if (!cell->neighbor(face)->active() + if (!cell->neighbor_or_periodic_neighbor(face)->active() || (neighbor->subdomain_id() != cell->subdomain_id())) { @@ -755,7 +759,9 @@ namespace DoFTools const typename DoFHandlerType::face_iterator cell_face = cell->face (face_n); - if (cell->at_boundary (face_n)) + const bool periodic_neighbor = cell->has_periodic_neighbor (face_n); + + if (cell->at_boundary (face_n) && (!periodic_neighbor)) { for (unsigned int i=0; ineighbor (face_n); + neighbor = cell->neighbor_or_periodic_neighbor (face_n); // If the cells are on the same level then only add to // the sparsity pattern if the current cell is 'greater' // in the total ordering. @@ -796,7 +802,8 @@ namespace DoFTools // if the neighbor is coarser than the current cell. if (neighbor->level () != cell->level () && - !cell->neighbor_is_coarser (face_n)) + ((!periodic_neighbor && !cell->neighbor_is_coarser (face_n)) + ||(periodic_neighbor && !cell->periodic_neighbor_is_coarser (face_n)))) continue; // (the neighbor is finer) const unsigned int @@ -810,7 +817,9 @@ namespace DoFTools { const typename DoFHandlerType::level_cell_iterator sub_neighbor - = cell->neighbor_child_on_subface (face_n, sub_nr); + = periodic_neighbor? + cell->periodic_neighbor_child_on_subface (face_n, sub_nr): + cell->neighbor_child_on_subface (face_n, sub_nr); sub_neighbor->get_dof_indices (dofs_on_other_cell); for (unsigned int i=0; iuser_flag_set ()) continue; - if (cell->at_boundary (face) ) + const bool periodic_neighbor = cell->has_periodic_neighbor (face); + + if (cell->at_boundary (face) && (!periodic_neighbor)) { for (unsigned int i=0; iget_fe().dofs_per_cell; ++i) for (unsigned int j=0; jget_fe().dofs_per_cell; ++j) @@ -1019,14 +1030,17 @@ namespace DoFTools else { typename dealii::hp::DoFHandler::level_cell_iterator - neighbor = cell->neighbor(face); + neighbor = cell->neighbor_or_periodic_neighbor(face); // Refinement edges are taken care of by coarser cells - if (cell->neighbor_is_coarser(face)) + if ((!periodic_neighbor && cell->neighbor_is_coarser(face)) || + (periodic_neighbor && cell->periodic_neighbor_is_coarser(face))) continue; const unsigned int - neighbor_face = cell->neighbor_of_neighbor(face); + neighbor_face = periodic_neighbor? + cell->periodic_neighbor_of_periodic_neighbor(face): + cell->neighbor_of_neighbor(face); if (cell_face->has_children()) { @@ -1036,7 +1050,9 @@ namespace DoFTools { const typename dealii::hp::DoFHandler::level_cell_iterator sub_neighbor - = cell->neighbor_child_on_subface (face, sub_nr); + = periodic_neighbor? + cell->periodic_neighbor_child_on_subface (face, sub_nr): + cell->neighbor_child_on_subface (face, sub_nr); dofs_on_other_cell.resize (sub_neighbor->get_fe().dofs_per_cell); sub_neighbor->get_dof_indices (dofs_on_other_cell); -- 2.39.5