17 #ifndef __deal2__tria_accessor_h
18 #define __deal2__tria_accessor_h
21 #include <deal.II/base/config.h>
23 #include <deal.II/base/geometry_info.h>
24 #include <deal.II/base/point.h>
25 #include <deal.II/grid/tria_iterator_base.h>
26 #include <deal.II/grid/tria_iterator_selector.h>
27 #include <deal.II/grid/cell_id.h>
39 template <
int dim,
int spacedim>
class Boundary;
48 struct Implementation;
53 struct Implementation;
95 void operator ++ ()
const
100 void operator -- ()
const
125 template <
int dim,
int spacedim>
class TriaAccessor<0, dim, spacedim>;
126 template <
int spacedim>
class TriaAccessor<0, 1, spacedim>;
182 <<
"You can only set the child index if the cell has no "
183 <<
"children, or clear it. The given "
184 <<
"index was " << arg1 <<
" (-1 means: clear children)");
233 <<
"You can only set the child index of an even numbered child."
234 <<
"The number of the child given was " << arg1 <<
".");
263 template <
int structdim,
int dim,
int spacedim=dim>
321 const int level = -1,
322 const int index = -1,
323 const AccessorData * = 0);
516 typename ::internal::TriaAccessor::PresentLevelType<structdim,dim>::type
present_level;
563 template <
int structdim,
int dim,
int spacedim=dim>
587 const int level = -1,
588 const int index = -1,
589 const AccessorData *local_data = 0);
613 template <
typename OtherAccessor>
632 void operator -- ()
const;
667 template <
int structdim,
int dim,
int spacedim>
683 const int level = -1,
684 const int index = -1,
708 template <
int structdim2,
int dim2,
int spacedim2>
717 template <
int structdim2,
int dim2,
int spacedim2>
789 typename ::internal::Triangulation::Iterators<dim,spacedim>::line_iterator
790 line (
const unsigned int i)
const;
801 unsigned int line_index (
const unsigned int i)
const;
807 typename ::internal::Triangulation::Iterators<dim,spacedim>::quad_iterator
808 quad (
const unsigned int i)
const;
819 unsigned int quad_index (
const unsigned int i)
const;
869 bool face_flip (
const unsigned int face)
const;
984 child (
const unsigned int i)
const;
1506 void set (const ::internal::Triangulation::TriaObject<structdim> &o)
const;
1519 const bool orientation)
const;
1541 const bool orientation)
const;
1554 const bool flip)
const;
1567 const bool rotation)
const;
1663 friend struct ::internal::Triangulation::Implementation;
1664 friend struct ::internal::TriaAccessor::Implementation;
1678 template<
int dim,
int spacedim>
1701 template <
int spacedim>
1780 const VertexKind vertex_kind,
1795 const AccessorData * = 0);
1802 template <
int structdim2,
int dim2,
int spacedim2>
1810 template <
int structdim2,
int dim2,
int spacedim2>
1842 static int level ();
1927 unsigned int vertex_index (
const unsigned int i = 0)
const;
1952 typename ::internal::Triangulation::Iterators<1,spacedim>::line_iterator
1953 static line (
const unsigned int);
1964 static unsigned int line_index (
const unsigned int i);
1971 typename ::internal::Triangulation::Iterators<1,spacedim>::quad_iterator
1972 quad (
const unsigned int i);
1983 static unsigned int quad_index (
const unsigned int i);
2031 static bool face_flip (
const unsigned int face);
2086 child (
const unsigned int);
2221 template <
int dim,
int spacedim=dim>
2248 const int level = -1,
2249 const int index = -1,
2250 const AccessorData *local_data = 0);
2278 template <
int structdim2,
int dim2,
int spacedim2>
2287 template <
int structdim2,
int dim2,
int spacedim2>
2308 child (
const unsigned int i)
const;
2315 face (
const unsigned int i)
const;
2424 const unsigned int subface_no)
const;
2438 neighbor (
const unsigned int i)
const;
2520 std::pair<unsigned int, unsigned int>
2868 bool operator < (const CellAccessor<dim, spacedim> &other)
const;
3037 DeclException0 (ExcRefineCellNotActive);
3041 DeclException0 (ExcCellFlaggedForRefinement);
3045 DeclException0 (ExcCellFlaggedForCoarsening);
3089 template<
int dim_,
int spacedim_ >
3127 friend struct ::internal::Triangulation::Implementation;
3136 template <
int structdim,
int dim,
int spacedim>
3137 template <
typename OtherAccessor>
3142 ExcMessage (
"You are attempting an illegal conversion between "
3143 "iterator/accessor types. The constructor you call "
3144 "only exists to make certain template constructs "
3145 "easier to write as dimension independent code but "
3146 "the conversion is not valid in the current context."));
3151 template <
int structdim,
int dim,
int spacedim>
3152 template <
int structdim2,
int dim2,
int spacedim2>
3157 ExcMessage (
"You are attempting an illegal conversion between "
3158 "iterator/accessor types. The constructor you call "
3159 "only exists to make certain template constructs "
3160 "easier to write as dimension independent code but "
3161 "the conversion is not valid in the current context."));
3166 template <
int dim,
int spacedim>
3167 template <
int structdim2,
int dim2,
int spacedim2>
3172 ExcMessage (
"You are attempting an illegal conversion between "
3173 "iterator/accessor types. The constructor you call "
3174 "only exists to make certain template constructs "
3175 "easier to write as dimension independent code but "
3176 "the conversion is not valid in the current context."));
3181 template <
int structdim,
int dim,
int spacedim>
3182 template <
int structdim2,
int dim2,
int spacedim2>
3187 ExcMessage (
"You are attempting an illegal conversion between "
3188 "iterator/accessor types. The constructor you call "
3189 "only exists to make certain template constructs "
3190 "easier to write as dimension independent code but "
3191 "the conversion is not valid in the current context."));
3196 template <
int dim,
int spacedim>
3197 template <
int structdim2,
int dim2,
int spacedim2>
3202 ExcMessage (
"You are attempting an illegal conversion between "
3203 "iterator/accessor types. The constructor you call "
3204 "only exists to make certain template constructs "
3205 "easier to write as dimension independent code but "
3206 "the conversion is not valid in the current context."));
3209 template <
int dim,
int spacedim>
3213 std::vector<unsigned char> id(this->level(), -1);
3214 unsigned int coarse_index;
3217 while (ptr.
level()>0)
3221 for (
unsigned int c=0; c<ptr.
parent()->n_children(); ++c)
3231 id[ptr.
level()-1] = v;
3237 coarse_index = ptr.
index();
3239 return CellId(coarse_index,
id);
3255 DEAL_II_NAMESPACE_CLOSE
3258 # include "tria_accessor.templates.h"
void set_subdomain_id(const types::subdomain_id new_subdomain_id) const
bool flag_for_line_refinement(const unsigned int line_no) const
void set_neighbor(const unsigned int i, const TriaIterator< CellAccessor< dim, spacedim > > &pointer) const
const Triangulation< dim, spacedim > & get_triangulation() const
TriaAccessor< dim, dim, spacedim >::AccessorData AccessorData
TriaIterator< CellAccessor< dim, spacedim > > neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
unsigned int user_index() const
void set_user_index(const unsigned int p) const
void set_user_pointer(void *p) const
void set(const ::internal::Triangulation::TriaObject< structdim > &o) const
void clear_children() const
void clear_user_flag() const
bool operator<(const TriaAccessorBase &other) const
TriaIterator< CellAccessor< dim, spacedim > > child(const unsigned int i) const
void clear_user_pointer() const
std::pair< unsigned int, unsigned int > neighbor_of_coarser_neighbor(const unsigned int neighbor) const
unsigned char material_id
void recursively_set_user_pointer(void *p) const
TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
TriaAccessorBase & operator=(const TriaAccessorBase &)
::ExceptionBase & ExcMessage(std::string arg1)
unsigned int neighbor_of_neighbor_internal(const unsigned int neighbor) const
bool face_orientation(const unsigned int face) const
unsigned int neighbor_face_no(const unsigned int neighbor) const
TriaIterator< TriaAccessor< dim-1, dim, spacedim > > face(const unsigned int i) const
unsigned int neighbor_of_neighbor(const unsigned int neighbor) const
void copy_from(const TriaAccessorBase &)
::internal::SubfaceCase< dim > subface_case(const unsigned int face_no) const
TriaIterator< CellAccessor< dim, spacedim > > neighbor(const unsigned int i) const
void set_face_flip(const unsigned int face, const bool flip) const
bool line_orientation(const unsigned int line) const
int neighbor_level(const unsigned int i) const
double minimum_vertex_distance() const
const Triangulation< 1, spacedim > * tria
const Boundary< dim, spacedim > & get_boundary() const
typename::internal::Triangulation::Iterators< dim, spacedim >::line_iterator line(const unsigned int i) const
types::subdomain_id level_subdomain_id() const
bool face_flip(const unsigned int face) const
void set_user_flag() const
bool has_children() const
bool is_locally_owned() const
void copy_from(const InvalidAccessor &)
bool is_translation_of(const TriaIterator< TriaAccessor< structdim, dim, spacedim > > &o) const
void recursively_clear_user_index() const
static const unsigned int space_dimension
unsigned int max_refinement_depth() const
void clear_user_data() const
void clear_used_flag() const
bool has_children() const
TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
void set_refinement_case(const RefinementCase< structdim > &ref_case) const
bool point_inside_codim(const Point< spacedim_ > &p) const
void set_used_flag() const
void operator=(const CellAccessor< dim, spacedim > &)
void set_parent(const unsigned int parent_index)
#define DeclException1(Exception1, type1, outsequence)
void set_direction_flag(const bool new_direction_flag) const
#define Assert(cond, exc)
types::boundary_id boundary_indicator() const
void set_boundary_indicator(const types::boundary_id) const
TriaAccessor(const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *local_data=0)
unsigned int vertex_index(const unsigned int i) const
bool user_flag_set() const
const Triangulation< dim, spacedim > * tria
#define DeclException0(Exception0)
types::material_id material_id() const
void recursively_clear_user_flag() const
bool is_locally_owned_on_level() const
::internal::Triangulation::TriaObjects<::internal::Triangulation::TriaObject< structdim > > & objects() const
bool is_artificial() const
bool coarsen_flag_set() const
TriaIterator< TriaAccessor< structdim, dim, spacedim > > isotropic_child(const unsigned int i) const
void set_all_boundary_indicators(const types::boundary_id) const
typename::internal::TriaAccessor::PresentLevelType< structdim, dim >::type present_level
TriaIterator< CellAccessor< dim, spacedim > > parent() const
void clear_refine_flag() const
void set_refine_flag(const RefinementCase< dim > ref_case=RefinementCase< dim >::isotropic_refinement) const
bool operator==(const TriaAccessorBase &) const
unsigned int face_index(const unsigned int i) const
unsigned int subdomain_id
void set_level_subdomain_id(const types::subdomain_id new_level_subdomain_id) const
void clear_coarsen_flag() const
typename::internal::Triangulation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i) const
bool point_inside(const Point< spacedim > &p) const
Point< spacedim > center() const
int neighbor_index(const unsigned int i) const
bool has_boundary_lines() const
void recursively_set_material_id(const types::material_id new_material_id) const
RefinementCase< dim > refine_flag_set() const
bool neighbor_is_coarser(const unsigned int neighbor) const
void set_face_orientation(const unsigned int face, const bool orientation) const
unsigned int number_of_children() const
TriaIterator< TriaAccessor< structdim, dim, spacedim > > child(const unsigned int i) const
unsigned int quad_index(const unsigned int i) const
void clear_refinement_case() const
bool flag_for_face_refinement(const unsigned int face_no, const RefinementCase< dim-1 > &face_refinement_case=RefinementCase< dim-1 >::isotropic_refinement) const
unsigned int n_children() const
IteratorState::IteratorStates state() const
void recursively_set_user_index(const unsigned int p) const
void set_coarsen_flag() const
static const unsigned int structure_dimension
CellAccessor(const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *local_data=0)
void set_face_rotation(const unsigned int face, const bool rotation) const
int isotropic_child_index(const unsigned int i) const
void set_children(const unsigned int i, const int index) const
Triangulation< dim, spacedim > Container
void * user_pointer() const
void recursively_clear_user_pointer() const
void recursively_set_subdomain_id(const types::subdomain_id new_subdomain_id) const
int child_index(const unsigned int i) const
InvalidAccessor(const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *local_data=0)
double extent_in_direction(const unsigned int axis) const
unsigned char boundary_id
types::subdomain_id subdomain_id() const
void recursively_set_user_flag() const
static const unsigned int dimension
unsigned int global_vertex_index
void clear_user_index() const
TriaAccessorBase(const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *=0)
::ExceptionBase & ExcInternalError()
void set_material_id(const types::material_id new_material_id) const
void operator=(const TriaAccessor &)
bool operator==(const InvalidAccessor &) const
Point< spacedim > barycenter() const
bool face_rotation(const unsigned int face) const
RefinementCase< structdim > refinement_case() const
void set_line_orientation(const unsigned int line, const bool orientation) const
bool direction_flag() const
bool operator!=(const TriaAccessorBase &) const
Point< spacedim > & vertex(const unsigned int i) const
unsigned int line_index(const unsigned int i) const