Reference documentation for deal.II version GIT relicensing-214-g6e74dec06b 2024-03-27 18:10:01+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | Private Member Functions | Friends | List of all members

#include <deal.II/grid/tria_accessor.h>

Inheritance diagram for CellAccessor< dim, spacedim >:
Inheritance graph
[legend]

Public Types

using AccessorData = typename TriaAccessor< dim, dim, spacedim >::AccessorData
 
using Container = Triangulation< dim, spacedim >
 
using LocalData = void *
 

Public Member Functions

types::global_cell_index global_active_cell_index () const
 
types::global_cell_index global_level_cell_index () const
 
bool point_inside (const Point< 1 > &p) const
 
bool point_inside (const Point< 2 > &p) const
 
bool point_inside (const Point< 3 > &p) const
 
bool point_inside (const Point< 2 > &p) const
 
bool point_inside (const Point< 3 > &p) const
 
bool point_inside (const Point< 3 > &p) const
 
bool used () const
 
void set_all_manifold_ids (const types::manifold_id manifold_ind) const
 
double extent_in_direction (const unsigned int axis) const
 
double extent_in_direction (const unsigned int axis) const
 
double extent_in_direction (const unsigned int axis) const
 
double extent_in_direction (const unsigned int axis) const
 
double extent_in_direction (const unsigned int axis) const
 
