Reference documentation for deal.II version 9.4.1
|
#include <deal.II/grid/reference_cell.h>
Public Member Functions | |
constexpr | ReferenceCell () |
Querying information about the kind of reference cells | |
bool | is_hyper_cube () const |
bool | is_simplex () const |
unsigned int | get_dimension () const |
Shape functions, mappings, quadratures defined on a reference cell | |
template<int dim> | |
double | d_linear_shape_function (const Point< dim > &xi, const unsigned int i) const |
template<int dim> | |
Tensor< 1, dim > | d_linear_shape_function_gradient (const Point< dim > &xi, const unsigned int i) const |
template<int dim, int spacedim = dim> | |
std::unique_ptr< Mapping< dim, spacedim > > | get_default_mapping (const unsigned int degree) const |
template<int dim, int spacedim = dim> | |
const Mapping< dim, spacedim > & | get_default_linear_mapping () const |
template<int dim> | |
Quadrature< dim > | get_gauss_type_quadrature (const unsigned n_points_1D) const |
template<int dim> | |
Quadrature< dim > | get_midpoint_quadrature () const |
template<int dim> | |
const Quadrature< dim > & | get_nodal_type_quadrature () const |
Querying the number of building blocks of a reference cell | |
unsigned int | n_vertices () const |
std_cxx20::ranges::iota_view< unsigned int, unsigned int > | vertex_indices () const |
template<int dim> | |
Point< dim > | vertex (const unsigned int v) const |
unsigned int | n_lines () const |
std_cxx20::ranges::iota_view< unsigned int, unsigned int > | line_indices () const |
unsigned int | n_faces () const |
std_cxx20::ranges::iota_view< unsigned int, unsigned int > | face_indices () const |
unsigned int | n_isotropic_children () const |
std_cxx20::ranges::iota_view< unsigned int, unsigned int > | isotropic_child_indices () const |
ReferenceCell | face_reference_cell (const unsigned int face_no) const |
double | volume () const |
template<int dim> | |
Point< dim > | barycenter () const |
template<int dim> | |
bool | contains_point (const Point< dim > &p, const double tolerance=0) const |
template<int dim> | |
Tensor< 1, dim > | unit_tangential_vectors (const unsigned int face_no, const unsigned int i) const |
template<int dim> | |
Tensor< 1, dim > | unit_normal_vectors (const unsigned int face_no) const |
template<typename T , std::size_t N> | |
unsigned char | compute_orientation (const std::array< T, N > &vertices_0, const std::array< T, N > &vertices_1) const |
template<typename T , std::size_t N> | |
std::array< T, N > | permute_according_orientation (const std::array< T, N > &vertices, const unsigned int orientation) const |
ArrayView< const unsigned int > | faces_for_given_vertex (const unsigned int vertex_index) const |
Relationships between objects in the cell and on faces | |
unsigned int | child_cell_on_face (const unsigned int face, const unsigned int subface, const unsigned char face_orientation=1) const |
std::array< unsigned int, 2 > | standard_vertex_to_face_and_vertex_index (const unsigned int vertex) const |
std::array< unsigned int, 2 > | standard_line_to_face_and_line_index (const unsigned int line) const |
unsigned int | face_to_cell_lines (const unsigned int face, const unsigned int line, const unsigned char face_orientation) const |
unsigned int | face_to_cell_vertices (const unsigned int face, const unsigned int vertex, const unsigned char face_orientation) const |
unsigned int | standard_to_real_face_vertex (const unsigned int vertex, const unsigned int face, const unsigned char face_orientation) const |
unsigned int | standard_to_real_face_line (const unsigned int line, const unsigned int face, const unsigned char face_orientation) const |
bool | standard_vs_true_line_orientation (const unsigned int line, const unsigned char face_orientation, const bool line_orientation) const |
Translating between deal.II indexing and formats used by other programs | |
unsigned int | exodusii_vertex_to_deal_vertex (const unsigned int vertex_n) const |
unsigned int | exodusii_face_to_deal_face (const unsigned int face_n) const |
unsigned int | unv_vertex_to_deal_vertex (const unsigned int vertex_n) const |
unsigned int | vtk_linear_type () const |
unsigned int | vtk_quadratic_type () const |
unsigned int | vtk_lagrange_type () const |
unsigned int | gmsh_element_type () const |
Other functions | |
std::string | to_string () const |
constexpr | operator std::uint8_t () const |
constexpr bool | operator== (const ReferenceCell &type) const |
constexpr bool | operator!= (const ReferenceCell &type) const |
template<class Archive > | |
void | serialize (Archive &archive, const unsigned int) |
Static Public Member Functions | |
static ReferenceCell | n_vertices_to_type (const int dim, const unsigned int n_vertices) |
Private Member Functions | |
constexpr | ReferenceCell (const std::uint8_t kind) |
Private Attributes | |
std::uint8_t | kind |
Friends | |
constexpr ReferenceCell | internal::ReferenceCell::make_reference_cell_from_int (const std::uint8_t) |
std::ostream & | operator<< (std::ostream &out, const ReferenceCell &reference_cell) |
std::istream & | operator>> (std::istream &in, ReferenceCell &reference_cell) |
A type that describes the kinds of reference cells that can be used. This includes quadrilaterals and hexahedra (i.e., "hypercubes"), triangles and tetrahedra (simplices), and the pyramids and wedges necessary when using mixed 3d meshes. This class then describes geometric, topological, and other kinds of information about these kinds of reference cells. This includes how many vertices or faces a certain kind of reference cell has (topological information), where these vertices lie, what the cell's volume or center of mass is (geometric information), and how to output these cells in various output formats or what appropriate quadrature rules are. The documentation of this class is separated into a number of sections to group the many member functions into different categories such as those mentioned above.
Objects of this type should not be created in user code, and as a consequence the class does not have a user-accessible constructor other than the default constructor (which creates an invalid object). Rather, there is a finite number of specific reference cell objects defined in the ReferenceCells namespace that completely enumerate all of the possible values. User codes should therefore rely exclusively on assigning ReferenceCell objects from these special objects, and comparing against those special objects.
The purposes and intents of this class are described in the reference cell glossary entry.
Definition at line 100 of file reference_cell.h.
|
inlineconstexpr |
Default constructor. Initialize this object as an invalid object. The end result is that the current object equals ReferenceCells::Invalid.
Generally, ReferenceCell objects are created by assignment from the special objects in namespace ReferenceCells, which is the only way to obtain a valid object.
Definition at line 830 of file reference_cell.h.
|
inlineconstexprprivate |
Constructor. This is the constructor used to create the different static
member variables of this class. It is private
but can be called by a function in an internal namespace that is a friend
of this class.
Definition at line 729 of file reference_cell.h.
|
inlinestatic |
Return the correct ReferenceCell for a given structural dimension and number of vertices. For example, if dim==2
and n_vertices==4
, this function will return ReferenceCells::Quadrilateral. But if dim==3
and n_vertices==4
, it will return ReferenceCells::Tetrahedron.
Definition at line 1844 of file reference_cell.h.
|
inline |
Return true
if the object is a ReferenceCells::Vertex, ReferenceCells::Line, ReferenceCells::Quadrilateral, or ReferenceCells::Hexahedron.
Definition at line 912 of file reference_cell.h.
|
inline |
Return true if the object is a Vertex, Line, Triangle, or Tetrahedron.
Definition at line 922 of file reference_cell.h.
|
inline |
Return the dimension of the reference cell represented by the current object.
Definition at line 932 of file reference_cell.h.
|
inline |
Compute the value of the \(i\)-th linear shape function at location \(\xi\) for the current reference-cell type.
Definition at line 1886 of file reference_cell.h.
|
inline |
Compute the gradient of the \(i\)-th linear shape function at location \(\xi\) for the current reference-cell type.
Definition at line 1980 of file reference_cell.h.
std::unique_ptr< Mapping< dim, spacedim > > ReferenceCell::get_default_mapping | ( | const unsigned int | degree | ) | const |
Return a default mapping of degree degree
matching the current reference cell. If this reference cell is a hypercube, then the returned mapping is a MappingQ; otherwise, it is an object of type MappingFE initialized with FE_SimplexP (if the reference cell is a triangle or tetrahedron), with FE_PyramidP (if the reference cell is a pyramid), or with FE_WedgeP (if the reference cell is a wedge).
Definition at line 111 of file reference_cell.cc.
const Mapping< dim, spacedim > & ReferenceCell::get_default_linear_mapping |
Return a default linear mapping matching the current reference cell. If this reference cell is a hypercube, then the returned mapping is a MappingQ1; otherwise, it is an object of type MappingFE initialized with FE_SimplexP (if the reference cell is a triangle or tetrahedron), with FE_PyramidP (if the reference cell is a pyramid), or with FE_WedgeP (if the reference cell is a wedge). In other words, the term "linear" in the name of the function has to be understood as \(d\)-linear (i.e., bilinear or trilinear) for some of the coordinate directions.
Definition at line 138 of file reference_cell.cc.
Quadrature< dim > ReferenceCell::get_gauss_type_quadrature | ( | const unsigned | n_points_1D | ) | const |
Return a Gauss-type quadrature matching the given reference cell (QGauss, QGaussSimplex, QGaussPyramid, QGaussWedge).
[in] | n_points_1D | The number of quadrature points in each direction (QGauss) or an indication of what polynomial degree needs to be integrated exactly for the other types. |
Definition at line 176 of file reference_cell.cc.
Quadrature< dim > ReferenceCell::get_midpoint_quadrature |
Return a quadrature object that has a single quadrature point at the barycenter of the cell with quadrature weight equal to the volume of the reference cell. This quadrature formula is exact for integrals of constant and linear integrands.
The object returned by this function generalizes what the QMidpoint class represents to other reference cells.
Definition at line 955 of file reference_cell.h.
const Quadrature< dim > & ReferenceCell::get_nodal_type_quadrature |
Return a quadrature rule whose quadrature points are the vertices of the given reference cell. For 1d line segments, this corresponds to the quadrature points of the trapezoidal rule, which by taking tensor products easily generalizes also to other hypercube elements (see also QTrapezoid). For all reference cell shapes, the quadrature points are ordered in the same order as the vertices of the reference cell.
Definition at line 198 of file reference_cell.cc.
|
inline |
Return the number of vertices that make up the reference cell in question. A vertex is a "corner" (a zero-dimensional object) of the reference cell.
Definition at line 964 of file reference_cell.h.
|
inline |
Return an object that can be thought of as an array containing all indices from zero to n_vertices().
Definition at line 1185 of file reference_cell.h.
Return the location of the v
th vertex of the reference cell that corresponds to the current object.
Because the ReferenceCell class does not have a dim
argument, it has to be explicitly specified in the call to this function.
Definition at line 1017 of file reference_cell.h.
|
inline |
Return the number of lines that make up the reference cell in question. A line is an "edge" (a one-dimensional object) of the reference cell.
Definition at line 990 of file reference_cell.h.
|
inline |
Return an object that can be thought of as an array containing all indices from zero to n_lines().
Definition at line 1193 of file reference_cell.h.
|
inline |
Return the number of faces that make up the reference cell in question. A face is a (dim-1)
-dimensional object bounding the reference cell.
Definition at line 1112 of file reference_cell.h.
|
inline |
Return an object that can be thought of as an array containing all indices from zero to n_faces().
Definition at line 1138 of file reference_cell.h.
|
inline |
Return the number of cells one would get by isotropically refining the current cell. Here, "isotropic refinement" means that we subdivide in each "direction" of a cell. For example, a square would be refined into four children by introducing new vertices along each edge and a new vertex in the cell center. For triangles, one would introduce new vertices at the center of each edge, and connect them to obtain four children. Similar constructions can be done for the other reference cell types.
Definition at line 1146 of file reference_cell.h.
|
inline |
Return an object that can be thought of as an array containing all indices from zero to n_isotropic_children().
Definition at line 1177 of file reference_cell.h.
|
inline |
Return the reference-cell type of face face_no
of the current object. For example, if the current object is ReferenceCells::Tetrahedron, then face_no
must be between in the interval \([0,4)\) and the function will always return ReferenceCells::Triangle. If the current object is ReferenceCells::Hexahedron, then face_no
must be between in the interval \([0,6)\) and the function will always return ReferenceCells::Quadrilateral. For wedges and pyramids, the returned object may be either ReferenceCells::Triangle or ReferenceCells::Quadrilateral, depending on the given index.
Definition at line 1201 of file reference_cell.h.
|
inline |
Return which child cells are adjacent to a certain face of the mother cell.
For example, in 2D the layout of a quadrilateral cell is as follows:
* 3 * 2-->--3 * | | * 0 ^ ^ 1 * | | * 0-->--1 * 2 *
Vertices and faces are indicated with their numbers, faces also with their directions.
Now, when refined, the layout is like this:
* *---*---* * | 2 | 3 | * *---*---* * | 0 | 1 | * *---*---* *
Thus, the child cells on face 0 are (ordered in the direction of the face) 0 and 2, on face 3 they are 2 and 3, etc.
For three spatial dimensions, the exact order of the children is laid down in the general documentation of this class.
The face_orientation
argument is meant exclusively for quadrilaterals and hexahedra at the moment. It determines how this function handles faces oriented in the standard and non-standard orientation. It represents a bit-code for the overall face_orientation
, face_flip
and face_rotation
and defaults to the standard orientation. The concept of face orientations is explained in this glossary entry.
Definition at line 1239 of file reference_cell.h.
|
inline |
For a given vertex in a cell, return a pair of a face index and a vertex index within this face.
Definition at line 1310 of file reference_cell.h.
|
inline |
For a given line in a cell, return a pair of a face index and a line index within this face.
Definition at line 1370 of file reference_cell.h.
|
inline |
Map face line number to cell line number.
Definition at line 1437 of file reference_cell.h.
|
inline |
Map face vertex number to cell vertex number.
Definition at line 1519 of file reference_cell.h.
|
inline |
Correct vertex index depending on face orientation.
Definition at line 1604 of file reference_cell.h.
|
inline |
Correct line index depending on face orientation.
Definition at line 1703 of file reference_cell.h.
|
inline |
Return whether the line with index line
is oriented in standard direction within a cell, given the face_orientation
of the face within the current cell, and line_orientation
flag for the line within that face. true
indicates that the line is oriented from vertex 0 to vertex 1, whereas it is the other way around otherwise. In 1d and 2d, this is always true
, but in 3d it may be different, see the respective discussion in the documentation of the GeometryInfo class.
Definition at line 2234 of file reference_cell.h.
|
inline |
Return the \(d\)-dimensional volume of the reference cell that corresponds to the current object, where \(d\) is the dimension of the space it lives in. For example, since the quadrilateral reference cell is \([0,1]^2\), its volume is one, whereas the volume of the reference triangle is 0.5 because it occupies the area \(\{0 \le x,y \le 1, x+y\le 1\}\).
For ReferenceCells::Vertex, the reference cell is a zero-dimensional point in a zero-dimensional space. As a consequence, one cannot meaningfully define a volume for it. The function returns one for this case, because this makes it possible to define useful quadrature rules based on the center of a reference cell and its volume.
Definition at line 2010 of file reference_cell.h.
Return the barycenter (i.e., the center of mass) of the reference cell that corresponds to the current object. The function is not called center()
because one can define the center of an object in a number of different ways whereas the barycenter of a reference cell \(K\) is unambiguously defined as
\[ \mathbf x_K = \frac{1}{V} \int_K \mathbf x \; dx \]
where \(V\) is the volume of the reference cell (see also the volume() function).
Definition at line 2037 of file reference_cell.h.
|
inline |
Return true if the given point is inside the reference cell of the present space dimension up to some tolerance. This function accepts an additional parameter (which defaults to zero) which specifies by how much the point position may actually be outside the true reference cell. This is useful because in practice we may often not be able to compute the coordinates of a point in reference coordinates exactly, but only up to numerical roundoff. For example, strictly speaking one would expect that for points on the boundary of the reference cell, the function would return true
if the tolerance was zero. But in practice, this may or may not actually be true; for example, the point \((1/3, 2/3)\) is on the boundary of the reference triangle because \(1/3+2/3 \le 1\), but since neither of its coordinates are exactly representable in floating point arithmetic, the floating point representations of \(1/3\) and \(2/3\) may or may not add up to anything that is less than or equal to one.
The tolerance parameter may be less than zero, indicating that the point should be safely inside the cell.
Definition at line 2066 of file reference_cell.h.
|
inline |
Definition at line 2128 of file reference_cell.h.
|
inline |
Return the unit normal vector of a face of the reference cell.
Definition at line 2204 of file reference_cell.h.
|
inline |
Determine the orientation of the current entity described by its vertices var_1
relative to an entity described by var_0
.
Definition at line 2342 of file reference_cell.h.
|
inline |
Inverse function of compute_orientation().
Definition at line 2437 of file reference_cell.h.
|
inline |
Return a vector of faces a given vertex_index
belongs to.
Definition at line 846 of file reference_cell.h.
Map an ExodusII vertex number to a deal.II vertex number.
Definition at line 243 of file reference_cell.cc.
Map an ExodusII face number to a deal.II face number.
Definition at line 289 of file reference_cell.cc.
Map a UNV vertex number to a deal.II vertex number.
Definition at line 339 of file reference_cell.cc.
unsigned int ReferenceCell::vtk_linear_type | ( | ) | const |
Return a VTK linear shape constant that corresponds to the reference cell.
Definition at line 375 of file reference_cell.cc.
unsigned int ReferenceCell::vtk_quadratic_type | ( | ) | const |
Return a VTK quadratic shape constant that corresponds to the reference cell.
Definition at line 404 of file reference_cell.cc.
unsigned int ReferenceCell::vtk_lagrange_type | ( | ) | const |
Return a VTK Lagrange shape constant that corresponds to the reference cell.
Definition at line 433 of file reference_cell.cc.
unsigned int ReferenceCell::gmsh_element_type | ( | ) | const |
Return the GMSH element type code that corresponds to the reference cell.
Definition at line 462 of file reference_cell.cc.
std::string ReferenceCell::to_string | ( | ) | const |
Return a text representation of the reference cell represented by the current object.
Definition at line 81 of file reference_cell.cc.
|
inlineconstexpr |
Conversion operator to an integer.
Definition at line 735 of file reference_cell.h.
|
inlineconstexpr |
Operator for equality comparison.
Definition at line 743 of file reference_cell.h.
|
inlineconstexpr |
Operator for inequality comparison.
Definition at line 751 of file reference_cell.h.
|
inline |
Write and read the data of this object from a stream for the purpose of serialization using the BOOST serialization library.
Definition at line 838 of file reference_cell.h.
|
friend |
A kind of constructor – not quite private because it can be called by anyone, but at least hidden in an internal namespace.
|
friend |
Output operator that writes the reference_cell
object to the stream in a text format in which the object is represented by an integer. The details of which integer value represents each kind of reference cell is unimportant and consequently not specified. If you want a string representation of what a ReferenceCell is, use ReferenceCell::to_string().
Definition at line 546 of file reference_cell.cc.
|
friend |
Input operator that reads the reference_cell
object from the stream in a text format in which the object is represented by an integer. Which specific integer value represents which reference cell is unspecified, but the function uses the same translation as the corresponding output operator<<
.
Definition at line 560 of file reference_cell.cc.
|
private |
The variable that stores what this object actually corresponds to.
Definition at line 683 of file reference_cell.h.