16 #ifndef dealii_tria_accessor_h 17 #define dealii_tria_accessor_h 20 #include <deal.II/base/config.h> 22 #include <deal.II/base/bounding_box.h> 23 #include <deal.II/base/exceptions.h> 24 #include <deal.II/base/geometry_info.h> 25 #include <deal.II/base/point.h> 27 #include <deal.II/grid/cell_id.h> 28 #include <deal.II/grid/tria_iterator_base.h> 29 #include <deal.II/grid/tria_iterator_selector.h> 34 DEAL_II_NAMESPACE_OPEN
36 template <
int dim,
int spacedim>
38 template <
typename Accessor>
40 template <
typename Accessor>
42 template <
typename Accessor>
45 template <
int dim,
int spacedim>
51 namespace TriangulationImplementation
60 namespace TriaAccessorImplementation
69 template <
int structdim,
int dim>
124 template <
int structdim,
int dim,
int spacedim>
126 template <
int dim,
int spacedim>
128 template <
int spacedim>
146 "The operation you are attempting can only be performed for " 147 "(cell, face, or edge) iterators that point to valid " 148 "objects. These objects need not necessarily be active, " 149 "i.e., have no children, but they need to be part of a " 150 "triangulation. (The objects pointed to by an iterator " 151 "may -- after coarsening -- also be objects that used " 152 "to be part of a triangulation, but are now no longer " 153 "used. Their memory location may have been retained " 154 "for re-use upon the next mesh refinement, but is " 155 "currently unused.)");
166 "The operation you are attempting can only be performed for " 167 "(cell, face, or edge) iterators that point to 'active' " 168 "objects. 'Active' objects are those that do not have " 169 "children (in the case of cells), or that are part of " 170 "an active cell (in the case of faces or edges). However, " 171 "the object on which you are trying the current " 172 "operation is not 'active' in this sense.");
179 "The operation you are attempting can only be performed for " 180 "(cell, face, or edge) iterators that have children, " 181 "but the object on which you are trying the current " 182 "operation does not have any.");
190 "The operation you are attempting can only be performed for " 191 "(cell, face, or edge) iterators that have a parent object, " 192 "but the object on which you are trying the current " 193 "operation does not have one -- i.e., it is on the " 194 "coarsest level of the triangulation.");
200 <<
"You can only set the child index if the cell does not " 201 <<
"currently have children registered; or you can clear it. " 202 <<
"The given index was " << arg1
203 <<
" (-1 means: clear children).");
207 template <
typename AccessorType>
210 <<
"You tried to dereference an iterator for which this " 211 <<
"is not possible. More information on this iterator: " 212 <<
"index=" << arg1.index() <<
", state=" 222 "Iterators can only be compared if they point to the same " 223 "triangulation, or if neither of them are associated " 224 "with a triangulation.");
257 <<
"You can only set the child index of an even numbered child." 258 <<
"The number of the child given was " << arg1 <<
".");
288 template <
int structdim,
int dim,
int spacedim = dim>
297 static const unsigned int space_dimension = spacedim;
304 static const unsigned int dimension = dim;
311 static const unsigned int structure_dimension = structdim;
325 const int level = -1,
326 const int index = -1,
495 get_triangulation()
const;
505 typename ::internal::TriaAccessorImplementation::
520 template <
typename Accessor>
522 template <
typename Accessor>
524 template <
typename Accessor>
551 template <
int structdim,
int dim,
int spacedim = dim>
569 const int level = -1,
570 const int index = -1,
586 template <
typename OtherAccessor>
623 has_children()
const;
635 vertex(
const unsigned int i)
const;
641 typename ::internal::TriangulationImplementation::
642 Iterators<dim, spacedim>::line_iterator
643 line(
const unsigned int i)
const;
649 typename ::internal::TriangulationImplementation::
650 Iterators<dim, spacedim>::quad_iterator
651 quad(
const unsigned int i)
const;
674 template <
int structdim,
int dim,
int spacedim>
688 const int level = -1,
689 const int index = -1,
704 template <
int structdim2,
int dim2,
int spacedim2>
711 template <
int structdim2,
int dim2,
int spacedim2>
747 vertex_iterator(
const unsigned int i)
const;
765 vertex_index(
const unsigned int i)
const;
805 vertex(
const unsigned int i)
const;
810 typename ::internal::TriangulationImplementation::
811 Iterators<dim, spacedim>::line_iterator
812 line(
const unsigned int i)
const;
821 line_index(
const unsigned int i)
const;
826 typename ::internal::TriangulationImplementation::
827 Iterators<dim, spacedim>::quad_iterator
828 quad(
const unsigned int i)
const;
837 quad_index(
const unsigned int i)
const;
861 face_orientation(
const unsigned int face)
const;
873 face_flip(
const unsigned int face)
const;
885 face_rotation(
const unsigned int face)
const;
898 line_orientation(
const unsigned int line)
const;
914 has_children()
const;
937 number_of_children()
const;
953 max_refinement_depth()
const;
959 child(
const unsigned int i)
const;
970 isotropic_child(
const unsigned int i)
const;
976 refinement_case()
const;
984 child_index(
const unsigned int i)
const;
992 isotropic_child_index(
const unsigned int i)
const;
1015 boundary_id()
const;
1090 at_boundary()
const;
1102 get_manifold()
const;
1125 manifold_id()
const;
1180 user_flag_set()
const;
1188 set_user_flag()
const;
1196 clear_user_flag()
const;
1204 recursively_set_user_flag()
const;
1212 recursively_clear_user_flag()
const;
1220 clear_user_data()
const;
1234 set_user_pointer(
void *p)
const;
1242 clear_user_pointer()
const;
1260 user_pointer()
const;
1284 recursively_set_user_pointer(
void *p)
const;
1293 recursively_clear_user_pointer()
const;
1305 set_user_index(
const unsigned int p)
const;
1313 clear_user_index()
const;
1347 recursively_set_user_index(
const unsigned int p)
const;
1358 recursively_clear_user_index()
const;
1406 std::pair<Point<spacedim>,
double>
1407 enclosing_ball()
const;
1413 bounding_box()
const;
1425 extent_in_direction(
const unsigned int axis)
const;
1431 minimum_vertex_distance()
const;
1472 real_to_unit_cell_affine_approximation(
const Point<spacedim> &point)
const;
1509 center(
const bool respect_manifold =
false,
1510 const bool interpolate_from_surrounding =
false)
const;
1596 set(const ::internal::TriangulationImplementation::TriaObject<structdim>
1607 set_line_orientation(
const unsigned int line,
const bool orientation)
const;
1620 set_face_orientation(
const unsigned int face,
const bool orientation)
const;
1629 set_face_flip(
const unsigned int face,
const bool flip)
const;
1638 set_face_rotation(
const unsigned int face,
const bool rotation)
const;
1644 set_used_flag()
const;
1650 clear_used_flag()
const;
1671 clear_refinement_case()
const;
1680 set_children(
const unsigned int i,
const int index)
const;
1687 clear_children()
const;
1693 friend struct ::internal::TriangulationImplementation::Implementation;
1694 friend struct ::internal::TriaAccessorImplementation::Implementation;
1718 template <
int dim,
int spacedim>
1727 static const unsigned int space_dimension = spacedim;
1734 static const unsigned int dimension = dim;
1741 static const unsigned int structure_dimension = 0;
1753 const unsigned int vertex_index);
1761 const int level = 0,
1762 const int index = 0,
1768 template <
int structdim2,
int dim2,
int spacedim2>
1774 template <
int structdim2,
int dim2,
int spacedim2>
1854 vertex_index(
const unsigned int i = 0)
const;
1862 vertex(
const unsigned int i = 0)
const;
1868 typename ::internal::TriangulationImplementation::
1869 Iterators<dim, spacedim>::line_iterator
static line(
const unsigned int);
1875 line_index(
const unsigned int i);
1880 static typename ::internal::TriangulationImplementation::
1881 Iterators<dim, spacedim>::quad_iterator
1882 quad(
const unsigned int i);
1888 quad_index(
const unsigned int i);
1916 extent_in_direction(
const unsigned int axis)
const;
1926 center(
const bool respect_manifold =
false,
1927 const bool interpolate_from_surrounding =
false)
const;
1954 face_orientation(
const unsigned int face);
1960 face_flip(
const unsigned int face);
1966 face_rotation(
const unsigned int face);
1972 line_orientation(
const unsigned int line);
2003 number_of_children();
2009 max_refinement_depth();
2015 child(
const unsigned int);
2021 isotropic_child(
const unsigned int);
2033 child_index(
const unsigned int i);
2039 isotropic_child_index(
const unsigned int i);
2072 template <
typename Accessor>
2074 template <
typename Accessor>
2076 template <
typename Accessor>
2099 template <
int spacedim>
2108 static const unsigned int space_dimension = spacedim;
2115 static const unsigned int dimension = 1;
2122 static const unsigned int structure_dimension = 0;
2162 const unsigned int vertex_index);
2177 template <
int structdim2,
int dim2,
int spacedim2>
2183 template <
int structdim2,
int dim2,
int spacedim2>
2273 vertex_index(
const unsigned int i = 0)
const;
2281 vertex(
const unsigned int i = 0)
const;
2294 typename ::internal::TriangulationImplementation::
2295 Iterators<1, spacedim>::line_iterator
static line(
const unsigned int);
2304 line_index(
const unsigned int i);
2309 static typename ::internal::TriangulationImplementation::
2310 Iterators<1, spacedim>::quad_iterator
2311 quad(
const unsigned int i);
2320 quad_index(
const unsigned int i);
2332 at_boundary()
const;
2349 boundary_id()
const;
2355 get_manifold()
const;
2364 manifold_id()
const;
2378 face_orientation(
const unsigned int face);
2384 face_flip(
const unsigned int face);
2390 face_rotation(
const unsigned int face);
2396 line_orientation(
const unsigned int line);
2427 number_of_children();
2433 max_refinement_depth();
2439 child(
const unsigned int);
2445 isotropic_child(
const unsigned int);
2457 child_index(
const unsigned int i);
2463 isotropic_child_index(
const unsigned int i);
2581 template <
int dim,
int spacedim = dim>
2606 const int level = -1,
2607 const int index = -1,
2627 template <
int structdim2,
int dim2,
int spacedim2>
2634 template <
int structdim2,
int dim2,
int spacedim2>
2665 child(
const unsigned int i)
const;
2671 face(
const unsigned int i)
const;
2683 face_index(
const unsigned int i)
const;
2733 neighbor_child_on_subface(
const unsigned int face_no,
2734 const unsigned int subface_no)
const;
2758 neighbor(
const unsigned int i)
const;
2765 neighbor_index(
const unsigned int i)
const;
2772 neighbor_level(
const unsigned int i)
const;
2786 neighbor_of_neighbor(
const unsigned int neighbor)
const;
2799 neighbor_is_coarser(
const unsigned int neighbor)
const;
2815 std::pair<unsigned int, unsigned int>
2816 neighbor_of_coarser_neighbor(
const unsigned int neighbor)
const;
2825 neighbor_face_no(
const unsigned int neighbor)
const;
2847 has_periodic_neighbor(
const unsigned int i)
const;
2866 periodic_neighbor(
const unsigned int i)
const;
2876 neighbor_or_periodic_neighbor(
const unsigned int i)
const;
2893 periodic_neighbor_child_on_subface(
const unsigned int face_no,
2894 const unsigned int subface_no)
const;
2906 std::pair<unsigned int, unsigned int>
2907 periodic_neighbor_of_coarser_periodic_neighbor(
const unsigned face_no)
const;
2915 periodic_neighbor_index(
const unsigned int i)
const;
2923 periodic_neighbor_level(
const unsigned int i)
const;
2940 periodic_neighbor_of_periodic_neighbor(
const unsigned int i)
const;
2948 periodic_neighbor_face_no(
const unsigned int i)
const;
2957 periodic_neighbor_is_coarser(
const unsigned int i)
const;
2976 at_boundary(
const unsigned int i)
const;
2987 at_boundary()
const;
2997 has_boundary_lines()
const;
3024 refine_flag_set()
const;
3051 clear_refine_flag()
const;
3061 flag_for_face_refinement(
3062 const unsigned int face_no,
3072 flag_for_line_refinement(
const unsigned int line_no)
const;
3083 subface_case(
const unsigned int face_no)
const;
3089 coarsen_flag_set()
const;
3096 set_coarsen_flag()
const;
3102 clear_coarsen_flag()
const;
3126 material_id()
const;
3178 subdomain_id()
const;
3203 level_subdomain_id()
const;
3210 set_level_subdomain_id(
3230 recursively_set_subdomain_id(
3250 direction_flag()
const;
3278 active_cell_index()
const;
3288 parent_index()
const;
3340 is_locally_owned()
const;
3347 is_locally_owned_on_level()
const;
3402 is_artificial()
const;
3429 set_neighbor(
const unsigned int i,
3483 neighbor_of_neighbor_internal(
const unsigned int neighbor)
const;
3490 template <
int dim_,
int spacedim_>
3502 set_active_cell_index(
const unsigned int active_cell_index);
3508 set_parent(
const unsigned int parent_index);
3517 set_direction_flag(
const bool new_direction_flag)
const;
3522 friend struct ::internal::TriangulationImplementation::Implementation;
3530 template <
int structdim,
int dim,
int spacedim>
3531 template <
typename OtherAccessor>
3533 const OtherAccessor &)
3536 ExcMessage(
"You are attempting an illegal conversion between " 3537 "iterator/accessor types. The constructor you call " 3538 "only exists to make certain template constructs " 3539 "easier to write as dimension independent code but " 3540 "the conversion is not valid in the current context."));
3545 template <
int structdim,
int dim,
int spacedim>
3546 template <
int structdim2,
int dim2,
int spacedim2>
3551 ExcMessage(
"You are attempting an illegal conversion between " 3552 "iterator/accessor types. The constructor you call " 3553 "only exists to make certain template constructs " 3554 "easier to write as dimension independent code but " 3555 "the conversion is not valid in the current context."));
3560 template <
int dim,
int spacedim>
3561 template <
int structdim2,
int dim2,
int spacedim2>
3566 ExcMessage(
"You are attempting an illegal conversion between " 3567 "iterator/accessor types. The constructor you call " 3568 "only exists to make certain template constructs " 3569 "easier to write as dimension independent code but " 3570 "the conversion is not valid in the current context."));
3575 template <
int structdim,
int dim,
int spacedim>
3576 template <
int structdim2,
int dim2,
int spacedim2>
3581 ExcMessage(
"You are attempting an illegal conversion between " 3582 "iterator/accessor types. The constructor you call " 3583 "only exists to make certain template constructs " 3584 "easier to write as dimension independent code but " 3585 "the conversion is not valid in the current context."));
3590 template <
int dim,
int spacedim>
3591 template <
int structdim2,
int dim2,
int spacedim2>
3596 ExcMessage(
"You are attempting an illegal conversion between " 3597 "iterator/accessor types. The constructor you call " 3598 "only exists to make certain template constructs " 3599 "easier to write as dimension independent code but " 3600 "the conversion is not valid in the current context."));
3632 DEAL_II_NAMESPACE_CLOSE
3635 #include "tria_accessor.templates.h"
static::ExceptionBase & ExcCellHasNoChildren()
InvalidAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
static::ExceptionBase & ExcNeighborIsNotCoarser()
const Triangulation< dim, spacedim > * tria
unsigned int global_vertex_index
static::ExceptionBase & ExcNeighborIsCoarser()
static::ExceptionBase & ExcDereferenceInvalidObject(AccessorType arg1)
const Triangulation< 1, spacedim > * tria
static::ExceptionBase & ExcFacesHaveNoLevel()
static::ExceptionBase & ExcCellNotUsed()
static::ExceptionBase & ExcMessage(std::string arg1)
unsigned int subdomain_id
#define DeclException1(Exception1, type1, outsequence)
const Triangulation< dim, spacedim > * tria
static::ExceptionBase & ExcSetOnlyEvenChildren(int arg1)
#define Assert(cond, exc)
#define DeclExceptionMsg(Exception, defaulttext)
#define DeclException0(Exception0)
static::ExceptionBase & ExcCellNotActive()
static::ExceptionBase & ExcNoPeriodicNeighbor()
static::ExceptionBase & ExcCantSetChildren(int arg1)
bool point_inside(const Point< spacedim > &p) const
CellAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
static::ExceptionBase & ExcCellHasNoParent()
TriaAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
Iterator reached end of container.
Iterator points to a valid object.
unsigned int global_vertex_index
void set_all_manifold_ids(const types::manifold_id) const
typename::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
static::ExceptionBase & ExcInternalError()
static::ExceptionBase & ExcCantCompareIterators()