15#ifndef dealii_tria_reference_cell_h
16#define dealii_tria_reference_cell_h
29#include <boost/container/small_vector.hpp>
38template <
int dim,
int spacedim>
171 const unsigned int i)
const;
182 template <
int dim,
int spacedim = dim>
183 std::unique_ptr<Mapping<dim, spacedim>>
197 template <
int dim,
int spacedim = dim>
276 vertex(
const unsigned int v)
const;
357 static constexpr unsigned char
370 static constexpr unsigned char
416 const unsigned int subface,
417 const unsigned char face_orientation =
429 std::array<unsigned int, 2>
441 std::array<unsigned int, 2>
454 const unsigned int vertex)
const;
461 const unsigned int line,
462 const unsigned char face_orientation)
const;
469 const unsigned int vertex,
470 const unsigned char face_orientation)
const;
494 const unsigned int vertex)
const;
501 const unsigned int face,
502 const unsigned char face_orientation)
const;
509 const unsigned int face,
510 const unsigned char face_orientation)
const;
524 const unsigned int face,
525 const unsigned char face_orientation,
526 const bool line_orientation)
const;
611 const unsigned int i)
const;
644 template <
typename T, std::
size_t N>
647 const std::array<T, N> &vertices_1)
const;
674 template <
typename T>
693 template <
typename T, std::
size_t N>
696 const unsigned int orientation)
const;
709 template <
typename T>
710 boost::container::small_vector<T, 8>
712 const unsigned char orientation)
const;
790 const std::array<unsigned, dim> &node_indices,
791 const std::array<unsigned, dim> &nodes_per_direction,
792 const bool legacy_format)
const;
825 constexpr operator std::uint8_t()
const;
844 template <
class Archive>
846 serialize(Archive &archive,
const unsigned int );
851 static constexpr std::size_t
884 {{{1, 0}}, {{0, 1}}}};
917 friend std::ostream &
920 friend std::istream &
952inline constexpr ReferenceCell::operator std::uint8_t()
const
980#ifndef DEAL_II_CXX14_CONSTEXPR_BUG
983 Assert((kind == std::numeric_limits<std::uint8_t>::max()) || (kind < 8),
1022 std::numeric_limits<std::uint8_t>::max());
1051template <
class Archive>
1060inline constexpr std::size_t
1081 {{{0, 2}}, {{0, 1}}, {{1, 2}}}};
1087 {{{0, 1, 2}}, {{0, 1, 3}}, {{0, 2, 3}}, {{1, 2, 3}}}};
1177 std::vector<double>({
volume()}));
1494inline constexpr unsigned char
1505inline constexpr unsigned char
1517 const unsigned int face,
1518 const unsigned int subface,
1519 const unsigned char combined_face_orientation)
const
1535 {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
1537 return subcells[face][subface];
1541 const auto [face_orientation, face_rotation, face_flip] =
1561 const auto [face_orientation, face_rotation, face_flip] =
1581inline std::array<unsigned int, 2>
1583 const unsigned int vertex)
const
1599 {{{0, 0}}, {{0, 1}}, {{1, 1}}, {{X, X}}, {{X, X}}, {{X, X}}}};
1611 {{{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 2}}, {{X, X}}, {{X, X}}}};
1618 {{{0, 0}}, {{0, 1}}, {{0, 2}}, {{0, 3}}, {{1, 2}}, {{X, X}}}};
1625 {{{0, 1}}, {{0, 0}}, {{0, 2}}, {{1, 0}}, {{1, 1}}, {{1, 2}}}};
1643inline std::array<unsigned int, 2>
1645 const unsigned int line)
const
1661 static const std::array<unsigned int, 2> table[6] = {
1662 {{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 1}}, {{1, 2}}, {{2, 1}}};
1668 static const std::array<unsigned int, 2> table[8] = {{{0, 0}},
1681 static const std::array<unsigned int, 2> table[9] = {{{0, 0}},
1707 const unsigned int vertex)
const
1720 {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
1721 return table[line][
vertex];
1726 {{{0, 2}}, {{1, 3}}, {{0, 1}}, {{2, 3}}}};
1727 return table[line][
vertex];
1732 {{{0, 1}}, {{1, 2}}, {{2, 0}}, {{0, 3}}, {{1, 3}}, {{2, 3}}}};
1733 return table[line][
vertex];
1745 return table[line][
vertex];
1758 return table[line][
vertex];
1776 return table[line][
vertex];
1789 const unsigned int face,
1790 const unsigned int line,
1791 const unsigned char combined_face_orientation)
const
1807 const auto [face_orientation, face_rotation, face_flip] =
1811 face, line, face_orientation, face_flip, face_rotation);
1819 const auto [face_orientation, face_rotation, face_flip] =
1823 face, line, face_orientation, face_flip, face_rotation);
1828 {{{0, 1, 2}}, {{0, 3, 4}}, {{2, 5, 3}}, {{1, 4, 5}}}};
1831 line, face, combined_face_orientation)];
1843 line, face, combined_face_orientation)];
1855 line, face, combined_face_orientation)];
1859 const auto [face_orientation, face_rotation, face_flip] =
1863 face, line, face_orientation, face_flip, face_rotation);
1876 const unsigned int face,
1877 const unsigned int vertex,
1878 const unsigned char combined_face_orientation)
const
1885 Assert(combined_face_orientation ==
1887 ExcMessage(
"In 1D, all faces must have the default orientation."));
1900 const auto [face_orientation, face_rotation, face_flip] =
1904 face,
vertex, face_orientation, face_flip, face_rotation);
1909 {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
1911 return table[face][combined_face_orientation ==
1918 const auto [face_orientation, face_rotation, face_flip] =
1922 face,
vertex, face_orientation, face_flip, face_rotation);
1927 {{{0, 1, 2}}, {{1, 0, 3}}, {{0, 2, 3}}, {{2, 1, 3}}}};
1930 vertex, face, combined_face_orientation)];
1943 vertex, face, combined_face_orientation)];
1956 vertex, face, combined_face_orientation)];
1960 const auto [face_orientation, face_rotation, face_flip] =
1964 face,
vertex, face_orientation, face_flip, face_rotation);
1978 const unsigned int vertex)
const
1980 return this->
template vertex<dim>(
1988 const unsigned int vertex,
1989 const unsigned int face,
1990 const unsigned char face_orientation)
const
2003 "In 1D, all faces must have the default orientation."));
2036 const unsigned int line,
2037 const unsigned int face,
2038 const unsigned char combined_face_orientation)
const
2060 return triangle_table[combined_face_orientation][line];
2064 const auto [face_orientation, face_rotation, face_flip] =
2074 return triangle_table[combined_face_orientation][line];
2079 const auto [face_orientation, face_rotation, face_flip] =
2088 return triangle_table[combined_face_orientation][line];
2100 return table[combined_face_orientation][line];
2190 ExcMessage(
"The combination of dim = " + std::to_string(dim) +
2191 " and n_vertices = " + std::to_string(
n_vertices) +
2192 " does not correspond to a known reference cell type."));
2201 const unsigned int i)
const
2218 return 1.0 - xi[
std::min(0, dim - 1)] -
2234 return 1.0 - xi[
std::min(0, dim - 1)] -
2249 const double Q14 = 0.25;
2251 const double r = xi[
std::min(0, dim - 1)];
2252 const double s = xi[
std::min(1, dim - 1)];
2253 const double t = xi[
std::min(2, dim - 1)];
2255 const double ratio =
2256 (std::fabs(t - 1.0) > 1.0e-14 ? (r * s * t) / (1.0 - t) : 0.0);
2259 return Q14 * ((1.0 - r) * (1.0 - s) - t + ratio);
2261 return Q14 * ((1.0 + r) * (1.0 - s) - t - ratio);
2263 return Q14 * ((1.0 - r) * (1.0 + s) - t - ratio);
2265 return Q14 * ((1.0 + r) * (1.0 + s) - t + ratio);
2290 const unsigned int i)
const
2369 return Point<dim>(1. / 4., 1. / 4., 1. / 4.);
2375 return Point<dim>(1. / 2., 1. / 2., 1. / 2.);
2394 constexpr unsigned int x_coordinate = 0;
2395 constexpr unsigned int y_coordinate = (dim >= 2 ? 1 : 0);
2396 constexpr unsigned int z_coordinate = (dim >= 3 ? 2 : 0);
2406 ExcMessage(
"Vertices are zero-dimensional objects and "
2407 "as a consequence have no coordinates. You "
2408 "cannot meaningfully ask whether a point is "
2409 "inside a vertex (within a certain tolerance) "
2410 "without coordinate values."));
2417 for (
unsigned int d = 0; d < dim; ++d)
2418 if ((p[d] < -tolerance) || (p[d] > 1 + tolerance))
2426 for (
unsigned int d = 0; d < dim; ++d)
2427 if (p[d] < -tolerance)
2439 for (
unsigned int d = 0; d < dim; ++d)
2441 return (sum <= 1 + tolerance *
std::sqrt(1. * dim));
2446 if (p[z_coordinate] < -tolerance)
2450 if (p[z_coordinate] > 1 + tolerance)
2457 const double distance_from_axis =
2458 std::max(std::fabs(p[x_coordinate]), std::fabs(p[y_coordinate]));
2462 return (distance_from_axis <= 1 + tolerance - p[z_coordinate]);
2471 if ((p[x_coordinate] < -tolerance) || (p[y_coordinate] < -tolerance))
2474 const double sum = p[x_coordinate] + p[y_coordinate];
2475 if (sum > 1 + tolerance *
std::sqrt(2.0))
2478 if (p[z_coordinate] < -tolerance)
2480 if (p[z_coordinate] > 1 + tolerance)
2497 const unsigned int i)
const
2513 static const std::array<Tensor<1, dim>, 3> table = {
2518 return table[face_no];
2532 +
std::pow(1.0 / 3.0, 1.0 / 4.0))}}}};
2534 return table[face_no][i];
2541 {{
Point<dim>(+1.0 / sqrt(2.0), 0, +1.0 / sqrt(2.0)),
2543 {{
Point<dim>(+1.0 / sqrt(2.0), 0, -1.0 / sqrt(2.0)),
2546 Point<dim>(0, +1.0 / sqrt(2.0), +1.0 / sqrt(2.0))}},
2548 Point<dim>(0, +1.0 / sqrt(2.0), -1.0 / sqrt(2.0))}}}};
2550 return table[face_no][i];
2563 return table[face_no][i];
2591 return cross_product_2d(unit_tangential_vectors<dim>(face_no, 0));
2595 return cross_product_3d(unit_tangential_vectors<dim>(face_no, 0),
2596 unit_tangential_vectors<dim>(face_no, 1));
2627 const unsigned int line,
2628 const unsigned int face,
2629 const unsigned char combined_face_orientation,
2630 const bool line_orientation)
const
2634 static constexpr ::ndarray<bool, 2, 8> bool_table{
2635 {{{
true,
true,
false,
true,
false,
false,
true,
false}},
2636 {{
true,
true,
true,
false,
false,
false,
false,
true}}}};
2638 return (line_orientation ==
2639 bool_table[line / 2][combined_face_orientation]);
2644 static constexpr ::ndarray<unsigned int, 4, 3> combined_lines{
2645 {{{0, 0, 0}}, {{X, 0, 1}}, {{X, 0, X}}, {{X, X, X}}}};
2647 const auto combined_line = combined_lines[face][line];
2649 Assert(combined_line != X,
2651 "This function can only be called for following face-line "
2652 "combinations: (0,0), (0,1), (0,2), (1,1), (1,2), (2,1),"));
2654 static constexpr ::ndarray<bool, 2, 6> bool_table{
2655 {{{
false,
true,
false,
true,
false,
true}},
2656 {{
true,
false,
true,
false,
true,
false}}}};
2658 return (line_orientation ==
2659 bool_table[combined_line][combined_face_orientation]);
2671 template <
typename T>
2704 for (
unsigned int i = 0; i < n_vertices; ++i)
2707 if (i + 1 != n_vertices)
2711 out <<
"] is not a valid permutation of [";
2713 for (
unsigned int i = 0; i < n_vertices; ++i)
2716 if (i + 1 != n_vertices)
2720 out <<
"]." << std::endl;
2750 <<
"The reference-cell type used on this cell (" << arg1.
to_string()
2751 <<
") does not match the reference-cell type of the finite element "
2752 <<
"associated with this cell (" << arg2.to_string() <<
"). "
2753 <<
"Did you accidentally use simplex elements on hypercube meshes "
2754 <<
"(or the other way around), or are you using a mixed mesh and "
2755 <<
"assigned a simplex element to a hypercube cell (or the other "
2756 <<
"way around) via the active_fe_index?");
2761template <
typename T, std::
size_t N>
2764 const std::array<T, N> &vertices_1)
const
2767 ExcMessage(
"The number of array elements must be equal to or "
2768 "greater than the number of vertices of the cell "
2769 "referenced by this object."));
2782template <
typename T>
2789 ExcMessage(
"The number of array elements must be equal to "
2790 "the number of vertices of the cell "
2791 "referenced by this object."));
2793 ExcMessage(
"The number of array elements must be equal to "
2794 "the number of vertices of the cell "
2795 "referenced by this object."));
2798 for (
unsigned char o = 0; o < table.size(); ++o)
2801 for (
unsigned int j = 0; j < table[o].size(); ++j)
2802 match = (match && vertices_0[j] == vertices_1[table[o][j]]);
2809 return std::numeric_limits<unsigned char>::max();
2818 if (vertices_0[0] == vertices_1[0])
2832 return std::numeric_limits<unsigned char>::max();
2837template <
typename T, std::
size_t N>
2838inline std::array<T, N>
2841 const unsigned int orientation)
const
2844 ExcMessage(
"The number of array elements must be equal to or "
2845 "greater than the number of vertices of the cell "
2846 "referenced by this object."));
2856 std::array<T, N> temp;
2857 std::copy(permutation.begin(), permutation.end(), temp.begin());
2864template <
typename T>
2865boost::container::small_vector<T, 8>
2868 const unsigned char orientation)
const
2871 boost::container::small_vector<T, 8> result(
n_vertices());
2873 auto permute = [&](
const auto &table) {
2875 for (
unsigned int j = 0; j < table[orientation].size(); ++j)
2876 result[j] =
vertices[table[orientation][j]];
2906 const unsigned char orientation)
const
2922 constexpr std::array<unsigned char, 6> inverses{{0, 1, 2, 5, 4, 3}};
2923 return inverses[orientation];
2928 constexpr std::array<unsigned char, 8> inverses{
2929 {0, 1, 2, 7, 4, 5, 6, 3}};
2930 return inverses[orientation];
2936 return std::numeric_limits<unsigned char>::max();
2943ReferenceCell::vtk_lexicographic_to_node_index<0>(
2944 const std::array<unsigned, 0> &node_indices,
2945 const std::array<unsigned, 0> &nodes_per_direction,
2946 const bool legacy_format)
const;
2950ReferenceCell::vtk_lexicographic_to_node_index<1>(
2951 const std::array<unsigned, 1> &node_indices,
2952 const std::array<unsigned, 1> &nodes_per_direction,
2953 const bool legacy_format)
const;
2957ReferenceCell::vtk_lexicographic_to_node_index<2>(
2958 const std::array<unsigned, 2> &node_indices,
2959 const std::array<unsigned, 2> &nodes_per_direction,
2960 const bool legacy_format)
const;
2964ReferenceCell::vtk_lexicographic_to_node_index<3>(
2965 const std::array<unsigned, 3> &node_indices,
2966 const std::array<unsigned, 3> &nodes_per_direction,
2967 const bool legacy_format)
const;
ArrayView< typename std::remove_reference< typename std::iterator_traits< Iterator >::reference >::type, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
Abstract base class for mapping classes.
static constexpr Point< dim, Number > unit_vector(const unsigned int i)
ArrayView< const unsigned int > faces_for_given_vertex(const unsigned int vertex_index) const
unsigned int child_cell_on_face(const unsigned int face, const unsigned int subface, const unsigned char face_orientation=default_combined_face_orientation()) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
unsigned int n_vertices() const
boost::container::small_vector< T, 8 > permute_by_combined_orientation(const ArrayView< const T > &vertices, const unsigned char orientation) const
Quadrature< dim > get_gauss_type_quadrature(const unsigned n_points_1d) const
bool is_hyper_cube() const
unsigned int vtk_linear_type() const
void serialize(Archive &archive, const unsigned int)
unsigned char get_inverse_combined_orientation(const unsigned char orientation) const
unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const unsigned char face_orientation) const
static constexpr unsigned char default_combined_face_orientation()
unsigned char compute_orientation(const std::array< T, N > &vertices_0, const std::array< T, N > &vertices_1) const
bool standard_vs_true_line_orientation(const unsigned int line, const unsigned int face, const unsigned char face_orientation, const bool line_orientation) const
std::array< unsigned int, 2 > standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const
Point< dim > vertex(const unsigned int v) const
unsigned int standard_to_real_face_vertex(const unsigned int vertex, const unsigned int face, const unsigned char face_orientation) const
const Quadrature< dim > & get_nodal_type_quadrature() const
Point< dim > face_vertex_location(const unsigned int face, const unsigned int vertex) const
static constexpr ndarray< unsigned int, 1, 1 > vertex_vertex_permutations
std_cxx20::ranges::iota_view< unsigned int, unsigned int > isotropic_child_indices() const
std::array< unsigned int, 2 > standard_line_to_face_and_line_index(const unsigned int line) const
double d_linear_shape_function(const Point< dim > &xi, const unsigned int i) const
Tensor< 1, dim > unit_tangential_vectors(const unsigned int face_no, const unsigned int i) const
static constexpr ndarray< unsigned int, 2, 2 > line_vertex_permutations
unsigned int vtk_lexicographic_to_node_index(const std::array< unsigned, dim > &node_indices, const std::array< unsigned, dim > &nodes_per_direction, const bool legacy_format) const
std::array< T, N > permute_according_orientation(const std::array< T, N > &vertices, const unsigned int orientation) const
bool contains_point(const Point< dim > &p, const double tolerance=0) const
unsigned int n_isotropic_children() const
unsigned int gmsh_element_type() const
unsigned int n_face_orientations(const unsigned int face_no) const
static constexpr ndarray< unsigned int, 8, 4 > quadrilateral_vertex_permutations
unsigned int n_faces() const
unsigned int line_to_cell_vertices(const unsigned int line, const unsigned int vertex) const
constexpr ReferenceCell()
unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const
unsigned int unv_vertex_to_deal_vertex(const unsigned int vertex_n) const
unsigned int vtk_vertex_to_deal_vertex(const unsigned int vertex_index) const
std::unique_ptr< Mapping< dim, spacedim > > get_default_mapping(const unsigned int degree) const
Quadrature< dim > get_midpoint_quadrature() const
unsigned int vtk_quadratic_type() const
static constexpr unsigned char reversed_combined_line_orientation()
unsigned int n_lines() const
unsigned int vtk_lagrange_type() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > line_indices() const
unsigned int get_dimension() const
ReferenceCell face_reference_cell(const unsigned int face_no) const
unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const
static ReferenceCell n_vertices_to_type(const int dim, const unsigned int n_vertices)
unsigned int standard_to_real_face_line(const unsigned int line, const unsigned int face, const unsigned char face_orientation) const
static constexpr std::size_t memory_consumption()
friend std::istream & operator>>(std::istream &in, ReferenceCell &reference_cell)
const Mapping< dim, spacedim > & get_default_linear_mapping() const
std::string to_string() const
unsigned int face_to_cell_lines(const unsigned int face, const unsigned int line, const unsigned char face_orientation) const
constexpr bool operator!=(const ReferenceCell &type) const
unsigned char get_combined_orientation(const ArrayView< const T > &vertices_0, const ArrayView< const T > &vertices_1) const
constexpr bool operator==(const ReferenceCell &type) const
Point< dim > barycenter() const
Tensor< 1, dim > unit_normal_vectors(const unsigned int face_no) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices() const
friend std::ostream & operator<<(std::ostream &out, const ReferenceCell &reference_cell)
Point< dim > closest_point(const Point< dim > &p) const
static constexpr ndarray< unsigned int, 6, 3 > triangle_vertex_permutations
Tensor< 1, dim > d_linear_shape_function_gradient(const Point< dim > &xi, const unsigned int i) const
virtual void print_info(std::ostream &out) const override
const ReferenceCell entity_type
virtual ~NoPermutation() noexcept override=default
const ArrayView< const T > vertices_0
const ArrayView< const T > vertices_1
NoPermutation(const ReferenceCell &entity_type, const ArrayView< const T > &vertices_0, const ArrayView< const T > &vertices_1)
#define DEAL_II_DEPRECATED
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNonMatchingReferenceCellTypes(ReferenceCell arg1, ReferenceCell arg2)
#define Assert(cond, exc)
#define DeclException2(Exception2, type1, type2, outsequence)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define DEAL_II_ASSERT_UNREACHABLE()
#define DEAL_II_NOT_IMPLEMENTED()
constexpr const ReferenceCell Tetrahedron
constexpr const ReferenceCell Quadrilateral
constexpr const ReferenceCell Wedge
constexpr const ReferenceCell Pyramid
constexpr const ReferenceCell Invalid
constexpr const ReferenceCell Triangle
constexpr const ReferenceCell Hexahedron
constexpr const ReferenceCell & get_hypercube()
constexpr const ReferenceCell Vertex
constexpr const ReferenceCell Line
constexpr const ReferenceCell & get_simplex()
constexpr ReferenceCell make_reference_cell_from_int(const std::uint8_t kind)
std::tuple< bool, bool, bool > split_face_orientation(const unsigned char combined_face_orientation)
static const unsigned int invalid_unsigned_int
boost::integer_range< IncrementableType > iota_view
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > sqrt(const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > pow(const ::VectorizedArray< Number, width > &, const Number p)
typename internal::ndarray::HelperArray< T, Ns... >::type ndarray
std::istream & operator>>(std::istream &in, ReferenceCell &reference_cell)
std::ostream & operator<<(std::ostream &out, const ReferenceCell &reference_cell)
static unsigned int child_cell_on_face(const RefinementCase< dim > &ref_case, const unsigned int face, const unsigned int subface, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false, const RefinementCase< dim - 1 > &face_refinement_case=RefinementCase< dim - 1 >::isotropic_refinement)
static unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
static std::array< unsigned int, 2 > standard_hex_line_to_quad_line_index(const unsigned int line)
static unsigned int standard_to_real_face_line(const unsigned int line, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
static std::array< unsigned int, 2 > standard_hex_vertex_to_quad_vertex_index(const unsigned int vertex)
static std::array< unsigned int, 2 > standard_quad_vertex_to_line_vertex_index(const unsigned int vertex)
static double d_linear_shape_function(const Point< dim > &xi, const unsigned int i)
static unsigned int face_to_cell_lines(const unsigned int face, const unsigned int line, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
static Tensor< 1, dim > d_linear_shape_function_gradient(const Point< dim > &xi, const unsigned int i)