From: Peter Munch Date: Tue, 1 Aug 2023 17:09:53 +0000 (+0200) Subject: Introduce Triangulation::as_dof_handler_level_iterator() X-Git-Tag: relicensing~615^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=cbbbbba621150e001429bc0a35c9619e872be7c4;p=dealii.git Introduce Triangulation::as_dof_handler_level_iterator() --- diff --git a/doc/news/changes/minor/20230801Munch_2 b/doc/news/changes/minor/20230801Munch_2 new file mode 100644 index 0000000000..8025fb633c --- /dev/null +++ b/doc/news/changes/minor/20230801Munch_2 @@ -0,0 +1,4 @@ +New: The new function Triangulation::as_dof_handler_level_iterator() +allows to create level iterators based on other cell iterators. +
+(Peter Munch, 2023/08/01) diff --git a/include/deal.II/grid/tria_accessor.h b/include/deal.II/grid/tria_accessor.h index 65e1a02047..8d9265b90c 100644 --- a/include/deal.II/grid/tria_accessor.h +++ b/include/deal.II/grid/tria_accessor.h @@ -3184,6 +3184,19 @@ public: TriaActiveIterator> as_dof_handler_iterator(const DoFHandler &dof_handler) const; + /** + * A function similar to as_dof_handler_iterator(). It converts + * a Triangulation active/level cell iterator to a + * DoFHandler active level cell iterator, or a DoFHandler level active/cell + * iterator to a level cell iterator of another DoFHandler. The + * @p iterator must be associated with the triangulation of the + * @p dof_handler. + */ + TriaIterator> + as_dof_handler_level_iterator( + const DoFHandler &dof_handler) const; + + /** * @} */ diff --git a/source/grid/tria_accessor.cc b/source/grid/tria_accessor.cc index 98602affe2..3267aae53e 100644 --- a/source/grid/tria_accessor.cc +++ b/source/grid/tria_accessor.cc @@ -2071,6 +2071,25 @@ CellAccessor::as_dof_handler_iterator( } + +template +TriaIterator> +CellAccessor::as_dof_handler_level_iterator( + const DoFHandler &dof_handler) const +{ + Assert(&this->get_triangulation() == &dof_handler.get_triangulation(), + ExcMessage("The triangulation associated with the iterator does not " + "match that of the DoFHandler.")); + + return typename DoFHandler::level_cell_iterator( + &dof_handler.get_triangulation(), + this->level(), + this->index(), + &dof_handler); +} + + + // For codim>0 we proceed as follows: // 1) project point onto manifold and // 2) transform to the unit cell with a Q1 mapping