Constructors
 CellAccessor (const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
 
 CellAccessor (const TriaAccessor< dim, dim, spacedim > &cell_accessor)
 
template<int structdim2, int dim2, int spacedim2>
 CellAccessor (const InvalidAccessor< structdim2, dim2, spacedim2 > &)
 
template<int structdim2, int dim2, int spacedim2>
 CellAccessor (const TriaAccessor< structdim2, dim2, spacedim2 > &)
 
 CellAccessor (const CellAccessor< dim, spacedim > &)=default
 
 CellAccessor (CellAccessor< dim, spacedim > &&)=default
 
 ~CellAccessor ()=default
 
CellAccessor< dim, spacedim > & operator= (const CellAccessor< dim, spacedim > &)=delete
 
CellAccessor< dim, spacedim > & operator= (CellAccessor< dim, spacedim > &&)=default
 
Converting iterators
TriaActiveIterator< DoFCellAccessor< dim, spacedim, false > > as_dof_handler_iterator (const DoFHandler< dim, spacedim > &dof_handler) const
 
TriaIterator< DoFCellAccessor< dim, spacedim, true > > as_dof_handler_level_iterator (const DoFHandler< dim, spacedim > &dof_handler) const
 
Dealing with periodic neighbors
bool has_periodic_neighbor (const unsigned int i) const
 
TriaIterator< CellAccessor< dim, spacedim > > periodic_neighbor (const unsigned int i) const
 
TriaIterator< CellAccessor< dim, spacedim > > neighbor_or_periodic_neighbor (const unsigned int i) const
 
TriaIterator< CellAccessor< dim, spacedim > > periodic_neighbor_child_on_subface (const unsigned int face_no, const unsigned int subface_no) const
 
std::pair< unsigned int, unsigned intperiodic_neighbor_of_coarser_periodic_neighbor (const unsigned face_no) const
 
int periodic_neighbor_index (const unsigned int i) const
 
int periodic_neighbor_level (const unsigned int i) const
 
unsigned int periodic_neighbor_of_periodic_neighbor (const unsigned int i) const
 
unsigned int periodic_neighbor_face_no (const unsigned int i) const
 
bool periodic_neighbor_is_coarser (const unsigned int i) const
 
Dealing with boundary indicators
bool at_boundary (const unsigned int i) const
 
bool at_boundary () const
 
bool has_boundary_lines () const
 
Dealing with refinement indicators
RefinementCase< dim > refine_flag_set () const
 
void set_refine_flag (const RefinementCase< dim > ref_case=RefinementCase< dim >::isotropic_refinement) const
 
void clear_refine_flag () const
 
bool flag_for_face_refinement (const unsigned int face_no, const RefinementCase< dim - 1 > &face_refinement_case=RefinementCase< dim - 1 >::isotropic_refinement) const
 
bool flag_for_line_refinement (const unsigned int line_no) const
 
::internal::SubfaceCase< dim > subface_case (const unsigned int face_no) const
 
bool coarsen_flag_set () const
 
void set_coarsen_flag () const
 
void clear_coarsen_flag () const
 
Dealing with material indicators
types::material_id material_id () const
 
void set_material_id (const types::material_id new_material_id) const
 
void recursively_set_material_id (const types::material_id new_material_id) const
 
Dealing with subdomain indicators
types::subdomain_id subdomain_id () const
 
void set_subdomain_id (const types::subdomain_id new_subdomain_id) const
 
types::subdomain_id level_subdomain_id () const
 
void set_level_subdomain_id (const types::subdomain_id new_level_subdomain_id) const
 
void recursively_set_subdomain_id (const types::subdomain_id new_subdomain_id) const
 
Dealing with codim 1 cell orientation
bool direction_flag () const
 
unsigned int active_cell_index () const
 
int parent_index () const
 
TriaIterator< CellAccessor< dim, spacedim > > parent () const
 
Other functions
bool is_active () const
 
bool is_locally_owned () const
 
bool is_locally_owned_on_level () const
 
bool is_ghost () const
 
bool is_ghost_on_level () const
 
bool is_artificial () const
 
bool is_artificial_on_level () const
 
bool point_inside (const Point< spacedim > &p) const
 
void set_neighbor (const unsigned int i, const TriaIterator< CellAccessor< dim, spacedim > > &pointer) const
 
CellId id () const
 
double diameter (const Mapping< dim, spacedim > &mapping) const
 
Accessing sub-objects
TriaIterator< TriaAccessor< 0, dim, spacedim > > vertex_iterator (const unsigned int i) const
 
unsigned int vertex_index (const unsigned int i) const
 
Point< spacedim > & vertex (const unsigned int i) const
 
typename::internal::TriangulationImplementation::Iterators< dim, spacedim >::line_iterator line (const unsigned int i) const
 
unsigned int line_index (const unsigned int i) const
 
typename::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator quad (const unsigned int i) const
 
unsigned int quad_index (const unsigned int i) const
 
Orientation of sub-objects
unsigned char combined_face_orientation (const unsigned int face) const
 
bool face_orientation (const unsigned int face) const
 
bool face_flip (const unsigned int face) const
 
bool face_rotation (const unsigned int face) const
 
bool line_orientation (const unsigned int line) const
 
Accessing children
bool has_children () const
 
unsigned int n_children () const
 
unsigned int n_active_descendants () const
 
unsigned int max_refinement_depth () const
 
unsigned int child_iterator_to_index (const TriaIterator< TriaAccessor< structdim, dim, spacedim > > &child) const
 
TriaIterator< TriaAccessor< structdim, dim, spacedim > > isotropic_child (const unsigned int i) const
 
RefinementCase< structdim > refinement_case () const
 
int child_index (const unsigned int i) const
 
int isotropic_child_index (const unsigned int i) const
 
Dealing with boundary indicators
types::boundary_id boundary_id () const
 
void set_boundary_id (const types::boundary_id) const
 
void set_all_boundary_ids (const types::boundary_id) const
 
const Manifold< dim, spacedim > & get_manifold () const
 
Dealing with manifold indicators
types::manifold_id manifold_id () const
 
void set_manifold_id (const types::manifold_id) const
 
void set_all_manifold_ids (const types::manifold_id) const
 
User data
bool user_flag_set () const
 
void set_user_flag () const
 
void clear_user_flag () const
 
void recursively_set_user_flag () const
 
void recursively_clear_user_flag () const
 
void clear_user_data () const
 
void set_user_pointer (void *p) const
 
void clear_user_pointer () const
 
void * user_pointer () const
 
void recursively_set_user_pointer (void *p) const
 
void recursively_clear_user_pointer () const
 
void set_user_index (const unsigned int p) const
 
void clear_user_index () const
 
unsigned int user_index () const
 
void recursively_set_user_index (const unsigned int p) const
 
void recursively_clear_user_index () const
 
Geometric information about an object
double diameter () const
 
std::pair< Point< spacedim >, double > enclosing_ball () const
 
BoundingBox< spacedim > bounding_box () const
 
double extent_in_direction (const unsigned int axis) const
 
double minimum_vertex_distance () const
 
Point< spacedim > intermediate_point (const Point< structdim > &coordinates) const
 
Point< structdim > real_to_unit_cell_affine_approximation (const Point< spacedim > &point) const
 
Point< spacedim > center (const bool respect_manifold=false, const bool interpolate_from_surrounding=false) const
 
Point< spacedim > barycenter () const
 
double measure () const
 
bool is_translation_of (const TriaIterator< TriaAccessor< structdim, dim, spacedim > > &o) const
 
ReferenceCell reference_cell () const
 
unsigned int n_vertices () const
 
unsigned int n_lines () const
 
unsigned int n_faces () const
 
std_cxx20::ranges::iota_view< unsigned int, unsigned intvertex_indices () const
 
std_cxx20::ranges::iota_view< unsigned int, unsigned intline_indices () const
 
std_cxx20::ranges::iota_view< unsigned int, unsigned intface_indices () const
 
Iterator address and state
int level () const
 
int index () const
 
IteratorState::IteratorStates state () const
 
const Triangulation< dim, spacedim > & get_triangulation () const
 

Static Public Member Functions

static ::ExceptionBaseExcRefineCellNotActive ()
 
static ::ExceptionBaseExcCellFlaggedForRefinement ()
 
static ::ExceptionBaseExcCellFlaggedForCoarsening ()
 

Static Public Attributes

static constexpr unsigned int space_dimension
 
static constexpr unsigned int dimension
 
static const unsigned int structure_dimension
 

Protected Member Functions

unsigned int neighbor_of_neighbor_internal (const unsigned int neighbor) const
 
template<int dim_, int spacedim_>
bool point_inside_codim (const Point< spacedim_ > &p) const
 
void copy_from (const TriaAccessorBase &)
 
bool operator< (const TriaAccessorBase &other) const
 
bool operator== (const TriaAccessorBase &) const
 
bool operator!= (const TriaAccessorBase &) const
 
::internal::TriangulationImplementation::TriaObjectsobjects () const
 
Advancement of iterators
void operator++ ()
 
void operator-- ()
 

Protected Attributes

typename::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
 
int present_index
 
const Triangulation< dim, spacedim > * tria
 

Private Member Functions

void set_active_cell_index (const unsigned int active_cell_index) const
 
void set_global_active_cell_index (const types::global_cell_index index) const
 
void set_global_level_cell_index (const types::global_cell_index index) const
 
void set_parent (const unsigned int parent_index)
 
void set_direction_flag (const bool new_direction_flag) const
 
void set_boundary_id_internal (const types::boundary_id id) const
 
void set_bounding_object_indices (const std::initializer_list< int > &new_indices) const
 
void set_bounding_object_indices (const std::initializer_list< unsigned int > &new_indices) const
 
void set_line_orientation (const unsigned int line, const bool orientation) const
 
void set_combined_face_orientation (const unsigned int face, const unsigned char combined_orientation) const
 
void set_used_flag () const
 
void clear_used_flag () const
 
void set_refinement_case (const RefinementCase< structdim > &ref_case) const
 
void clear_refinement_case () const
 
void set_children (const unsigned int i, const int index) const
 
void clear_children () const
 

Friends

class Triangulation< dim, spacedim >
 
class parallel::TriangulationBase< dim, spacedim >
 
struct ::internal::TriangulationImplementation::Implementation
 
struct ::internal::TriangulationImplementation::ImplementationMixedMesh
 

Accessing sub-objects and neighbors

TriaIterator< CellAccessor< dim, spacedim > > child (const unsigned int i) const
 
boost::container::small_vector< TriaIterator< CellAccessor< dim, spacedim > >, GeometryInfo< dim >::max_children_per_cell > child_iterators () const
 
TriaIterator< TriaAccessor< dim - 1, dim, spacedim > > face (const unsigned int i) const
 
unsigned int face_iterator_to_index (const TriaIterator< TriaAccessor< dim - 1, dim, spacedim > > &face) const
 
boost::container::small_vector< TriaIterator< TriaAccessor< dim - 1, dim, spacedim > >, GeometryInfo< dim >::faces_per_cell > face_iterators () const
 
unsigned int face_index (const unsigned int i) const
 
TriaIterator< CellAccessor< dim, spacedim > > neighbor_child_on_subface (const unsigned int face_no, const unsigned int subface_no) const
 
TriaIterator< CellAccessor< dim, spacedim > > neighbor (const unsigned int face_no) const
 
int neighbor_index (const unsigned int face_no) const
 
int neighbor_level (const unsigned int face_no) const
 
unsigned int neighbor_of_neighbor (const unsigned int face_no) const
 
bool neighbor_is_coarser (const unsigned int face_no) const
 
std::pair< unsigned int, unsigned intneighbor_of_coarser_neighbor (const unsigned int neighbor) const
 
unsigned int neighbor_face_no (const unsigned int neighbor) const
 
static bool is_level_cell ()
 

Detailed Description

template<int dim, int spacedim = dim>
class CellAccessor< dim, spacedim >

This class allows access to a cell: a line in one dimension, a quad in two dimension, etc.

The following refers to any dimension:

This class allows access to a cell, which is a line in 1d and a quad in 2d. Cells have more functionality than lines or quads by themselves, for example they can be flagged for refinement, they have neighbors, they have the possibility to check whether they are at the boundary etc. This class offers access to all this data.

Definition at line 3071 of file tria_accessor.h.

Member Typedef Documentation

◆ AccessorData

template<int dim, int spacedim = dim>
using CellAccessor< dim, spacedim >::AccessorData = typename TriaAccessor<dim, dim, spacedim>::AccessorData

Propagate the AccessorData type into the present class.

Definition at line 3077 of file tria_accessor.h.

◆ Container

template<int dim, int spacedim = dim>
using CellAccessor< dim, spacedim >::Container = Triangulation<dim, spacedim>

Define the type of the container this is part of.

Definition at line 3082 of file tria_accessor.h.

◆ LocalData

using TriaAccessorBase< structdim, dim, spacedim >::LocalData = void *
inherited

Data type to be used for passing parameters from iterators to the accessor classes in a unified way, no matter what the type of number of these parameters is.

Definition at line 448 of file tria_accessor.h.

Constructor & Destructor Documentation

◆ CellAccessor() [1/6]

template<int dim, int spacedim = dim>
CellAccessor< dim, spacedim >::CellAccessor ( const Triangulation< dim, spacedim > *  parent = nullptr,
const int  level = -1,
const int  index = -1,
const AccessorData local_data = nullptr 
)

Constructor.

◆ CellAccessor() [2/6]

template<int dim, int spacedim = dim>
CellAccessor< dim, spacedim >::CellAccessor ( const TriaAccessor< dim, dim, spacedim > &  cell_accessor)

Copy constructor.

◆ CellAccessor() [3/6]

template<int dim, int spacedim>
template<int structdim2, int dim2, int spacedim2>
CellAccessor< dim, spacedim >::CellAccessor ( const InvalidAccessor< structdim2, dim2, spacedim2 > &  )

Conversion constructor. This constructor exists to make certain constructs simpler to write in dimension independent code. For example, it allows assigning a face iterator to a line iterator, an operation that is useful in 2d but doesn't make any sense in 3d. The constructor here exists for the purpose of making the code conform to C++ but it will unconditionally abort; in other words, assigning a face iterator to a line iterator is better put into an if-statement that checks that the dimension is two, and assign to a quad iterator in 3d (an operator that, without this constructor would be illegal if we happen to compile for 2d).

Definition at line 4235 of file tria_accessor.h.

◆ CellAccessor() [4/6]

template<int dim, int spacedim>
template<int structdim2, int dim2, int spacedim2>
CellAccessor< dim, spacedim >::CellAccessor ( const TriaAccessor< structdim2, dim2, spacedim2 > &  )

Another conversion operator between objects that don't make sense, just like the previous one.

Definition at line 4265 of file tria_accessor.h.

◆ CellAccessor() [5/6]

template<int dim, int spacedim = dim>
CellAccessor< dim, spacedim >::CellAccessor ( const CellAccessor< dim, spacedim > &  )
default

Copy constructor.

◆ CellAccessor() [6/6]

template<int dim, int spacedim = dim>
CellAccessor< dim, spacedim >::CellAccessor ( CellAccessor< dim, spacedim > &&  )
default

Move constructor.

◆ ~CellAccessor()

template<int dim, int spacedim = dim>
CellAccessor< dim, spacedim >::~CellAccessor ( )
default

Destructor.

Member Function Documentation

◆ operator=() [1/2]

template<int dim, int spacedim = dim>
CellAccessor< dim, spacedim > & CellAccessor< dim, spacedim >::operator= ( const CellAccessor< dim, spacedim > &  )
delete

Copy operator. These operators are usually used in a context like iterator a,b; *a=*b;. Presumably, the intent here is to copy the object pointed to by b to the object pointed to by a. However, the result of dereferencing an iterator is not an object but an accessor; consequently, this operation is not useful for iterators on triangulations. Consequently, this operator is declared as deleted and can not be used.

◆ operator=() [2/2]

template<int dim, int spacedim = dim>
CellAccessor< dim, spacedim > & CellAccessor< dim, spacedim >::operator= ( CellAccessor< dim, spacedim > &&  )
default

Move assignment operator.

◆ as_dof_handler_iterator()

template<int dim, int spacedim>
TriaActiveIterator< DoFCellAccessor< dim, spacedim, false > > CellAccessor< dim, spacedim >::as_dof_handler_iterator ( const DoFHandler< dim, spacedim > &  dof_handler) const

A function that converts a Triangulation active cell iterator to a DoFHandler active cell iterator, or a DoFHandler active cell iterator to an active cell iterator of another DoFHandler. The iterator must be associated with the triangulation of the dof_handler.

Parameters
dof_handlerThe DoFHandler for the output active cell iterator.
Returns
An active cell iterator for the dof_handler, matching the cell referenced by the input iterator. The type of the returned object is a DoFHandler::active_cell_iterator.

Definition at line 2057 of file tria_accessor.cc.

◆ as_dof_handler_level_iterator()

template<int dim, int spacedim>
TriaIterator< DoFCellAccessor< dim, spacedim, true > > CellAccessor< dim, spacedim >::as_dof_handler_level_iterator ( const DoFHandler< dim, spacedim > &  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 iterator must be associated with the triangulation of the dof_handler.

Definition at line 2078 of file tria_accessor.cc.

◆ child()

template<int dim, int spacedim = dim>
TriaIterator< CellAccessor< dim, spacedim > > CellAccessor< dim, spacedim >::child ( const unsigned int  i) const

Return a pointer to the ith child. Overloaded version which returns a more reasonable iterator class.

◆ child_iterators()

template<int dim, int spacedim = dim>
boost::container::small_vector< TriaIterator< CellAccessor< dim, spacedim > >, GeometryInfo< dim >::max_children_per_cell > CellAccessor< dim, spacedim >::child_iterators ( ) const

Return an array of iterators to all children of this cell.

◆ face()

template<int dim, int spacedim = dim>
TriaIterator< TriaAccessor< dim - 1, dim, spacedim > > CellAccessor< dim, spacedim >::face ( const unsigned int  i) const

Return an iterator to the ith face of this cell.

◆ face_iterator_to_index()

template<int dim, int spacedim = dim>
unsigned int CellAccessor< dim, spacedim >::face_iterator_to_index ( const TriaIterator< TriaAccessor< dim - 1, dim, spacedim > > &  face) const

Return the face number of face on the current cell. This is the inverse function of TriaAccessor::face().

◆ face_iterators()

template<int dim, int spacedim = dim>
boost::container::small_vector< TriaIterator< TriaAccessor< dim - 1, dim, spacedim > >, GeometryInfo< dim >::faces_per_cell > CellAccessor< dim, spacedim >::face_iterators ( ) const

Return an array of iterators to all faces of this cell.

◆ face_index()

template<int dim, int spacedim = dim>
unsigned int CellAccessor< dim, spacedim >::face_index ( const unsigned int  i) const

Return the (global) index of the ith face of this cell.

Note
Despite the name, the index returned here is only global in the sense that it is specific to a particular Triangulation object or, in the case the triangulation is actually of type parallel::distributed::Triangulation, specific to that part of the distributed triangulation stored on the current processor.

◆ neighbor_child_on_subface()

template<int dim, int spacedim>
TriaIterator< CellAccessor< dim, spacedim > > CellAccessor< dim, spacedim >::neighbor_child_on_subface ( const unsigned int  face_no,
const unsigned int  subface_no 
) const

Return an iterator to that cell that neighbors the present cell on the given face and subface number.

To succeed, the present cell must not be further refined, and the neighbor on the given face must be further refined exactly once; the returned cell is then a child of that neighbor.

The function may not be called in 1d, since there we have no subfaces. The implementation of this function is rather straightforward in 2d, by first determining which face of the neighbor cell the present cell is bordering on (this is what the neighbor_of_neighbor function does), and then asking GeometryInfo::child_cell_on_subface for the index of the child.

However, the situation is more complicated in 3d, since there faces may have more than one orientation, and we have to use face_orientation, face_flip and face_rotation for both this and the neighbor cell to figure out which cell we want to have.

This can lead to surprising results: if we are sitting on a cell and are asking for a cell behind subface sf, then this means that we are considering the subface for the face in the natural direction for the present cell. However, if the face as seen from this cell has face_orientation()==false, then the child of the face that separates the present cell from the neighboring cell's child is not necessarily the sf-th child of the face of this cell. This is so because the subface_no on a cell corresponds to the subface with respect to the intrinsic ordering of the present cell, whereas children of face iterators are computed with respect to the intrinsic ordering of faces; these two orderings are only identical if the face orientation is true, and reversed otherwise.

Similarly, effects of face_flip()==true and face_rotation()==true(), both of which indicate a non-standard face have to be considered.

Fortunately, this is only very rarely of concern, since usually one simply wishes to loop over all finer neighbors at a given face of an active cell. Only in the process of refinement of a Triangulation we want to set neighbor information for both our child cells and the neighbor's children. Since we can respect orientation of faces from our current cell in that case, we do NOT respect face_orientation, face_flip and face_rotation of the present cell within this function, i.e. the returned neighbor's child is behind subface subface concerning the intrinsic ordering of the given face.

Definition at line 3026 of file tria_accessor.cc.

◆ neighbor()

template<int dim, int spacedim = dim>
TriaIterator< CellAccessor< dim, spacedim > > CellAccessor< dim, spacedim >::neighbor ( const unsigned int  face_no) const

Return an iterator to the neighboring cell on the other side of the face with number face_no. If the neighbor does not exist, i.e., if the face with number face_no of the current object is at the boundary, then an invalid iterator is returned.

Consequently, the index face_no must be less than n_faces().

The neighbor of a cell has at most the same level as this cell. For example, consider the following situation:

Here, if you are on the top right cell and you ask for its left neighbor (which is, according to the conventions spelled out in the GeometryInfo class, its zeroth neighbor), then you will get the mother cell of the four small cells at the top left. In other words, the cell you get as neighbor has the same refinement level as the one you're on right now (the top right one) and it may have children.

On the other hand, if you were at the top right cell of the four small cells at the top left, and you asked for the right neighbor (which is associated with index face_no=1), then you would get the large cell at the top right which in this case has a lower refinement level and no children of its own.

◆ neighbor_index()

template<int dim, int spacedim = dim>
int CellAccessor< dim, spacedim >::neighbor_index ( const unsigned int  face_no) const

Return the cell index of the neighboring cell on the other side of the face with index face_no. If the neighbor does not exist, this function returns -1.

This function is equivalent to cell->neighbor(face_no)->index(). See neighbor() for more details.

◆ neighbor_level()

template<int dim, int spacedim = dim>
int CellAccessor< dim, spacedim >::neighbor_level ( const unsigned int  face_no) const

Return the level of the neighboring cell on the other side of the face with number face_no. If the neighbor does not exist, this function returns -1.

This function is equivalent to cell->neighbor(face_no)->level(). See neighbor() for more details.

◆ neighbor_of_neighbor()

template<int dim, int spacedim>
unsigned int CellAccessor< dim, spacedim >::neighbor_of_neighbor ( const unsigned int  face_no) const

Return the how-many'th neighbor this cell is of cell->neighbor(face_no), i.e. return other_face_no such that cell->neighbor(face_no)->neighbor(other_face_no)==cell. This function is the right one if you want to know how to get back from a neighbor to the present cell.

Note that this operation is only useful if the neighbor is not coarser than the present cell. If the neighbor is coarser this function throws an exception. Use the neighbor_of_coarser_neighbor function in that case.

Definition at line 2494 of file tria_accessor.cc.

◆ neighbor_is_coarser()

template<int dim, int spacedim>
bool CellAccessor< dim, spacedim >::neighbor_is_coarser ( const unsigned int  face_no) const

Return, whether the neighbor is coarser then the present cell. This is important in case of anisotropic refinement where this information does not depend on the levels of the cells.

Note, that in an anisotropic setting, a cell can only be coarser than another one at a given face, not on a general basis. The face of the finer cell is contained in the corresponding face of the coarser cell, the finer face is either a child or a grandchild of the coarser face.

Definition at line 2508 of file tria_accessor.cc.

◆ neighbor_of_coarser_neighbor()

template<int dim, int spacedim>
std::pair< unsigned int, unsigned int > CellAccessor< dim, spacedim >::neighbor_of_coarser_neighbor ( const unsigned int  neighbor) const

This function is a generalization of the neighbor_of_neighbor function for the case of a coarser neighbor. It returns a pair of numbers, face_no and subface_no, with the following property, if the neighbor is not refined: cell->neighbor(neighbor)->neighbor_child_on_subface(face_no, subface_no)==cell. In 3d, a coarser neighbor can still be refined. In that case subface_no denotes the child index of the neighbors face that relates to our face: cell->neighbor(neighbor)->face(face_no)->child(subface_no)==cell->face(neighbor). This case in 3d and how it can happen is discussed in the introduction of the step-30 tutorial program.

This function is impossible for dim==1.

Definition at line 2519 of file tria_accessor.cc.

◆ neighbor_face_no()

template<int dim, int spacedim = dim>
unsigned int CellAccessor< dim, spacedim >::neighbor_face_no ( const unsigned int  neighbor) const

This function is a generalization of the neighbor_of_neighbor and the neighbor_of_coarser_neighbor functions. It checks whether the neighbor is coarser or not and calls the respective function. In both cases, only the face_no is returned.

◆ is_level_cell()

template<int dim, int spacedim = dim>
static bool CellAccessor< dim, spacedim >::is_level_cell ( )
static

Compatibility interface with DoFCellAccessor. Always returns false.

◆ has_periodic_neighbor()

template<int dim, int spacedim>
bool CellAccessor< dim, spacedim >::has_periodic_neighbor ( const unsigned int  i) const

If the cell has a periodic neighbor at its ith face, this function returns true, otherwise, the returned value is false.

Definition at line 2695 of file tria_accessor.cc.

◆ periodic_neighbor()

template<int dim, int spacedim>
TriaIterator< CellAccessor< dim, spacedim > > CellAccessor< dim, spacedim >::periodic_neighbor ( const unsigned int  i) const

For a cell with its ith face at a periodic boundary, see the entry for periodic boundaries, this function returns an iterator to the cell on the other side of the periodic boundary. If there is no periodic boundary at the ith face, an exception will be thrown. In order to avoid running into an exception, check the result of has_periodic_neighbor() for the ith face prior to using this function. The behavior of periodic_neighbor() is similar to neighbor(), in the sense that the returned cell has at most the same level of refinement as the current cell. On distributed meshes, by calling Triangulation::add_periodicity(), we can make sure that the element on the other side of the periodic boundary exists in this rank as a ghost cell or a locally owned cell.

Definition at line 2723 of file tria_accessor.cc.

◆ neighbor_or_periodic_neighbor()

template<int dim, int spacedim>
TriaIterator< CellAccessor< dim, spacedim > > CellAccessor< dim, spacedim >::neighbor_or_periodic_neighbor ( const unsigned int  i) const

For a cell whose ith face is not at a boundary, this function returns the same result as neighbor(). If the ith face is at a periodic boundary this function returns the same result as periodic_neighbor(). If neither of the aforementioned conditions are met, i.e. the ith face is on a nonperiodic boundary, an exception will be thrown.

Definition at line 2750 of file tria_accessor.cc.

◆ periodic_neighbor_child_on_subface()

template<int dim, int spacedim>
TriaIterator< CellAccessor< dim, spacedim > > CellAccessor< dim, spacedim >::periodic_neighbor_child_on_subface ( const unsigned int  face_no,
const unsigned int  subface_no 
) const

Return an iterator to the periodic neighbor of the cell at a given face and subface number. The general guidelines for using this function is similar to the function neighbor_child_on_subface(). The implementation of this function is consistent with periodic_neighbor_of_coarser_periodic_neighbor(). For instance, assume that we are sitting on a cell named cell1, whose neighbor behind the ith face is one level coarser. Let us name this coarser neighbor cell2. Then, by calling periodic_neighbor_of_coarser_periodic_neighbor(), from cell1, we get a face_num and a subface_num. Now, if we call periodic_neighbor_child_on_subface() from cell2, with the above face_num and subface_num, we get an iterator to cell1.

Definition at line 2767 of file tria_accessor.cc.

◆ periodic_neighbor_of_coarser_periodic_neighbor()

template<int dim, int spacedim = dim>
std::pair< unsigned int, unsigned int > CellAccessor< dim, spacedim >::periodic_neighbor_of_coarser_periodic_neighbor ( const unsigned  face_no) const

This function is a generalization of periodic_neighbor_of_periodic_neighbor() for those cells which have a coarser periodic neighbor. The returned pair of numbers can be used in periodic_neighbor_child_on_subface() to get back to the current cell. In other words, the following assertion should be true, for a cell with coarser periodic neighbor: cell->periodic_neighbor(i)->periodic_neighbor_child_on_subface(face_no, subface_no)==cell

Definition at line 2819 of file tria_accessor.cc.

◆ periodic_neighbor_index()

template<int dim, int spacedim>
int CellAccessor< dim, spacedim >::periodic_neighbor_index ( const unsigned int  i) const

This function returns the index of the periodic neighbor at the ith face of the current cell. If there is no periodic neighbor at the given face, the returned value is -1.

Definition at line 2880 of file tria_accessor.cc.

◆ periodic_neighbor_level()

template<int dim, int spacedim>
int CellAccessor< dim, spacedim >::periodic_neighbor_level ( const unsigned int  i) const

This function returns the level of the periodic neighbor at the ith face of the current cell. If there is no periodic neighbor at the given face, the returned value is -1.

Definition at line 2890 of file tria_accessor.cc.

◆ periodic_neighbor_of_periodic_neighbor()

template<int dim, int spacedim>
unsigned int CellAccessor< dim, spacedim >::periodic_neighbor_of_periodic_neighbor ( const unsigned int  i) const

For a cell with a periodic neighbor at its ith face, this function returns the face number of that periodic neighbor such that the current cell is the periodic neighbor of that neighbor. In other words the following assertion holds for those cells which have a periodic neighbor with the same or a higher level of refinement as the current cell: {cell->periodic_neighbor(i)-> periodic_neighbor(cell->periodic_neighbor_of_periodic_neighbor(i))==cell} For the cells with a coarser periodic neighbor, one should use periodic_neighbor_of_coarser_periodic_neighbor() and periodic_neighbor_child_on_subface() to get back to the current cell.

Definition at line 2900 of file tria_accessor.cc.

◆ periodic_neighbor_face_no()

template<int dim, int spacedim>
unsigned int CellAccessor< dim, spacedim >::periodic_neighbor_face_no ( const unsigned int  i) const

If a cell has a periodic neighbor at its ith face, this function returns the face number of the periodic neighbor, which is connected to the ith face of this cell.

Definition at line 2910 of file tria_accessor.cc.

◆ periodic_neighbor_is_coarser()

template<int dim, int spacedim>
bool CellAccessor< dim, spacedim >::periodic_neighbor_is_coarser ( const unsigned int  i) const

This function returns true if the element on the other side of the periodic boundary is coarser and returns false otherwise. The implementation allows this function to work in the case of anisotropic refinement.

Definition at line 2941 of file tria_accessor.cc.

◆ at_boundary() [1/2]

template<int dim, int spacedim>
bool CellAccessor< dim, spacedim >::at_boundary ( const unsigned int  i) const

Return whether the ith vertex or face (depending on the dimension) is part of the boundary. This is true, if the ith neighbor does not exist.

Definition at line 2996 of file tria_accessor.cc.

◆ at_boundary() [2/2]

template<int dim, int spacedim>
bool CellAccessor< dim, spacedim >::at_boundary ( ) const

Return whether the cell is at the boundary. Being at the boundary is defined by one face being on the boundary. Note that this does not catch cases where only one vertex of a 2d or 3d subobject is at the boundary, or where only one line of a hex is at the boundary while the interiors of all faces are in the interior of the domain. For the latter case, the has_boundary_lines function is the right one to ask.

Definition at line 2146 of file tria_accessor.cc.

◆ has_boundary_lines()

template<int dim, int spacedim>
bool CellAccessor< dim, spacedim >::has_boundary_lines ( ) const

This is a slight variation to the at_boundary function: for 1 and 2 dimensions, it is equivalent, for three dimensions it returns whether at least one of the 12 lines of the hexahedron is at a boundary. This, of course, includes the case where a whole face is at the boundary, but also some other cases.

Definition at line 3008 of file tria_accessor.cc.

◆ refine_flag_set()

template<int dim, int spacedim = dim>
RefinementCase< dim > CellAccessor< dim, spacedim >::refine_flag_set ( ) const

Return the RefinementCase<dim> this cell was flagged to be refined with. The return value of this function can be compared to a bool to check if this cell is flagged for any kind of refinement. For example, if you have previously called cell->set_refine_flag() for a cell, then you will enter the 'if' block in the following snippet:

if (cell->refine_flag_set())
{
// yes, this cell is marked for refinement.
}

◆ set_refine_flag()

template<int dim, int spacedim = dim>
void CellAccessor< dim, spacedim >::set_refine_flag ( const RefinementCase< dim >  ref_case = RefinementCase< dim >::isotropic_refinement) const

Flag the cell pointed to for refinement. This function is only allowed for active cells. Keeping the default value for ref_case will mark this cell for isotropic refinement.

If you choose anisotropic refinement, for example by passing as argument one of the flags RefinementCase::cut_x, RefinementCase::cut_y, RefinementCase::cut_z, or a combination of these, then keep in mind that refining in x-, y-, or z-direction happens with regard to the local coordinate system of the cell. In other words, these flags determine which edges and faces of the cell will be cut into new edges and faces. On the other hand, this process is independent of how the cell is oriented within the global coordinate system, and you should not assume any particular orientation of the cell's local coordinate system within the global coordinate system of the space it lives in.

◆ clear_refine_flag()

template<int dim, int spacedim = dim>
void CellAccessor< dim, spacedim >::clear_refine_flag ( ) const

Clear the refinement flag.

◆ flag_for_face_refinement()

template<int dim, int spacedim = dim>
bool CellAccessor< dim, spacedim >::flag_for_face_refinement ( const unsigned int  face_no,
const RefinementCase< dim - 1 > &  face_refinement_case = RefinementCase< dim - 1 >::isotropic_refinement 
) const

Modify the refinement flag of the cell to ensure (at least) the given refinement case face_refinement_case at face face_no, taking into account orientation, flip and rotation of the face. Return, whether the refinement flag had to be modified. This function is only allowed for active cells.

◆ flag_for_line_refinement()

template<int dim, int spacedim = dim>
bool CellAccessor< dim, spacedim >::flag_for_line_refinement ( const unsigned int  line_no) const

Modify the refinement flag of the cell to ensure that line face_no will be refined. Return, whether the refinement flag had to be modified. This function is only allowed for active cells.

◆ subface_case()

template<int dim, int spacedim = dim>
::internal::SubfaceCase< dim > CellAccessor< dim, spacedim >::subface_case ( const unsigned int  face_no) const

Return the SubfaceCase of face face_no. Note that this is not identical to asking cell->face(face_no)->refinement_case() since the latter returns a RefinementCase<dim-1> and thus only considers one (anisotropic) refinement, whereas this function considers the complete refinement situation including possible refinement of the face's children. This function may only be called for active cells in 2d and 3d.

◆ coarsen_flag_set()

template<int dim, int spacedim = dim>
bool CellAccessor< dim, spacedim >::coarsen_flag_set ( ) const

Return whether the coarsen flag is set or not.

◆ set_coarsen_flag()

template<int dim, int spacedim = dim>
void CellAccessor< dim, spacedim >::set_coarsen_flag ( ) const

Flag the cell pointed to for coarsening. This function is only allowed for active cells.

◆ clear_coarsen_flag()

template<int dim, int spacedim = dim>
void CellAccessor< dim, spacedim >::clear_coarsen_flag ( ) const

Clear the coarsen flag.

◆ material_id()

template<int dim, int spacedim>
types::material_id CellAccessor< dim, spacedim >::material_id ( ) const

Return the material id of this cell.

For a typical use of this function, see the step-28 tutorial program.

See the glossary for more information.

Definition at line 2159 of file tria_accessor.cc.

◆ set_material_id()

template<int dim, int spacedim>
void CellAccessor< dim, spacedim >::set_material_id ( const types::material_id  new_material_id) const

Set the material id of this cell.

For a typical use of this function, see the step-28 tutorial program.

See the glossary for more information.

Definition at line 2171 of file tria_accessor.cc.

◆ recursively_set_material_id()

template<int dim, int spacedim>
void CellAccessor< dim, spacedim >::recursively_set_material_id ( const types::material_id  new_material_id) const

Set the material id of this cell and all its children (and grand-children, and so on) to the given value.

See the glossary for more information.

Definition at line 2185 of file tria_accessor.cc.

◆ subdomain_id()

template<int dim, int spacedim = dim>
types::subdomain_id CellAccessor< dim, spacedim >::subdomain_id ( ) const

Return the subdomain id of this cell.

See the glossary for more information.

Note
The subdomain of a cell is a property only defined for active cells, i.e., cells that are not further refined. Consequently, you can only call this function if the cell it refers to has no children. For multigrid methods in parallel, it is also important to know which processor owns non-active cells, and for this you can call level_subdomain_id().

◆ set_subdomain_id()

template<int dim, int spacedim>
void CellAccessor< dim, spacedim >::set_subdomain_id ( const types::subdomain_id  new_subdomain_id) const

Set the subdomain id of this cell.

See the glossary for more information. This function should not be called if you use a parallel::distributed::Triangulation object.

Note
The subdomain of a cell is a property only defined for active cells, i.e., cells that are not further refined. Consequently, you can only call this function if the cell it refers to has no children. For multigrid methods in parallel, it is also important to know which processor owns non-active cells, and for this you can call level_subdomain_id().

Definition at line 2199 of file tria_accessor.cc.

◆ level_subdomain_id()

template<int dim, int spacedim = dim>
types::subdomain_id CellAccessor< dim, spacedim >::level_subdomain_id ( ) const

Get the level subdomain id of this cell. This is used for parallel multigrid where not only the global mesh (consisting of the active cells) is partitioned among processors, but also the individual levels of the hierarchy of recursively refined cells that make up the mesh. In other words, the level subdomain id is a property that is also defined for non-active cells if a multigrid hierarchy is used.

◆ set_level_subdomain_id()

template<int dim, int spacedim>
void CellAccessor< dim, spacedim >::set_level_subdomain_id ( const types::subdomain_id  new_level_subdomain_id) const

Set the level subdomain id of this cell. This is used for parallel multigrid.

Definition at line 2213 of file tria_accessor.cc.

◆ recursively_set_subdomain_id()

template<int dim, int spacedim>
void CellAccessor< dim, spacedim >::recursively_set_subdomain_id ( const types::subdomain_id  new_subdomain_id) const

Set the subdomain id of this cell (if it is active) or all its terminal children (and grand-children, and so on, as long as they have no children of their own) to the given value. Since the subdomain id is a concept that is only defined for cells that are active (i.e., have no children of their own), this function only sets the subdomain ids for all children and grand children of this cell that are actually active, skipping intermediate child cells.

See the glossary for more information. This function should not be called if you use a parallel::distributed::Triangulation object since there the subdomain id is implicitly defined by which processor you're on.

Definition at line 2343 of file tria_accessor.cc.

◆ global_active_cell_index()

template<int dim, int spacedim = dim>
types::global_cell_index CellAccessor< dim, spacedim >::global_active_cell_index ( ) const

Return a globally unique cell index for the current cell, assuming it is not artificial. The value is identical to active_cell_index() if the cell is part of a serial triangulation.

In the context of parallel triangulations, locally-owned cells are enumerated contiguously within each subdomain of the mesh. This ensures that the index returned by this function can be used as the index into vectors with a total of Triangulation::n_globally_active_cells() entries, and for which every process stores a contiguous part. If such a cell-data vector has been set up with parallel::TriangulationBase::global_active_cell_index_partitioner(), the index returned by this function can then be used to access the correct vector entry.

◆ global_level_cell_index()

template<int dim, int spacedim = dim>
types::global_cell_index CellAccessor< dim, spacedim >::global_level_cell_index ( ) const

Return a globally unique index for a non-artificial level cell.

Note
Similar to global_active_cell_index(), with the difference that the cell-data vector has been set up with parallel::TriangulationBase::global_level_cell_index_partitioner().

◆ direction_flag()

template<int dim, int spacedim>
bool CellAccessor< dim, spacedim >::direction_flag ( ) const

Return the orientation of this cell. This function always returns true if dim==spacedim. It can return true or false if dim==spacedim-1. The function cannot be called (and will abort with an error) if called for dim<spacedim-1.

For the meaning of this flag, see GlossDirectionFlag.

Definition at line 2224 of file tria_accessor.cc.

◆ active_cell_index()

template<int dim, int spacedim = dim>
unsigned int CellAccessor< dim, spacedim >::active_cell_index ( ) const

Return the how many-th active cell the current cell is (assuming the current cell is indeed active). This is useful, for example, if you are accessing the elements of a vector with as many entries as there are active cells. Such vectors are used for estimating the error on each cell of a triangulation, for specifying refinement criteria passed to the functions in GridRefinement, and for generating cell-wise output.

The function throws an exception if the current cell is not active.

Note
If the triangulation this function is called on is of type parallel::distributed::Triangulation, then active cells may be locally owned, ghost cells, or artificial (see GlossLocallyOwnedCell, GlossGhostCell, and GlossArtificialCell). This function counts over all of them, including ghost and artificial active cells. This implies that the index returned by this function uniquely identifies a cell within the triangulation on a single processor, but does not uniquely identify the cell among the (parts of the) triangulation that is shared among processors. If you would like to identify active cells across processors, you need to consider the CellId of a cell returned by CellAccessor::id().

◆ parent_index()

template<int dim, int spacedim>
int CellAccessor< dim, spacedim >::parent_index ( ) const

Return the index of the parent of this cell within the level of the triangulation to which the parent cell belongs. The level of the parent is of course one lower than that of the present cell. If the parent does not exist (i.e., if the object is at the coarsest level of the mesh hierarchy), an exception is generated.

Definition at line 2281 of file tria_accessor.cc.

◆ parent()

template<int dim, int spacedim>
TriaIterator< CellAccessor< dim, spacedim > > CellAccessor< dim, spacedim >::parent ( ) const

Return an iterator to the parent. If the parent does not exist (i.e., if the object is at the coarsest level of the mesh hierarchy), an exception is generated.

Definition at line 2329 of file tria_accessor.cc.

◆ is_active()

template<int dim, int spacedim = dim>
bool CellAccessor< dim, spacedim >::is_active ( ) const

Test that the cell has no children (this is the criterion for whether a cell is called "active").

See the glossary for more information.

◆ is_locally_owned()

template<int dim, int spacedim = dim>
bool CellAccessor< dim, spacedim >::is_locally_owned ( ) const

Return whether this cell is owned by the current processor or is owned by another processor. The function always returns true if applied to an object of type Triangulation, but may yield false if the triangulation is of type parallel::distributed::Triangulation.

See the glossary and the Parallel computing with multiple processors using distributed memory module for more information.

Postcondition
The returned value is equal to !is_ghost() && !is_artificial().
Note
Whether a cell is a ghost cell, artificial, or is locally owned or is a property that only pertains to cells that are active. Consequently, you can only call this function if the cell it refers to has no children.

◆ is_locally_owned_on_level()

template<int dim, int spacedim = dim>
bool CellAccessor< dim, spacedim >::is_locally_owned_on_level ( ) const

Return true if either the Triangulation is not distributed or if level_subdomain_id() is equal to the id of the current processor.

◆ is_ghost()

template<int dim, int spacedim = dim>
bool CellAccessor< dim, spacedim >::is_ghost ( ) const

Return true if:

  1. This cell exists in the global mesh (i.e., it is not artificial), and
  2. This cell is owned by another processor (i.e., has a subdomain_id different from Triangulation::locally_owned_subdomain())

In all other cases the returned value is false. In particular, only parallel Triangulations (i.e., Triangulations inheriting from parallel::TriangulationBase) can have ghost cells, so for a serial Triangulation the returned value is false.

See the glossary and the Parallel computing with multiple processors using distributed memory module for more information.

Postcondition
The returned value is equal to !is_locally_owned() && !is_artificial().
Note
For parallel::distributed::Triangulation and parallel::fullydistributed::Triangulation, ghost cells are always adjacent to locally owned cells. For parallel::shared::Triangulation they may not be, dependent on whether or not the triangulation uses artificial cells - see parallel::shared::Triangulation::Triangulation() for more information.
Whether a cell is a ghost cell, artificial, or is locally owned is a property that only pertains to cells that are active. Consequently, you can only call this function if the cell it refers to has no children.

◆ is_ghost_on_level()

template<int dim, int spacedim = dim>
bool CellAccessor< dim, spacedim >::is_ghost_on_level ( ) const

Return true if either the Triangulation is not distributed or if the cell is not artificial and the level_subdomain_id() is not equal to the id of the current processor.

◆ is_artificial()

template<int dim, int spacedim = dim>
bool CellAccessor< dim, spacedim >::is_artificial ( ) const

Return whether this cell is artificial, i.e. it isn't one of the cells owned by the current processor, and it also doesn't border on one. As a consequence, it exists in the mesh to ensure that each processor has all coarse mesh cells and that the 2:1 ratio of neighboring cells is maintained, but it is not one of the cells we should work on on the current processor. In particular, there is no guarantee that this cell isn't, in fact, further refined on one of the other processors.

This function only makes sense if the triangulation used is of kind parallel::distributed::Triangulation. In all other cases, the returned value is always false.

See the glossary and the Parallel computing with multiple processors using distributed memory module for more information.

Postcondition
The returned value is equal to !is_ghost() && !is_locally_owned().
Note
Whether a cell is a ghost cell, artificial, or is locally owned is a property that only pertains to cells that are active. Consequently, you can only call this function if the cell it refers to has no children.

◆ is_artificial_on_level()

template<int dim, int spacedim = dim>
bool CellAccessor< dim, spacedim >::is_artificial_on_level ( ) const

Similar to is_artificial() but checking the conditions on the levels.

Postcondition
The returned value is equal to !is_ghost_on_level() && !is_locally_owned_on_level().

◆ point_inside() [1/7]

template<int dim, int spacedim = dim>
bool CellAccessor< dim, spacedim >::point_inside ( const Point< spacedim > &  p) const

Test whether the point p is inside this cell. Points on the boundary are counted as being inside the cell.

Note that this function assumes that the mapping between unit cell and real cell is (bi-, tri-)linear, i.e. that faces in 2d and edges in 3d are straight lines. If you have higher order transformations, results may be different as to whether a point is in- or outside the cell in real space.

In case of codim>0, the point is first projected to the manifold where the cell is embedded and then check if this projection is inside the cell.

◆ set_neighbor()

template<int dim, int spacedim>
void CellAccessor< dim, spacedim >::set_neighbor ( const unsigned int  i,
const TriaIterator< CellAccessor< dim, spacedim > > &  pointer 
) const

Set the neighbor i of this cell to the cell pointed to by pointer.

This function shouldn't really be public (but needs to for various reasons in order not to make a long list of functions friends): it modifies internal data structures and may leave things. Do not use it from application codes.

Definition at line 2357 of file tria_accessor.cc.

◆ id()

template<int dim, int spacedim>
CellId CellAccessor< dim, spacedim >::id ( ) const

Return a unique ID for the current cell. This ID is constructed from the path in the hierarchy from the coarse father cell and works correctly in parallel computations using objects of type parallel::distributed::Triangulation. This function is therefore useful in providing a unique identifier for cells (active or not) that also works for parallel triangulations. See the documentation of the CellId class for more information.

Note
This operation takes O(level) time to compute. In most practical cases, the number of levels of a triangulation will depend logarithmically on the number of cells in the triangulation.

Definition at line 2387 of file tria_accessor.cc.

◆ diameter() [1/2]

template<int dim, int spacedim = dim>
double CellAccessor< dim, spacedim >::diameter ( const Mapping< dim, spacedim > &  mapping) const

The same as TriaAccessor::diameter() but also taking a Mapping class.

◆ neighbor_of_neighbor_internal()

template<int dim, int spacedim>
unsigned int CellAccessor< dim, spacedim >::neighbor_of_neighbor_internal ( const unsigned int  neighbor) const
protected

This function assumes that the neighbor is not coarser than the current cell. In this case it returns the neighbor_of_neighbor() value. If, however, the neighbor is coarser this function returns an invalid_unsigned_int.

This function is not for public use. Use the function neighbor_of_neighbor() instead which throws an exception if called for a coarser neighbor. If neighbor is indeed coarser (you get to know this by e.g. the neighbor_is_coarser() function) then the neighbor_of_coarser_neighbor() function should be call. If you'd like to know only the face_no which is required to get back from the neighbor to the present cell then simply use the neighbor_face_no() function which can be used for coarser as well as non-coarser neighbors.

Definition at line 2428 of file tria_accessor.cc.

◆ point_inside_codim()

template<int dim, int spacedim>
template<int dim_, int spacedim_>
bool CellAccessor< dim, spacedim >::point_inside_codim ( const Point< spacedim_ > &  p) const
protected

As for any codim>0 we can use a similar code and c++ does not allow partial templates. we use this auxiliary function that is then called from point_inside.

Definition at line 2101 of file tria_accessor.cc.

◆ set_active_cell_index()

template<int dim, int spacedim>
void CellAccessor< dim, spacedim >::set_active_cell_index ( const unsigned int  active_cell_index) const
private

Set the active cell index of a cell. This is done at the end of refinement.

Definition at line 2296 of file tria_accessor.cc.

◆ set_global_active_cell_index()

template<int dim, int spacedim>
void CellAccessor< dim, spacedim >::set_global_active_cell_index ( const types::global_cell_index  index) const
private

Set global active cell index for a cell.

Definition at line 2307 of file tria_accessor.cc.

◆ set_global_level_cell_index()

template<int dim, int spacedim>
void CellAccessor< dim, spacedim >::set_global_level_cell_index ( const types::global_cell_index  index) const
private

Set global level cell index for a level cell.

Definition at line 2318 of file tria_accessor.cc.

◆ set_parent()

template<int dim, int spacedim>
void CellAccessor< dim, spacedim >::set_parent ( const unsigned int  parent_index)
private

Set the parent of a cell.

Definition at line 2269 of file tria_accessor.cc.

◆ set_direction_flag()

template<int dim, int spacedim>
void CellAccessor< dim, spacedim >::set_direction_flag ( const bool  new_direction_flag) const
private

Set the orientation of this cell. This function can only be called if the argument is true if dim==spacedim. It can be called with either true or false if dim==spacedim-1. The function cannot be called (and will abort with an error) if called for dim<spacedim-1.

For the meaning of this flag, see GlossDirectionFlag.

Definition at line 2244 of file tria_accessor.cc.

◆ point_inside() [2/7]

bool CellAccessor< 1 >::point_inside ( const Point< 1 > &  p) const

Definition at line 1920 of file tria_accessor.cc.

◆ point_inside() [3/7]

bool CellAccessor< 2 >::point_inside ( const Point< 2 > &  p) const

Definition at line 1933 of file tria_accessor.cc.

◆ point_inside() [4/7]

bool CellAccessor< 3 >::point_inside ( const Point< 3 > &  p) const

Definition at line 1998 of file tria_accessor.cc.

◆ point_inside() [5/7]

bool CellAccessor< 1, 2 >::point_inside ( const Point< 2 > &  p) const

Definition at line 2119 of file tria_accessor.cc.

◆ point_inside() [6/7]

bool CellAccessor< 1, 3 >::point_inside ( const Point< 3 > &  p) const

Definition at line 2127 of file tria_accessor.cc.

◆ point_inside() [7/7]

bool CellAccessor< 2, 3 >::point_inside ( const Point< 3 > &  p) const

Definition at line 2135 of file tria_accessor.cc.

◆ used()

template<int structdim, int dim, int spacedim>
bool TriaAccessor< structdim, dim, spacedim >::used ( ) const
inherited

Test for the element being used or not. The return value is true for all iterators that are either normal iterators or active iterators, only raw iterators can return false. Since raw iterators are only used in the interiors of the library, you will not usually need this function.

◆ vertex_iterator()

template<int structdim, int dim, int spacedim>
TriaIterator< TriaAccessor< 0, dim, spacedim > > TriaAccessor< structdim, dim, spacedim >::vertex_iterator ( const unsigned int  i) const
inherited

Pointer to the ith vertex bounding this object. Throw an exception if dim=1.

◆ vertex_index()

template<int structdim, int dim, int spacedim>
unsigned int TriaAccessor< structdim, dim, spacedim >::vertex_index ( const unsigned int  i) const
inherited

Return the global index of i-th vertex of the current object. The convention regarding the numbering of vertices is laid down in the documentation of the GeometryInfo class.

Note that the returned value is only the index of the geometrical vertex. It has nothing to do with possible degrees of freedom associated with it. For this, see the DoFAccessor::vertex_dof_index functions.

Note
Despite the name, the index returned here is only global in the sense that it is specific to a particular Triangulation object or, in the case the triangulation is actually of type parallel::distributed::Triangulation, specific to that part of the distributed triangulation stored on the current processor.

◆ vertex()

template<int structdim, int dim, int spacedim>
Point< spacedim > & TriaAccessor< structdim, dim, spacedim >::vertex ( const unsigned int  i) const
inherited

Return a reference to the ith vertex. The reference is not const, i.e., it is possible to call this function on the left hand side of an assignment, thereby moving the vertex of a cell within the triangulation. Of course, doing so requires that you ensure that the new location of the vertex remains useful – for example, avoiding inverted or otherwise distorted (see also this glossary entry).

Note
When a cell is refined, its children inherit the position of the vertex positions of those vertices they share with the mother cell (plus the locations of the new vertices on edges, faces, and cell interiors that are created for the new child cells). If the vertex of a cell is moved, this implies that its children will also use these new locations. On the other hand, imagine a 2d situation where you have one cell that is refined (with four children) and then you move the central vertex connecting all four children. If you coarsen these four children again to the mother cell, then the location of the moved vertex is lost and if, in a later step, you refine the mother cell again, the then again new vertex will be placed again at the same position as the first time around – i.e., not at the location you had previously moved it to.
The behavior described above is relevant if you have a parallel::distributed::Triangulation object. There, refining a mesh always involves a re-partitioning. In other words, vertices of locally owned cells (see this glossary entry) that you may have moved to a different location on one processor may be moved to a different processor upon mesh refinement (even if these particular cells were not refined) which will re-create their position based on the position of the coarse cells they previously had, not based on the position these vertices had on the processor that previously owned them. In other words, in parallel computations, you will probably have to move nodes explicitly after every mesh refinement because vertex positions may or may not be preserved across the re-partitioning that accompanies mesh refinement.

◆ line()

template<int structdim, int dim, int spacedim>
typename::internal::TriangulationImplementation::Iterators< dim, spacedim >::line_iterator TriaAccessor< structdim, dim, spacedim >::line ( const unsigned int  i) const
inherited

Pointer to the ith line bounding this object.

◆ line_index()

template<int structdim, int dim, int spacedim>
unsigned int TriaAccessor< structdim, dim, spacedim >::line_index ( const unsigned int  i) const
inherited

Line index of the ith line bounding this object.

Implemented only for structdim>1, otherwise an exception generated.

◆ quad()

template<int structdim, int dim, int spacedim>
typename::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator TriaAccessor< structdim, dim, spacedim >::quad ( const unsigned int  i) const
inherited

Pointer to the ith quad bounding this object.

◆ quad_index()

template<int structdim, int dim, int spacedim>
unsigned int TriaAccessor< structdim, dim, spacedim >::quad_index ( const unsigned int  i) const
inherited

Quad index of the ith quad bounding this object.

Implemented only for structdim>2, otherwise an exception generated.

◆ has_children()

template<int structdim, int dim, int spacedim>
bool TriaAccessor< structdim, dim, spacedim >::has_children ( ) const
inherited

Test whether the object has children.

◆ n_children()

template<int structdim, int dim, int spacedim>
unsigned int TriaAccessor< structdim, dim, spacedim >::n_children ( ) const
inherited

Return the number of immediate children of this object. The number of children of an unrefined cell is zero.

◆ n_active_descendants()

template<int structdim, int dim, int spacedim>
unsigned int TriaAccessor< structdim, dim, spacedim >::n_active_descendants ( ) const
inherited

Compute and return the number of active descendants of this objects. For example, if all of the eight children of a hex are further refined isotropically exactly once, the returned number will be 64, not 80.

If the present cell is not refined, one is returned.

If one considers a triangulation as a forest where the root of each tree are the coarse mesh cells and nodes have descendants (the children of a cell), then this function returns the number of terminal nodes in the sub-tree originating from the current object; consequently, if the current object is not further refined, the answer is one.

◆ max_refinement_depth()

template<int structdim, int dim, int spacedim>
unsigned int TriaAccessor< structdim, dim, spacedim >::max_refinement_depth ( ) const
inherited

Return the number of times that this object is refined. Note that not all its children are refined that often (which is why we prepend max_), the returned number is rather the maximum number of refinement in any branch of children of this object.

For example, if this object is refined, and one of its children is refined exactly one more time, then max_refinement_depth should return 2.

If this object is not refined (i.e. it is active), then the return value is zero.

◆ child_iterator_to_index()

template<int structdim, int dim, int spacedim>
unsigned int TriaAccessor< structdim, dim, spacedim >::child_iterator_to_index ( const TriaIterator< TriaAccessor< structdim, dim, spacedim > > &  child) const
inherited

Return the child number of child on the current cell. This is the inverse function of TriaAccessor::child().

◆ isotropic_child()

template<int structdim, int dim, int spacedim>
TriaIterator< TriaAccessor< structdim, dim, spacedim > > TriaAccessor< structdim, dim, spacedim >::isotropic_child ( const unsigned int  i) const
inherited

Return an iterator to that object that is identical to the ith child for isotropic refinement. If the current object is refined isotropically, then the returned object is the ith child. If the current object is refined anisotropically, the returned child may in fact be a grandchild of the object, or may not exist at all (in which case an exception is generated).

◆ refinement_case()

template<int structdim, int dim, int spacedim>
RefinementCase< structdim > TriaAccessor< structdim, dim, spacedim >::refinement_case ( ) const
inherited

Return the RefinementCase of this cell.

◆ child_index()

template<int structdim, int dim, int spacedim>
int TriaAccessor< structdim, dim, spacedim >::child_index ( const unsigned int  i) const
inherited

Index of the ith child. The level of the child is one higher than that of the present cell, if the children of a cell are accessed. The children of faces have no level. If the child does not exist, -1 is returned.

◆ isotropic_child_index()

template<int structdim, int dim, int spacedim>
int TriaAccessor< structdim, dim, spacedim >::isotropic_child_index ( const unsigned int  i) const
inherited

Index of the ith isotropic child. See the isotropic_child() function for a definition of this concept. If the child does not exist, -1 is returned.

◆ boundary_id()

template<int structdim, int dim, int spacedim>
types::boundary_id TriaAccessor< structdim, dim, spacedim >::boundary_id ( ) const
inherited

Return the boundary indicator of this object.

If the return value is the special value numbers::internal_face_boundary_id, then this object is in the interior of the domain.

See also
Glossary entry on boundary indicators

◆ set_boundary_id()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::set_boundary_id ( const types::boundary_id  ) const
inherited

Set the boundary indicator of the current object. The same applies as for the boundary_id() function.

This function only sets the boundary object of the current object itself, not the indicators of the ones that bound it. For example, in 3d, if this function is called on a face, then the boundary indicator of the 4 edges that bound the face remain unchanged. If you want to set the boundary indicators of face and edges at the same time, use the set_all_boundary_ids() function. You can see the result of not using the correct function in the results section of step-49.

Warning
You should never set the boundary indicator of an interior face (a face not at the boundary of the domain), or set the boundary indicator of an exterior face to numbers::internal_face_boundary_id (this value is reserved for another purpose). Algorithms may not work or produce very confusing results if boundary cells have a boundary indicator of numbers::internal_face_boundary_id or if interior cells have boundary indicators other than numbers::internal_face_boundary_id. Unfortunately, the current object has no means of finding out whether it really is at the boundary of the domain and so cannot determine whether the value you are trying to set makes sense under the current circumstances.
See also
Glossary entry on boundary indicators

◆ set_all_boundary_ids()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::set_all_boundary_ids ( const types::boundary_id  ) const
inherited

Do as set_boundary_id() but also set the boundary indicators of the objects that bound the current object. For example, in 3d, if set_boundary_id() is called on a face, then the boundary indicator of the 4 edges that bound the face remain unchanged. In contrast, if you call the current function, the boundary indicators of face and edges are all set to the given value.

This function is useful if you set boundary indicators of faces in 3d (in 2d, the function does the same as set_boundary_id()) and you do so because you want a curved boundary object to represent the part of the boundary that corresponds to the current face. In that case, the Triangulation class needs to figure out where to put new vertices upon mesh refinement, and higher order Mapping objects also need to figure out where new interpolation points for a curved boundary approximation should be. In either case, the two classes first determine where interpolation points on the edges of a boundary face should be, asking the boundary object, before asking the boundary object for the interpolation points corresponding to the interior of the boundary face. For this to work properly, it is not sufficient to have set the boundary indicator for the face alone, but you also need to set the boundary indicators of the edges that bound the face. This function does all of this at once. You can see the result of not using the correct function in the results section of step-49.

See also
Glossary entry on boundary indicators

◆ get_manifold()

template<int structdim, int dim, int spacedim>
const Manifold< dim, spacedim > & TriaAccessor< structdim, dim, spacedim >::get_manifold ( ) const
inherited

Return a constant reference to the manifold object used for this object.

As explained in the Manifold description for triangulations module, the process involved in finding the appropriate manifold description involves querying both the manifold or boundary indicators. See there for more information.

◆ manifold_id()

template<int structdim, int dim, int spacedim>
types::manifold_id TriaAccessor< structdim, dim, spacedim >::manifold_id ( ) const
inherited

Return the manifold indicator of this object.

If the return value is the special value numbers::flat_manifold_id, then this object is associated with a standard Cartesian Manifold Description.

See also
Glossary entry on manifold indicators

◆ set_all_manifold_ids()

void TriaAccessor< 3, 3, 3 >::set_all_manifold_ids ( const types::manifold_id  manifold_ind) const
inherited

Definition at line 1651 of file tria_accessor.cc.

◆ user_flag_set()

template<int structdim, int dim, int spacedim>
bool TriaAccessor< structdim, dim, spacedim >::user_flag_set ( ) const
inherited

Read the user flag. See GlossUserFlags for more information.

◆ set_user_flag()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::set_user_flag ( ) const
inherited

Set the user flag. See GlossUserFlags for more information.

◆ clear_user_flag()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::clear_user_flag ( ) const
inherited

Clear the user flag. See GlossUserFlags for more information.

◆ recursively_set_user_flag()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::recursively_set_user_flag ( ) const
inherited

Set the user flag for this and all descendants. See GlossUserFlags for more information.

◆ recursively_clear_user_flag()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::recursively_clear_user_flag ( ) const
inherited

Clear the user flag for this and all descendants. See GlossUserFlags for more information.

◆ clear_user_data()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::clear_user_data ( ) const
inherited

Reset the user data to zero, independent if pointer or index. See GlossUserData for more information.

◆ set_user_pointer()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::set_user_pointer ( void *  p) const
inherited

Set the user pointer to p.

Note
User pointers and user indices are mutually exclusive. Therefore, you can only use one of them, unless you call Triangulation::clear_user_data() in between.

See GlossUserData for more information.

◆ clear_user_pointer()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::clear_user_pointer ( ) const
inherited

Reset the user pointer to a nullptr pointer. See GlossUserData for more information.

◆ user_pointer()

template<int structdim, int dim, int spacedim>
void * TriaAccessor< structdim, dim, spacedim >::user_pointer ( ) const
inherited

Access the value of the user pointer. It is in the responsibility of the user to make sure that the pointer points to something useful. You should use the new style cast operator to maintain a minimum of type safety, e.g.

Note
User pointers and user indices are mutually exclusive. Therefore, you can only use one of them, unless you call Triangulation::clear_user_data() in between. A a=static_cast<A>(cell->user_pointer());.

See GlossUserData for more information.

◆ recursively_set_user_pointer()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::recursively_set_user_pointer ( void *  p) const
inherited

Set the user pointer of this object and all its children to the given value. This is useful for example if all cells of a certain subdomain, or all faces of a certain part of the boundary should have user pointers pointing to objects describing this part of the domain or boundary.

Note that the user pointer is not inherited under mesh refinement, so after mesh refinement there might be cells or faces that don't have user pointers pointing to the describing object. In this case, simply loop over all the elements of the coarsest level that has this information, and use this function to recursively set the user pointer of all finer levels of the triangulation.

Note
User pointers and user indices are mutually exclusive. Therefore, you can only use one of them, unless you call Triangulation::clear_user_data() in between.

See GlossUserData for more information.

◆ recursively_clear_user_pointer()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::recursively_clear_user_pointer ( ) const
inherited

Clear the user pointer of this object and all of its descendants. The same holds as said for the recursively_set_user_pointer() function. See GlossUserData for more information.

◆ set_user_index()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::set_user_index ( const unsigned int  p) const
inherited

Set the user index to p.

Note
User pointers and user indices are mutually exclusive. Therefore, you can only use one of them, unless you call Triangulation::clear_user_data() in between. See GlossUserData for more information.

◆ clear_user_index()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::clear_user_index ( ) const
inherited

Reset the user index to 0. See GlossUserData for more information.

◆ user_index()

template<int structdim, int dim, int spacedim>
unsigned int TriaAccessor< structdim, dim, spacedim >::user_index ( ) const
inherited

Access the value of the user index.

Note
User pointers and user indices are mutually exclusive. Therefore, you can only use one of them, unless you call Triangulation::clear_user_data() in between.

See GlossUserData for more information.

◆ recursively_set_user_index()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::recursively_set_user_index ( const unsigned int  p) const
inherited

Set the user index of this object and all its children.

Note that the user index is not inherited under mesh refinement, so after mesh refinement there might be cells or faces that don't have the expected user indices. In this case, simply loop over all the elements of the coarsest level that has this information, and use this function to recursively set the user index of all finer levels of the triangulation.

Note
User pointers and user indices are mutually exclusive. Therefore, you can only use one of them, unless you call Triangulation::clear_user_data() in between.

See GlossUserData for more information.

◆ recursively_clear_user_index()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::recursively_clear_user_index ( ) const
inherited

Clear the user index of this object and all of its descendants. The same holds as said for the recursively_set_user_index() function.

See GlossUserData for more information.

◆ diameter() [2/2]

template<int structdim, int dim, int spacedim>
double TriaAccessor< structdim, dim, spacedim >::diameter ( ) const
inherited

Diameter of the object.

The diameter of an object is computed to be the largest diagonal of the current object. If this object is a quadrilateral, then there are two such diagonal, and if it is a hexahedron, then there are four diagonals that connect "opposite" points. For triangles and tetrahedra, the function simply returns the length of the longest edge.

The situation is more difficult for wedges and pyramids: For wedges, we return the length of the longest diagonal of the three quadrilateral faces or the longest edge length of the two triangular faces. For pyramids, the same principle is applied.

In all of these cases, this definition of "diameter" is not necessarily the true diameter in the sense of the largest distance between points inside the object. Indeed, one can often construct objects for which it is not, though these are generally quite deformed compared to the reference shape. Furthermore, for objects that may use higher order mappings, one may have bulging faces that also create trouble for computing an exact representation of the diameter of the object. That said, the definition used above is completely sufficient for most computations.

◆ enclosing_ball()

template<int structdim, int dim, int spacedim>
std::pair< Point< spacedim >, double > TriaAccessor< structdim, dim, spacedim >::enclosing_ball ( ) const
inherited

Return a pair of Point and double corresponding to the center and the radius of a reasonably small enclosing ball of the object.

The function implements Ritter's O(n) algorithm to get a reasonably small enclosing ball around the vertices of the object. The initial guess for the enclosing ball is taken to be the ball which contains the largest diagonal of the object as its diameter. Starting from such an initial guess, the algorithm tests whether all the vertices of the object (except the vertices of the largest diagonal) are geometrically within the ball. If any vertex (v) is found to be geometrically outside the ball, a new iterate of the ball is constructed by shifting its center and increasing the radius so as to geometrically enclose both the previous ball and the vertex (v). The algorithm terminates when all the vertices are geometrically inside the ball.

If a vertex (v) is geometrically inside a particular iterate of the ball, then it will continue to be so in the subsequent iterates of the ball (this is true by construction).

Note
This function assumes d-linear mapping from the reference cell.

see this and [Ritter 1990]

◆ bounding_box()

template<int structdim, int dim, int spacedim>
BoundingBox< spacedim > TriaAccessor< structdim, dim, spacedim >::bounding_box ( ) const
inherited

Return the smallest bounding box that encloses the object.

Notice that this method is not aware of any mapping you may be using to do your computations. If you are using a mapping object that modifies the position of the vertices, like MappingQEulerian, or MappingFEField, then you should call the function Mapping::get_bounding_box() instead.

◆ extent_in_direction() [1/6]

template<int structdim, int dim, int spacedim>
double TriaAccessor< structdim, dim, spacedim >::extent_in_direction ( const unsigned int  axis) const
inherited

Length of an object in the direction of the given axis, specified in the local coordinate system. See the documentation of GeometryInfo for the meaning and enumeration of the local axes.

Note that the "length" of an object can be interpreted in a variety of ways. Here, we choose it as the maximal length of any of the edges of the object that are parallel to the chosen axis on the reference cell.

◆ extent_in_direction() [2/6]

double TriaAccessor< 1, 1, 1 >::extent_in_direction ( const unsigned int  axis) const
inherited

Definition at line 1578 of file tria_accessor.cc.

◆ extent_in_direction() [3/6]

double TriaAccessor< 1, 1, 2 >::extent_in_direction ( const unsigned int  axis) const
inherited

Definition at line 1589 of file tria_accessor.cc.

◆ extent_in_direction() [4/6]

double TriaAccessor< 2, 2, 2 >::extent_in_direction ( const unsigned int  axis) const
inherited

Definition at line 1600 of file tria_accessor.cc.

◆ extent_in_direction() [5/6]

double TriaAccessor< 2, 2, 3 >::extent_in_direction ( const unsigned int  axis) const
inherited

Definition at line 1614 of file tria_accessor.cc.

◆ extent_in_direction() [6/6]

double TriaAccessor< 3, 3, 3 >::extent_in_direction ( const unsigned int  axis) const
inherited

Definition at line 1629 of file tria_accessor.cc.

◆ minimum_vertex_distance()

template<int structdim, int dim, int spacedim>
double TriaAccessor< structdim, dim, spacedim >::minimum_vertex_distance ( ) const
inherited

Return the minimal distance between any two vertices.

◆ intermediate_point()

template<int structdim, int dim, int spacedim>
Point< spacedim > TriaAccessor< structdim, dim, spacedim >::intermediate_point ( const Point< structdim > &  coordinates) const
inherited

Return a point belonging to the Manifold<dim,spacedim> where this object lives, given its parametric coordinates on the reference structdim cell. This function queries the underlying manifold object, and can be used to obtain the exact geometrical location of arbitrary points on this object.

Notice that the argument coordinates are the coordinates on the reference cell, given in reference coordinates. In other words, the argument provides a weighting between the different vertices. For example, for lines, calling this function with argument Point<1>(.5), is equivalent to asking the line for its center.

Definition at line 1672 of file tria_accessor.cc.

◆ real_to_unit_cell_affine_approximation()

template<int structdim, int dim, int spacedim>
Point< structdim > TriaAccessor< structdim, dim, spacedim >::real_to_unit_cell_affine_approximation ( const Point< spacedim > &  point) const
inherited

This function computes a fast approximate transformation from the real to the unit cell by inversion of an affine approximation of the \(d\)-linear function from the reference \(d\)-dimensional cell.

The affine approximation of the unit to real cell mapping is found by a least squares fit of an affine function to the \(2^d\) vertices of the present object. For any valid mesh cell whose geometry is not degenerate, this operation results in a unique affine mapping. Thus, this function will return a finite result for all given input points, even in cases where the actual transformation by an actual bi-/trilinear or higher order mapping might be singular. Besides only approximating the mapping from the vertex points, this function also ignores the attached manifold descriptions. The result is only exact in case the transformation from the unit to the real cell is indeed affine, such as in one dimension or for Cartesian and affine (parallelogram) meshes in 2d/3d.

For exact transformations to the unit cell, use Mapping::transform_real_to_unit_cell().

Note
If dim<spacedim we first project p onto the plane.

Definition at line 1694 of file tria_accessor.cc.

◆ center()

template<int structdim, int dim, int spacedim>
Point< spacedim > TriaAccessor< structdim, dim, spacedim >::center ( const bool  respect_manifold = false,
const bool  interpolate_from_surrounding = false 
) const
inherited

Center of the object. The center of an object is defined to be the average of the locations of the vertices, which is also where a \(Q_1\) mapping would map the center of the reference cell. However, you can also ask this function to instead return the average of the vertices as computed by the underlying Manifold object associated with the current object, by setting to true the optional parameter respect_manifold. Manifolds would then typically pull back the coordinates of the vertices to a reference domain (not necessarily the reference cell), compute the average there, and then push forward the coordinates of the averaged point to the physical space again; the resulting point is guaranteed to lie within the manifold, even if the manifold is curved.

When the object uses a different manifold description as its surrounding, like when part of the bounding objects of this TriaAccessor use a non-flat manifold description but the object itself is flat, the result given by the TriaAccessor::center() function may not be accurate enough, even when parameter respect_manifold is set to true. If you find this to be case, than you can further refine the computation of the center by setting to true the second additional parameter interpolate_from_surrounding. This computes the location of the center by a so-called transfinite interpolation from the center of all the bounding objects. For a 2d object, it puts a weight of 1/2 on each of the four surrounding lines and a weight -1/4 on the four vertices. This corresponds to a linear interpolation between the descriptions of the four faces, subtracting the contribution of the vertices that is added twice when coming through both lines adjacent to the vertex. In 3d, the weights for faces are 1/2, the weights for lines are -1/4, and the weights for vertices are 1/8. For further information, also confer to the TransfiniteInterpolationManifold class that is able to not only apply this beneficial description to a single cell but all children of a coarse cell.

Definition at line 1713 of file tria_accessor.cc.

◆ barycenter()

template<int structdim, int dim, int spacedim>
Point< spacedim > TriaAccessor< structdim, dim, spacedim >::barycenter ( ) const
inherited

Return the barycenter (also called centroid) of the object. The barycenter for an object \(K\) of dimension \(d\) in \(D\) space dimensions is given by the \(D\)-dimensional vector \(\mathbf x_K\) defined by

\[ \mathbf x_K = \frac{1}{|K|} \int_K \mathbf x \; \textrm{d}x \]

where the measure of the object is given by

\[ |K| = \int_K \mathbf 1 \; \textrm{d}x. \]

This function assumes that \(K\) is mapped by a \(d\)-linear function from the reference \(d\)-dimensional cell. Then the integrals above can be pulled back to the reference cell and evaluated exactly (if through lengthy and, compared to the center() function, expensive computations).

◆ measure()

template<int structdim, int dim, int spacedim>
double TriaAccessor< structdim, dim, spacedim >::measure ( ) const
inherited

Compute the dim-dimensional measure of the object. For a dim-dimensional cell in dim-dimensional space, this equals its volume. On the other hand, for a 2d cell in 3d space, or if the current object pointed to is a 2d face of a 3d cell in 3d space, then the function computes the area the object occupies. For a one-dimensional object, return its length.

The function only computes the measure of cells, faces or edges assumed to be represented by (bi-/tri-)linear mappings. In other words, it only takes into account the locations of the vertices that bound the current object but not how the interior of the object may actually be mapped. In most simple cases, this is exactly what you want. However, for objects that are not "straight", e.g. 2d cells embedded in 3d space as part of a triangulation of a curved domain, two-dimensional faces of 3d cells that are not just parallelograms, or for faces that are at the boundary of a domain that is not just bounded by straight line segments or planes, this function only computes the dim-dimensional measure of a (bi-/tri-)linear interpolation of the "real" object as defined by the manifold or boundary object describing the real geometry of the object in question. If you want to consider the "real" geometry, you will need to compute the measure by integrating a function equal to one over the object, which after applying quadrature equals the summing the JxW values returned by the FEValues or FEFaceValues object you will want to use for the integral.

Note
There is no analytic formula for the area of a bilinear face (i.e., something with a quadrilateral reference cell) in 3D. This function uses the quadrature defined by QGauss<2>(4) to approximate the measure in this case.

◆ is_translation_of()

template<int structdim, int dim, int spacedim>
bool TriaAccessor< structdim, dim, spacedim >::is_translation_of ( const TriaIterator< TriaAccessor< structdim, dim, spacedim > > &  o) const
inherited

Return true if the current object is a translation of the given argument.

Note
For the purpose of a triangulation, cells, faces, etc are only characterized by their vertices. The current function therefore only compares the locations of vertices. For many practical applications, however, it is not only the vertices that determine whether one cell is a translation of another, but also how the cell is mapped from the reference cell to its location in real space. For example, if we are using higher order mappings, then not only do the vertices have to be translations of each other, but also the points along edges. In these questions, therefore, it would be appropriate to ask the mapping, not the current function, whether two objects are translations of each other.

◆ reference_cell()

template<int structdim, int dim, int spacedim>
ReferenceCell TriaAccessor< structdim, dim, spacedim >::reference_cell ( ) const
inherited

Reference cell type of the current object.

◆ n_vertices()

template<int structdim, int dim, int spacedim>
unsigned int TriaAccessor< structdim, dim, spacedim >::n_vertices ( ) const
inherited

Number of vertices.

◆ n_lines()

template<int structdim, int dim, int spacedim>
unsigned int TriaAccessor< structdim, dim, spacedim >::n_lines ( ) const
inherited

Number of lines.

◆ n_faces()

template<int structdim, int dim, int spacedim>
unsigned int TriaAccessor< structdim, dim, spacedim >::n_faces ( ) const
inherited

Return the number of faces for a cell. This function is only implemented for cells (i.e., structdim==dim) to avoid the question of what exactly is meant in a construct such as cell->face(f)->n_faces(). If you want to ask how many bounding lines a face of a 3d cell has, use cell->face(f)->n_lines(); if you want to ask about the number of vertices of a face of a 2d cell, use cell->face(f)->n_vertices().

◆ vertex_indices()

template<int structdim, int dim, int spacedim>
std_cxx20::ranges::iota_view< unsigned int, unsigned int > TriaAccessor< structdim, dim, spacedim >::vertex_indices ( ) const
inherited

Return an object that can be thought of as an array containing all indices from zero to n_vertices().

◆ line_indices()

template<int structdim, int dim, int spacedim>
std_cxx20::ranges::iota_view< unsigned int, unsigned int > TriaAccessor< structdim, dim, spacedim >::line_indices ( ) const
inherited

Return an object that can be thought of as an array containing all indices from zero to n_lines().

◆ face_indices()

template<int structdim, int dim, int spacedim>
std_cxx20::ranges::iota_view< unsigned int, unsigned int > TriaAccessor< structdim, dim, spacedim >::face_indices ( ) const
inherited

Return an object that can be thought of as an array containing all indices from zero to n_faces().

Note
Only implemented for cells (structdim==dim).

◆ set_boundary_id_internal()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::set_boundary_id_internal ( const types::boundary_id  id) const
privateinherited

Like set_boundary_id but without checking for internal faces or invalid ids.

◆ set_bounding_object_indices() [1/2]

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::set_bounding_object_indices ( const std::initializer_list< int > &  new_indices) const
privateinherited

Set the indices of those objects that bound the current object. For example, if the current object represents a cell, then the argument denotes the indices of the faces that bound the cell. If the current object represents a line, the argument denotes the indices of the vertices that bound it. And so on.

◆ set_bounding_object_indices() [2/2]

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::set_bounding_object_indices ( const std::initializer_list< unsigned int > &  new_indices) const
privateinherited

The same as above but for unsigned int.

◆ set_line_orientation()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::set_line_orientation ( const unsigned int  line,
const bool  orientation 
) const
privateinherited

Set the flag indicating, what line_orientation() will return.

It is only possible to set the line_orientation of faces in 3d (i.e. structdim==2 && dim==3).

◆ set_used_flag()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::set_used_flag ( ) const
privateinherited

Set the used flag. Only for internal use in the library.

◆ clear_used_flag()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::clear_used_flag ( ) const
privateinherited

Clear the used flag. Only for internal use in the library.

◆ set_refinement_case()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::set_refinement_case ( const RefinementCase< structdim > &  ref_case) const
privateinherited

Set the RefinementCase<dim> this TriaObject is refined with. Not defined for structdim=1 as lines are always refined resulting in 2 children lines (isotropic refinement).

You should know quite exactly what you are doing if you touch this function. It is exclusively for internal use in the library.

◆ clear_refinement_case()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::clear_refinement_case ( ) const
privateinherited

Clear the RefinementCase<dim> of this TriaObject, i.e. reset it to RefinementCase<dim>::no_refinement.

You should know quite exactly what you are doing if you touch this function. It is exclusively for internal use in the library.

◆ set_children()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::set_children ( const unsigned int  i,
const int  index 
) const
privateinherited

Set the index of the ith child. Since the children come at least in pairs, we need to store the index of only every second child, i.e. of the even numbered children. Make sure, that the index of child i=0 is set first. Calling this function for odd numbered children is not allowed.

◆ clear_children()

template<int structdim, int dim, int spacedim>
void TriaAccessor< structdim, dim, spacedim >::clear_children ( ) const
privateinherited

Clear the child field, i.e. set it to a value which indicates that this cell has no children.

◆ copy_from()

void TriaAccessorBase< structdim, dim, spacedim >::copy_from ( const TriaAccessorBase< structdim, dim, spacedim > &  )
protectedinherited

Copy operator. Since this is only called from iterators, do not return anything, since the iterator will return itself.

This method is protected, since it is only to be called from the iterator class.

◆ operator<()

bool TriaAccessorBase< structdim, dim, spacedim >::operator< ( const TriaAccessorBase< structdim, dim, spacedim > &  other) const
protectedinherited

Comparison operator for accessors. This operator is used when comparing iterators into objects of a triangulation, for example when putting them into a std::map.

If structure_dimension is less than dimension, we simply compare the index of such an object because faces and edges do not have levels. If structure_dimension equals dimension, we compare the level first, and the index only if levels are equal.

◆ operator==()

bool TriaAccessorBase< structdim, dim, spacedim >::operator== ( const TriaAccessorBase< structdim, dim, spacedim > &  ) const
protectedinherited

Compare for equality.

◆ operator!=()

bool TriaAccessorBase< structdim, dim, spacedim >::operator!= ( const TriaAccessorBase< structdim, dim, spacedim > &  ) const
protectedinherited

Compare for inequality.

◆ operator++()

void TriaAccessorBase< structdim, dim, spacedim >::operator++ ( )
protectedinherited

This operator advances the iterator to the next element.

For dim=1 only: The next element is next on this level if there are more. If the present element is the last on this level, the first on the next level is accessed.

◆ operator--()

void TriaAccessorBase< structdim, dim, spacedim >::operator-- ( )
protectedinherited

This operator moves the iterator to the previous element.

For dim=1 only: The previous element is previous on this level if index>0. If the present element is the first on this level, the last on the previous level is accessed.

◆ objects()

::internal::TriangulationImplementation::TriaObjects & TriaAccessorBase< structdim, dim, spacedim >::objects ( ) const
protectedinherited

Access to the other objects of a Triangulation with same dimension.

◆ level()

int TriaAccessorBase< structdim, dim, spacedim >::level ( ) const
inherited

For cells, this function returns the level within the mesh hierarchy at which this cell is located. For all other objects, the function returns zero.

Note
Within a Triangulation object, cells are uniquely identified by a pair (level, index) where the former is the cell's refinement level and the latter is the index of the cell within this refinement level (the former being what this function returns). Consequently, there may be multiple cells on different refinement levels but with the same index within their level. Contrary to this, if the current object corresponds to a face or edge, then the object is uniquely identified solely by its index as faces and edges do not have a refinement level. For these objects, the current function always returns zero as the level.

◆ index()

int TriaAccessorBase< structdim, dim, spacedim >::index ( ) const
inherited

Return the index of the element presently pointed to on the present level.

Within a Triangulation object, cells are uniquely identified by a pair (level, index) where the former is the cell's refinement level and the latter is the index of the cell within this refinement level (the latter being what this function returns). Consequently, there may be multiple cells on different refinement levels but with the same index within their level. Contrary to this, if the current object corresponds to a face or edge, then the object is uniquely identified solely by its index as faces and edges do not have a refinement level.

Note
The indices objects returned by this function are not a contiguous set of numbers on each level: going from cell to cell, some of the indices in a level may be unused.
If the triangulation is actually of type parallel::distributed::Triangulation then the indices are relatively only to that part of the distributed triangulation that is stored on the current processor. In other words, cells living in the partitions of the triangulation stored on different processors may have the same index even if they refer to the same cell, and the may have different indices even if they do refer to the same cell (e.g., if a cell is owned by one processor but is a ghost cell on another).

◆ state()

IteratorState::IteratorStates TriaAccessorBase< structdim, dim, spacedim >::state ( ) const
inherited

Return the state of the iterator. For the different states an accessor can be in, refer to the TriaRawIterator documentation.

◆ get_triangulation()

const Triangulation< dim, spacedim > & TriaAccessorBase< structdim, dim, spacedim >::get_triangulation ( ) const
inherited

Return a reference to the triangulation which the object pointed to by this class belongs to.

Friends And Related Symbol Documentation

◆ Triangulation< dim, spacedim >

template<int dim, int spacedim = dim>
friend class Triangulation< dim, spacedim >
friend

Definition at line 4187 of file tria_accessor.h.

◆ parallel::TriangulationBase< dim, spacedim >

template<int dim, int spacedim = dim>
friend class parallel::TriangulationBase< dim, spacedim >
friend

Definition at line 4187 of file tria_accessor.h.

◆ ::internal::TriangulationImplementation::Implementation

template<int dim, int spacedim = dim>
friend struct ::internal::TriangulationImplementation::Implementation
friend

Definition at line 4193 of file tria_accessor.h.

◆ ::internal::TriangulationImplementation::ImplementationMixedMesh

template<int dim, int spacedim = dim>
friend struct ::internal::TriangulationImplementation::ImplementationMixedMesh
friend

Definition at line 4194 of file tria_accessor.h.

Member Data Documentation

◆ space_dimension

constexpr unsigned int TriaAccessorBase< structdim, dim, spacedim >::space_dimension
staticconstexprinherited

Dimension of the space the object represented by this accessor lives in. For example, if this accessor represents a quadrilateral that is part of a two-dimensional surface in four-dimensional space, then this value is four.

Definition at line 316 of file tria_accessor.h.

◆ dimension

constexpr unsigned int TriaAccessorBase< structdim, dim, spacedim >::dimension
staticconstexprinherited

Dimensionality of the object that the thing represented by this accessor is part of. For example, if this accessor represents a line that is part of a hexahedron, then this value will be three.

Definition at line 323 of file tria_accessor.h.

◆ structure_dimension

const unsigned int TriaAccessorBase< structdim, dim, spacedim >::structure_dimension
staticinherited

Dimensionality of the current object represented by this accessor. For example, if it is line (irrespective of whether it is part of a 2d or 3d subobject), then this value equals 1.

Definition at line 330 of file tria_accessor.h.

◆ present_level

typename::internal::TriaAccessorImplementation::PresentLevelType<structdim,dim>::type TriaAccessorBase< structdim, dim, spacedim >::present_level
protectedinherited

The level if this is a cell (structdim==dim). Else, contains zero.

Definition at line 528 of file tria_accessor.h.

◆ present_index

int TriaAccessorBase< structdim, dim, spacedim >::present_index
protectedinherited

Used to store the index of the element presently pointed to on the level presently used.

Definition at line 534 of file tria_accessor.h.

◆ tria

const Triangulation<dim, spacedim>* TriaAccessorBase< structdim, dim, spacedim >::tria
protectedinherited

Pointer to the triangulation which we act on.

Definition at line 539 of file tria_accessor.h.


The documentation for this class was generated from the following files: