From: Wolfgang Bangerth Date: Thu, 11 Nov 2010 23:35:45 +0000 (+0000) Subject: Implement face iterators also in 1d -- i.e. you can write cell->face()->at_boundary... X-Git-Tag: v8.0.0~4958 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=8fbad252da04e5dcf7efac431f3c43d37192d64e;p=dealii.git Implement face iterators also in 1d -- i.e. you can write cell->face()->at_boundary() and similar queries also for 1d triangulations. git-svn-id: https://svn.dealii.org/trunk@22699 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/dofs/dof_accessor.h b/deal.II/include/deal.II/dofs/dof_accessor.h index 1ab4b54a07..d443ccf68c 100644 --- a/deal.II/include/deal.II/dofs/dof_accessor.h +++ b/deal.II/include/deal.II/dofs/dof_accessor.h @@ -107,7 +107,7 @@ namespace internal * See the full documentation for * more information. */ - typedef CellAccessor BaseClass; + typedef dealii::CellAccessor BaseClass; }; } } diff --git a/deal.II/include/deal.II/grid/tria_accessor.h b/deal.II/include/deal.II/grid/tria_accessor.h index 020ad28511..50c5d5dcc7 100644 --- a/deal.II/include/deal.II/grid/tria_accessor.h +++ b/deal.II/include/deal.II/grid/tria_accessor.h @@ -1602,28 +1602,98 @@ class TriaAccessor<0, dim, spacedim> * @author Wolfgang Bangerth, 2010 */ template -class TriaAccessor<0, 1, spacedim> : public TriaAccessorBase<0,1, spacedim> +class TriaAccessor<0, 1, spacedim> { public: /** - * Propagate typedef from - * base class to this class. + * Dimension of the space the + * object represented by this + * accessor lives in. For + * example, if this accessor + * represents a quad that is + * part of a two-dimensional + * surface in four-dimensional + * space, then this value is + * four. */ - typedef typename TriaAccessorBase<0,1,spacedim>::AccessorData AccessorData; + static const unsigned int space_dimension = spacedim; + + /** + * Dimensionality of the object + * that the thing represented by + * this accessopr is part of. For + * example, if this accessor + * represents a line that is part + * of a hexahedron, then this + * value will be three. + */ + static const unsigned int dimension = 1; + + /** + * Dimensionality of the current + * object represented by this + * accessor. For example, if it + * is line (irrespective of + * whether it is part of a quad + * or hex, and what dimension we + * are in), then this value + * equals 1. + */ + static const unsigned int structure_dimension = 0; + + /** + * Pointer to internal data. + */ + typedef void AccessorData; + + /** + * Whether the vertex represented + * here is at the left end of the + * domain, the right end, or in + * the interior. + */ + enum VertexKind + { + left_vertex, + interior_vertex, + right_vertex + }; /** * Constructor. Should never be * called and thus produces an * error. */ - TriaAccessor (const Triangulation<1,spacedim> *parent = 0, - const int level = -1, - const int index = -1, - const AccessorData *local_data = 0) - : - TriaAccessorBase<0,1, spacedim> (parent, level, index, local_data) + TriaAccessor (const Triangulation<1,spacedim> * tria, + const VertexKind vertex_kind, + const unsigned int vertex_index) + : + tria (tria), + vertex_kind (vertex_kind), + global_vertex_index (vertex_index) {} + /** + * Constructor. This constructor + * exists in order to maintain + * interface compatibility with + * the other accessor + * classes. However, it doesn't + * do anything useful here and so + * may not actually be called. + */ + TriaAccessor (const Triangulation<1,spacedim> * = 0, + const int = 0, + const int = 0, + const AccessorData * = 0) + : + tria (0), + vertex_kind (interior_vertex), + global_vertex_index (numbers::invalid_unsigned_int) + { + Assert (false, ExcInternalError()); + } + /** * Constructor. Should never be * called and thus produces an @@ -1631,6 +1701,10 @@ class TriaAccessor<0, 1, spacedim> : public TriaAccessorBase<0,1, spacedim> */ template TriaAccessor (const TriaAccessor &) + : + tria (0), + vertex_kind (interior_vertex), + global_vertex_index (numbers::invalid_unsigned_int) { Assert(false, ExcImpossibleInDim(0)); } @@ -1642,10 +1716,48 @@ class TriaAccessor<0, 1, spacedim> : public TriaAccessorBase<0,1, spacedim> */ template TriaAccessor (const InvalidAccessor &) + : + tria (0), + vertex_kind (interior_vertex), + global_vertex_index (numbers::invalid_unsigned_int) { Assert(false, ExcImpossibleInDim(0)); } + /** + * Return the state of the + * iterator. Since an iterator + * to points can not be + * incremented or decremented, + * its state remains constant, + * and in particular equal to + * IteratorState::valid. + */ + static IteratorState::IteratorStates state () + { + return IteratorState::valid; + } + + /** + * Level of this object. Vertices + * have no level, so this + * function always returns zero. + */ + static int level () + { + return 0; + } + + /** + * Index of this object. Vertices + * have no index, so this + * function always returns -1. + */ + static int index () + { + return -1; + } + /** * @name Advancement of iterators */ @@ -1716,10 +1828,10 @@ class TriaAccessor<0, 1, spacedim> : public TriaAccessorBase<0,1, spacedim> * @p DoFAccessor::vertex_dof_index * functions. */ - unsigned int vertex_index (const unsigned int i) const + unsigned int vertex_index (const unsigned int i = 0) const { - Assert(false, ExcNotImplemented()); - return numbers::invalid_unsigned_int; + Assert(i==0, ExcIndexRange(i, 0, 1)); + return global_vertex_index; } /** @@ -1730,11 +1842,10 @@ class TriaAccessor<0, 1, spacedim> : public TriaAccessorBase<0,1, spacedim> * object refers. Otherwise, it * throws an exception. */ - Point & vertex (const unsigned int i) const + Point & vertex (const unsigned int i = 0) const { - Assert(false, ExcNotImplemented()); - static const Point p; - return p; + Assert(i==0, ExcIndexRange(i, 0, 1)); + return const_cast &> (this->tria->vertices[global_vertex_index]); } @@ -1743,11 +1854,51 @@ class TriaAccessor<0, 1, spacedim> : public TriaAccessorBase<0,1, spacedim> * bounding this object. Will * point to an invalid object. */ -// typename internal::Triangulation::Iterators<1,spacedim>::line_iterator -// line (const unsigned int i) const -// { -// return TriaIterator >(); -// } + typename internal::Triangulation::Iterators<1,spacedim>::line_iterator + line (const unsigned int i) const + { + return typename internal::Triangulation::Iterators<1,spacedim>::line_iterator(); + } + + /** + * Return whether this point is + * at the boundary of the + * one-dimensional triangulation + * we deal with here. + */ + bool at_boundary () const + { + return vertex_kind != interior_vertex; + } + + + /** + * Boundary indicator of this + * object. The convention for one + * dimensional triangulations is + * that left end vertices have + * boundary indicator zero, and + * right end vertices have + * boundary indicator one. + * + * If the return value is the special + * value 255, then this object is in the + * interior of the domain. + */ + unsigned char boundary_indicator () const + { + switch (vertex_kind) + { + case left_vertex: + return 0; + case right_vertex: + return 1; + default: + return 255; + } + } + + /** @@ -1770,8 +1921,12 @@ class TriaAccessor<0, 1, spacedim> : public TriaAccessorBase<0,1, spacedim> * Pointer to the @p ith quad * bounding this object. */ -// typename internal::Triangulation::Iterators<1,spacedim>::quad_iterator -// quad (const unsigned int i) const; + typename internal::Triangulation::Iterators<1,spacedim>::quad_iterator + quad (const unsigned int i) const + { + return typename internal::Triangulation::Iterators<1,spacedim>::quad_iterator(); + } + /** * Quad index of the @p ith @@ -1925,15 +2080,6 @@ class TriaAccessor<0, 1, spacedim> : public TriaAccessorBase<0,1, spacedim> { return -1; } - /** - * @brief Do nothing and throw and error - */ - inline - unsigned char boundary_indicator () const - { - Assert(false, ExcImpossibleInDim(0)); - return 0; - } /** * @brief Do nothing and throw and error @@ -1942,7 +2088,9 @@ class TriaAccessor<0, 1, spacedim> : public TriaAccessorBase<0,1, spacedim> void set_boundary_indicator (const unsigned char) const { - Assert(false, ExcImpossibleInDim(0)); + Assert(false, + ExcMessage ("The boundary indicator of vertices is determined by their " + "location within a triangulation and can not be set explicitly")); } /** @@ -1950,18 +2098,24 @@ class TriaAccessor<0, 1, spacedim> : public TriaAccessorBase<0,1, spacedim> */ void set_all_boundary_indicators (const unsigned char) const { - Assert(false, ExcImpossibleInDim(0)); + Assert(false, + ExcMessage ("The boundary indicator of vertices is determined by their " + "location within a triangulation and can not be set explicitly")); } /** - * @brief Do nothing and throw and error + * Return whether the vertex + * pointed to here is used. */ bool used () const { - Assert(false, ExcImpossibleInDim(0)); - return false; + return tria->vertex_used(global_vertex_index); } + protected: + const Triangulation<1,spacedim> * tria; + const VertexKind vertex_kind; + const unsigned int global_vertex_index; }; diff --git a/deal.II/include/deal.II/grid/tria_accessor.templates.h b/deal.II/include/deal.II/grid/tria_accessor.templates.h index 3665478044..d8cf88b9eb 100644 --- a/deal.II/include/deal.II/grid/tria_accessor.templates.h +++ b/deal.II/include/deal.II/grid/tria_accessor.templates.h @@ -1313,7 +1313,7 @@ TriaAccessor:: parent_index () const { Assert (this->present_level > 0, TriaAccessorExceptions::ExcCellHasNoParent ()); - + // the parent of two consecutive cells // is stored only once, since it is // the same @@ -1987,28 +1987,61 @@ CellAccessor (const CellAccessor &cell_accessor) -template -inline -TriaIterator > -CellAccessor::face (const unsigned int i) const +namespace internal { - switch (dim) + namespace CellAccessor + { + template + inline + dealii::TriaIterator > + get_face (const dealii::CellAccessor<1,spacedim> &cell, + const unsigned int i) { - case 1: - { - Assert (false, ExcImpossibleInDim(1)); - return TriaIterator >(); - } + dealii::TriaAccessor<0, 1, spacedim> + a (&cell.get_triangulation(), + ((i == 0) && cell.at_boundary(0) + ? + dealii::TriaAccessor<0, 1, spacedim>::left_vertex + : + ((i == 1) && cell.at_boundary(1) + ? + dealii::TriaAccessor<0, 1, spacedim>::right_vertex + : + dealii::TriaAccessor<0, 1, spacedim>::interior_vertex)), + cell.vertex_index(i)); + return dealii::TriaIterator >(a); + } - case 2: - return this->line(i); - case 3: - return this->quad(i); + template + inline + dealii::TriaIterator > + get_face (const dealii::CellAccessor<2,spacedim> &cell, + const unsigned int i) + { + return cell.line(i); + } - default: - return TriaIterator >(); + + template + inline + dealii::TriaIterator > + get_face (const dealii::CellAccessor<3,spacedim> &cell, + const unsigned int i) + { + return cell.quad(i); } + } +} + + + +template +inline +TriaIterator > +CellAccessor::face (const unsigned int i) const +{ + return internal::CellAccessor::get_face (*this, i); } diff --git a/tests/deal.II/vertex_as_face_01.cc b/tests/deal.II/vertex_as_face_01.cc new file mode 100644 index 0000000000..7cfc4dff87 --- /dev/null +++ b/tests/deal.II/vertex_as_face_01.cc @@ -0,0 +1,61 @@ +//---------------------------- vertex_as_face_01.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2010 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- vertex_as_face_01.cc --------------------------- + +// verify that we can do things like cell->face() in 1d as well. here: +// test vertex location + + +#include "../tests.h" +#include +#include +#include +#include + +#include + + +template +void test () +{ + Triangulation<1,spacedim> tria; + GridGenerator::hyper_cube (tria); + + deallog << "Coarse mesh:" << std::endl; + deallog << "Left vertex=" << tria.begin_active()->face(0)->vertex(0) << std::endl; + deallog << "Right vertex=" << tria.begin_active()->face(1)->vertex(0) << std::endl; + + tria.refine_global (2); + + for (typename Triangulation<1,spacedim>::active_cell_iterator + cell = tria.begin_active(); + cell != tria.end(); ++cell) + { + deallog << "Cell: " << cell << std::endl; + deallog << "Left vertex=" << cell->face(0)->vertex(0) << std::endl; + deallog << "Right vertex=" << cell->face(1)->vertex(0) << std::endl; + } +} + + + +int main () +{ + std::ofstream logfile("vertex_as_face_01/output"); + deallog.attach(logfile); + deallog.depth_console(0); + + test<1> (); + test<2> (); + + return 0; +} diff --git a/tests/deal.II/vertex_as_face_01/cmp/generic b/tests/deal.II/vertex_as_face_01/cmp/generic new file mode 100644 index 0000000000..ba076b7506 --- /dev/null +++ b/tests/deal.II/vertex_as_face_01/cmp/generic @@ -0,0 +1,31 @@ + +DEAL::Coarse mesh: +DEAL::Left vertex=0.00000 +DEAL::Right vertex=1.00000 +DEAL::Cell: 2.0 +DEAL::Left vertex=0.00000 +DEAL::Right vertex=0.250000 +DEAL::Cell: 2.1 +DEAL::Left vertex=0.250000 +DEAL::Right vertex=0.500000 +DEAL::Cell: 2.2 +DEAL::Left vertex=0.500000 +DEAL::Right vertex=0.750000 +DEAL::Cell: 2.3 +DEAL::Left vertex=0.750000 +DEAL::Right vertex=1.00000 +DEAL::Coarse mesh: +DEAL::Left vertex=0.00000 0.00000 +DEAL::Right vertex=1.00000 0.00000 +DEAL::Cell: 2.0 +DEAL::Left vertex=0.00000 0.00000 +DEAL::Right vertex=0.250000 0.00000 +DEAL::Cell: 2.1 +DEAL::Left vertex=0.250000 0.00000 +DEAL::Right vertex=0.500000 0.00000 +DEAL::Cell: 2.2 +DEAL::Left vertex=0.500000 0.00000 +DEAL::Right vertex=0.750000 0.00000 +DEAL::Cell: 2.3 +DEAL::Left vertex=0.750000 0.00000 +DEAL::Right vertex=1.00000 0.00000 diff --git a/tests/deal.II/vertex_as_face_02.cc b/tests/deal.II/vertex_as_face_02.cc new file mode 100644 index 0000000000..aa0315d688 --- /dev/null +++ b/tests/deal.II/vertex_as_face_02.cc @@ -0,0 +1,61 @@ +//---------------------------- vertex_as_face_02.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2010 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- vertex_as_face_02.cc --------------------------- + +// verify that we can do things like cell->face() in 1d as well. here: +// test vertex_index + + +#include "../tests.h" +#include +#include +#include +#include + +#include + + +template +void test () +{ + Triangulation<1,spacedim> tria; + GridGenerator::hyper_cube (tria); + + deallog << "Coarse mesh:" << std::endl; + deallog << "Left vertex=" << tria.begin_active()->face(0)->vertex_index(0) << std::endl; + deallog << "Right vertex=" << tria.begin_active()->face(1)->vertex_index(0) << std::endl; + + tria.refine_global (2); + + for (typename Triangulation<1,spacedim>::active_cell_iterator + cell = tria.begin_active(); + cell != tria.end(); ++cell) + { + deallog << "Cell: " << cell << std::endl; + deallog << "Left vertex=" << cell->face(0)->vertex_index(0) << std::endl; + deallog << "Right vertex=" << cell->face(1)->vertex_index(0) << std::endl; + } +} + + + +int main () +{ + std::ofstream logfile("vertex_as_face_02/output"); + deallog.attach(logfile); + deallog.depth_console(0); + + test<1> (); + test<2> (); + + return 0; +} diff --git a/tests/deal.II/vertex_as_face_02/cmp/generic b/tests/deal.II/vertex_as_face_02/cmp/generic new file mode 100644 index 0000000000..e7b18abd8b --- /dev/null +++ b/tests/deal.II/vertex_as_face_02/cmp/generic @@ -0,0 +1,31 @@ + +DEAL::Coarse mesh: +DEAL::Left vertex=0 +DEAL::Right vertex=1 +DEAL::Cell: 2.0 +DEAL::Left vertex=0 +DEAL::Right vertex=3 +DEAL::Cell: 2.1 +DEAL::Left vertex=3 +DEAL::Right vertex=2 +DEAL::Cell: 2.2 +DEAL::Left vertex=2 +DEAL::Right vertex=4 +DEAL::Cell: 2.3 +DEAL::Left vertex=4 +DEAL::Right vertex=1 +DEAL::Coarse mesh: +DEAL::Left vertex=0 +DEAL::Right vertex=1 +DEAL::Cell: 2.0 +DEAL::Left vertex=0 +DEAL::Right vertex=3 +DEAL::Cell: 2.1 +DEAL::Left vertex=3 +DEAL::Right vertex=2 +DEAL::Cell: 2.2 +DEAL::Left vertex=2 +DEAL::Right vertex=4 +DEAL::Cell: 2.3 +DEAL::Left vertex=4 +DEAL::Right vertex=1 diff --git a/tests/deal.II/vertex_as_face_03.cc b/tests/deal.II/vertex_as_face_03.cc new file mode 100644 index 0000000000..6b186d5400 --- /dev/null +++ b/tests/deal.II/vertex_as_face_03.cc @@ -0,0 +1,61 @@ +//---------------------------- vertex_as_face_03.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2010 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- vertex_as_face_03.cc --------------------------- + +// verify that we can do things like cell->face() in 1d as well. here: +// test cell->face()->at_boundary() + + +#include "../tests.h" +#include +#include +#include +#include + +#include + + +template +void test () +{ + Triangulation<1,spacedim> tria; + GridGenerator::hyper_cube (tria); + + deallog << "Coarse mesh:" << std::endl; + deallog << "Left vertex=" << tria.begin_active()->face(0)->at_boundary() << std::endl; + deallog << "Right vertex=" << tria.begin_active()->face(1)->at_boundary() << std::endl; + + tria.refine_global (2); + + for (typename Triangulation<1,spacedim>::active_cell_iterator + cell = tria.begin_active(); + cell != tria.end(); ++cell) + { + deallog << "Cell: " << cell << std::endl; + deallog << "Left vertex=" << cell->face(0)->at_boundary() << std::endl; + deallog << "Right vertex=" << cell->face(1)->at_boundary() << std::endl; + } +} + + + +int main () +{ + std::ofstream logfile("vertex_as_face_03/output"); + deallog.attach(logfile); + deallog.depth_console(0); + + test<1> (); + test<2> (); + + return 0; +} diff --git a/tests/deal.II/vertex_as_face_03/cmp/generic b/tests/deal.II/vertex_as_face_03/cmp/generic new file mode 100644 index 0000000000..15796d162c --- /dev/null +++ b/tests/deal.II/vertex_as_face_03/cmp/generic @@ -0,0 +1,31 @@ + +DEAL::Coarse mesh: +DEAL::Left vertex=1 +DEAL::Right vertex=1 +DEAL::Cell: 2.0 +DEAL::Left vertex=1 +DEAL::Right vertex=0 +DEAL::Cell: 2.1 +DEAL::Left vertex=0 +DEAL::Right vertex=0 +DEAL::Cell: 2.2 +DEAL::Left vertex=0 +DEAL::Right vertex=0 +DEAL::Cell: 2.3 +DEAL::Left vertex=0 +DEAL::Right vertex=1 +DEAL::Coarse mesh: +DEAL::Left vertex=1 +DEAL::Right vertex=1 +DEAL::Cell: 2.0 +DEAL::Left vertex=1 +DEAL::Right vertex=0 +DEAL::Cell: 2.1 +DEAL::Left vertex=0 +DEAL::Right vertex=0 +DEAL::Cell: 2.2 +DEAL::Left vertex=0 +DEAL::Right vertex=0 +DEAL::Cell: 2.3 +DEAL::Left vertex=0 +DEAL::Right vertex=1 diff --git a/tests/deal.II/vertex_as_face_04.cc b/tests/deal.II/vertex_as_face_04.cc new file mode 100644 index 0000000000..619de73346 --- /dev/null +++ b/tests/deal.II/vertex_as_face_04.cc @@ -0,0 +1,61 @@ +//---------------------------- vertex_as_face_04.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2010 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- vertex_as_face_04.cc --------------------------- + +// verify that we can do things like cell->face() in 1d as well. here: +// test vertex location + + +#include "../tests.h" +#include +#include +#include +#include + +#include + + +template +void test () +{ + Triangulation<1,spacedim> tria; + GridGenerator::hyper_cube (tria); + + deallog << "Coarse mesh:" << std::endl; + deallog << "Left vertex=" << (int)tria.begin_active()->face(0)->boundary_indicator() << std::endl; + deallog << "Right vertex=" << (int)tria.begin_active()->face(1)->boundary_indicator() << std::endl; + + tria.refine_global (2); + + for (typename Triangulation<1,spacedim>::active_cell_iterator + cell = tria.begin_active(); + cell != tria.end(); ++cell) + { + deallog << "Cell: " << cell << std::endl; + deallog << "Left vertex=" << (int)cell->face(0)->boundary_indicator() << std::endl; + deallog << "Right vertex=" << (int)cell->face(1)->boundary_indicator() << std::endl; + } +} + + + +int main () +{ + std::ofstream logfile("vertex_as_face_04/output"); + deallog.attach(logfile); + deallog.depth_console(0); + + test<1> (); + test<2> (); + + return 0; +} diff --git a/tests/deal.II/vertex_as_face_04/cmp/generic b/tests/deal.II/vertex_as_face_04/cmp/generic new file mode 100644 index 0000000000..9fe3657a33 --- /dev/null +++ b/tests/deal.II/vertex_as_face_04/cmp/generic @@ -0,0 +1,31 @@ + +DEAL::Coarse mesh: +DEAL::Left vertex=0 +DEAL::Right vertex=1 +DEAL::Cell: 2.0 +DEAL::Left vertex=0 +DEAL::Right vertex=255 +DEAL::Cell: 2.1 +DEAL::Left vertex=255 +DEAL::Right vertex=255 +DEAL::Cell: 2.2 +DEAL::Left vertex=255 +DEAL::Right vertex=255 +DEAL::Cell: 2.3 +DEAL::Left vertex=255 +DEAL::Right vertex=1 +DEAL::Coarse mesh: +DEAL::Left vertex=0 +DEAL::Right vertex=1 +DEAL::Cell: 2.0 +DEAL::Left vertex=0 +DEAL::Right vertex=255 +DEAL::Cell: 2.1 +DEAL::Left vertex=255 +DEAL::Right vertex=255 +DEAL::Cell: 2.2 +DEAL::Left vertex=255 +DEAL::Right vertex=255 +DEAL::Cell: 2.3 +DEAL::Left vertex=255 +DEAL::Right vertex=1