57 VTK_QUADRATIC_EDGE = 21,
58 VTK_QUADRATIC_TRIANGLE = 22,
59 VTK_QUADRATIC_QUAD = 23,
60 VTK_QUADRATIC_TETRA = 24,
61 VTK_QUADRATIC_HEXAHEDRON = 25,
62 VTK_QUADRATIC_WEDGE = 26,
63 VTK_QUADRATIC_PYRAMID = 27,
65 VTK_LAGRANGE_CURVE = 68,
66 VTK_LAGRANGE_TRIANGLE = 69,
67 VTK_LAGRANGE_QUADRILATERAL = 70,
68 VTK_LAGRANGE_TETRAHEDRON = 71,
69 VTK_LAGRANGE_HEXAHEDRON = 72,
70 VTK_LAGRANGE_WEDGE = 73,
71 VTK_LAGRANGE_PYRAMID = 74,
120 template <
int dim,
int spacedim>
121 std::unique_ptr<Mapping<dim, spacedim>>
127 return std::make_unique<MappingQ<dim, spacedim>>(degree);
132 return std::make_unique<MappingFE<dim, spacedim>>(
135 return std::make_unique<MappingFE<dim, spacedim>>(
142 return std::make_unique<MappingQ<dim, spacedim>>(degree);
147 template <
int dim,
int spacedim>
265 constexpr std::array<unsigned int, 4> exodus_to_deal{{0, 1, 3, 2}};
266 return exodus_to_deal[vertex_n];
272 constexpr std::array<unsigned int, 8> exodus_to_deal{
273 {0, 1, 3, 2, 4, 5, 7, 6}};
274 return exodus_to_deal[vertex_n];
278 constexpr std::array<unsigned int, 6> exodus_to_deal{
280 return exodus_to_deal[vertex_n];
284 constexpr std::array<unsigned int, 5> exodus_to_deal{{0, 1, 3, 2, 4}};
285 return exodus_to_deal[vertex_n];
310 constexpr std::array<unsigned int, 4> exodus_to_deal{{2, 1, 3, 0}};
311 return exodus_to_deal[face_n];
315 constexpr std::array<unsigned int, 4> exodus_to_deal{{1, 3, 2, 0}};
316 return exodus_to_deal[face_n];
320 constexpr std::array<unsigned int, 6> exodus_to_deal{
322 return exodus_to_deal[face_n];
326 constexpr std::array<unsigned int, 6> exodus_to_deal{{3, 4, 2, 0, 1}};
327 return exodus_to_deal[face_n];
331 constexpr std::array<unsigned int, 5> exodus_to_deal{{3, 2, 4, 1, 0}};
332 return exodus_to_deal[face_n];
362 constexpr std::array<unsigned int, 4> unv_to_deal{{1, 0, 2, 3}};
363 return unv_to_deal[vertex_n];
367 constexpr std::array<unsigned int, 8> unv_to_deal{
368 {6, 7, 5, 4, 2, 3, 1, 0}};
369 return unv_to_deal[vertex_n];
385 return VTKCellType::VTK_VERTEX;
387 return VTKCellType::VTK_LINE;
389 return VTKCellType::VTK_TRIANGLE;
391 return VTKCellType::VTK_QUAD;
393 return VTKCellType::VTK_TETRA;
395 return VTKCellType::VTK_PYRAMID;
397 return VTKCellType::VTK_WEDGE;
399 return VTKCellType::VTK_HEXAHEDRON;
401 return VTKCellType::VTK_INVALID;
406 return VTKCellType::VTK_INVALID;
417 return VTKCellType::VTK_VERTEX;
419 return VTKCellType::VTK_QUADRATIC_EDGE;
421 return VTKCellType::VTK_QUADRATIC_TRIANGLE;
423 return VTKCellType::VTK_QUADRATIC_QUAD;
425 return VTKCellType::VTK_QUADRATIC_TETRA;
427 return VTKCellType::VTK_QUADRATIC_PYRAMID;
429 return VTKCellType::VTK_QUADRATIC_WEDGE;
431 return VTKCellType::VTK_QUADRATIC_HEXAHEDRON;
433 return VTKCellType::VTK_INVALID;
438 return VTKCellType::VTK_INVALID;
449 return VTKCellType::VTK_VERTEX;
451 return VTKCellType::VTK_LAGRANGE_CURVE;
453 return VTKCellType::VTK_LAGRANGE_TRIANGLE;
455 return VTKCellType::VTK_LAGRANGE_QUADRILATERAL;
457 return VTKCellType::VTK_LAGRANGE_TETRAHEDRON;
459 return VTKCellType::VTK_LAGRANGE_PYRAMID;
461 return VTKCellType::VTK_LAGRANGE_WEDGE;
463 return VTKCellType::VTK_LAGRANGE_HEXAHEDRON;
465 return VTKCellType::VTK_INVALID;
470 return VTKCellType::VTK_INVALID;
477 ReferenceCell::vtk_lexicographic_to_node_index<0>(
478 const std::array<unsigned, 0> &,
479 const std::array<unsigned, 0> &,
490 ReferenceCell::vtk_lexicographic_to_node_index<1>(
491 const std::array<unsigned, 1> &,
492 const std::array<unsigned, 1> &,
507 ReferenceCell::vtk_lexicographic_to_node_index<2>(
508 const std::array<unsigned, 2> &node_indices,
509 const std::array<unsigned, 2> &nodes_per_direction,
514 const unsigned int i = node_indices[0];
515 const unsigned int j = node_indices[1];
517 const bool ibdy = (i == 0 || i == nodes_per_direction[0]);
518 const bool jbdy = (j == 0 || j == nodes_per_direction[1]);
520 const int nbdy = (ibdy ? 1 : 0) + (jbdy ? 1 : 0);
524 return (i != 0u ? (j != 0u ? 2 : 1) : (j != 0u ? 3 : 0));
534 nodes_per_direction[0] - 1 + nodes_per_direction[1] - 1 :
542 (i != 0u ? nodes_per_direction[0] - 1 :
543 2 * (nodes_per_direction[0] - 1) +
544 nodes_per_direction[1] - 1) +
549 offset += 2 * (nodes_per_direction[0] - 1 + nodes_per_direction[1] - 1);
551 return offset + (i - 1) + (nodes_per_direction[0] - 1) * ((j - 1));
568 ReferenceCell::vtk_lexicographic_to_node_index<3>(
569 const std::array<unsigned, 3> &node_indices,
570 const std::array<unsigned, 3> &nodes_per_direction,
571 const bool legacy_format)
const
575 const unsigned int i = node_indices[0];
576 const unsigned int j = node_indices[1];
577 const unsigned int k = node_indices[2];
579 const bool ibdy = (i == 0 || i == nodes_per_direction[0]);
580 const bool jbdy = (j == 0 || j == nodes_per_direction[1]);
581 const bool kbdy = (k == 0 || k == nodes_per_direction[2]);
583 const int nbdy = (ibdy ? 1 : 0) + (jbdy ? 1 : 0) + (kbdy ? 1 : 0);
587 return (i != 0u ? (j != 0u ? 2 : 1) : (j != 0u ? 3 : 0)) +
598 nodes_per_direction[0] - 1 + nodes_per_direction[1] - 1 :
600 (k != 0u ? 2 * (nodes_per_direction[0] - 1 +
601 nodes_per_direction[1] - 1) :
608 (i != 0u ? nodes_per_direction[0] - 1 :
609 2 * (nodes_per_direction[0] - 1) +
610 nodes_per_direction[1] - 1) +
611 (k != 0u ? 2 * (nodes_per_direction[0] - 1 +
612 nodes_per_direction[1] - 1) :
618 4 * (nodes_per_direction[0] - 1) + 4 * (nodes_per_direction[1] - 1);
621 (nodes_per_direction[2] - 1) *
622 (i != 0u ? (j != 0u ? 3 : 1) : (j != 0u ? 2 : 0)) +
626 (nodes_per_direction[2] - 1) *
627 (i != 0u ? (j != 0u ? 2 : 1) : (j != 0u ? 3 : 0)) +
631 offset += 4 * (nodes_per_direction[0] - 1 + nodes_per_direction[1] - 1 +
632 nodes_per_direction[2] - 1);
637 return (j - 1) + ((nodes_per_direction[1] - 1) * (k - 1)) +
638 (i != 0u ? (nodes_per_direction[1] - 1) *
639 (nodes_per_direction[2] - 1) :
643 offset += 2 * (nodes_per_direction[1] - 1) * (nodes_per_direction[2] - 1);
646 return (i - 1) + ((nodes_per_direction[0] - 1) * (k - 1)) +
647 (j != 0u ? (nodes_per_direction[2] - 1) *
648 (nodes_per_direction[0] - 1) :
652 offset += 2 * (nodes_per_direction[2] - 1) * (nodes_per_direction[0] - 1);
654 return (i - 1) + ((nodes_per_direction[0] - 1) * (j - 1)) +
656 (nodes_per_direction[0] - 1) * (nodes_per_direction[1] - 1) :
662 offset += 2 * ((nodes_per_direction[1] - 1) * (nodes_per_direction[2] - 1) +
663 (nodes_per_direction[2] - 1) * (nodes_per_direction[0] - 1) +
664 (nodes_per_direction[0] - 1) * (nodes_per_direction[1] - 1));
665 return offset + (i - 1) +
666 (nodes_per_direction[0] - 1) *
667 ((j - 1) + (nodes_per_direction[1] - 1) * ((k - 1)));
693 static constexpr std::array<unsigned int, 4> index_translation_table =
695 return index_translation_table[vertex_index];
701 static constexpr std::array<unsigned int, 5> index_translation_table =
703 return index_translation_table[vertex_index];
709 static constexpr std::array<unsigned int, 8> index_translation_table =
710 {{0, 1, 3, 2, 4, 5, 7, 6}};
711 return index_translation_table[vertex_index];
816 std::pair<Point<dim>,
double>
827 const double t = ((x1 - x0) * (p - x0)) / ((x1 - x0).norm_square());
834 const auto p2 = x0 + t * (x1 - x0);
835 return std::make_pair(p2, p2.distance_square(p));
843 std::pair<Point<dim>,
double>
844 project_to_quad(
const std::array<
Point<dim>, 3> & ,
850 std::numeric_limits<double>::signaling_NaN());
871 std::pair<Point<3>,
double>
883 std::array<Point<3>, 3> shifted_vertices =
vertices;
885 for (
Point<3> &shifted_vertex : shifted_vertices)
886 shifted_vertex +=
shift;
895 const Point<3> vertex = shifted_vertices[0];
905 e0 = shifted_vertices[1] - shifted_vertices[0];
906 e1 = shifted_vertices[2] - shifted_vertices[0];
911 e0 = shifted_vertices[1] - shifted_vertices[0];
912 e1 = shifted_vertices[2] - shifted_vertices[0];
917 const double c0 = e0 * (shifted_p - vertex) / e0.
norm_square();
918 const double c1 = e1 * (shifted_p - vertex) / e1.norm_square();
919 const Point<3> projected_shifted_p = vertex + c0 * e0 + c1 * e1;
921 bool in_quad =
false;
925 for (
unsigned int i = 0; i < 3; ++i)
926 shifted_vertex_matrix[i] = shifted_vertices[i];
929 bool is_convex_combination =
true;
930 for (
unsigned int i = 0; i < 3; ++i)
931 is_convex_combination = is_convex_combination &&
932 (0.0 <= combination_coordinates[i]) &&
933 (combination_coordinates[i] <= 1.0);
934 in_quad = is_convex_combination;
937 in_quad = (0.0 <= c0 && c0 <= 1.0 && 0.0 <= c1 && c1 <= 1.0);
940 return std::make_pair(projected_shifted_p -
shift,
961 return project_to_line(vertex<dim>(0), vertex<dim>(1), p).first;
966 unsigned int closest_vertex_no = 0;
967 double closest_vertex_distance_square = vertex<dim>(0).
distance_square(p);
968 for (
unsigned int i = 1; i <
n_vertices(); ++i)
970 const double new_vertex_distance_square =
971 vertex<dim>(i).distance_square(p);
972 if (new_vertex_distance_square < closest_vertex_distance_square)
974 closest_vertex_no = i;
975 closest_vertex_distance_square = new_vertex_distance_square;
982 for (
const unsigned int face_no :
988 auto pair = project_to_line(
v0,
v1, p);
989 if (pair.second < min_distance_square)
992 min_distance_square = pair.second;
1009 const std::array<unsigned int, 5> all_pyramid_faces{{0, 1, 2, 3, 4}};
1013 for (
const unsigned int face_no : faces)
1018 std::array<Point<dim>, 3>
vertices;
1019 for (
unsigned int vertex_no = 0; vertex_no < 3; ++vertex_no)
1023 auto pair = project_to_quad(
vertices, p, face_cell);
1024 if (pair.second < min_distance_square)
1026 result = pair.first;
1027 min_distance_square = pair.second;
1031 for (
const unsigned int face_no :
1035 for (
const unsigned int face_line_no : face_cell.line_indices())
1037 const auto cell_line_no =
1045 auto pair = project_to_line(
v0,
v1, p);
1046 if (pair.second < min_distance_square)
1048 result = pair.first;
1049 min_distance_square = pair.second;
1064 constexpr
unsigned int x_index = 0;
1065 constexpr
unsigned int y_index = (dim >= 2 ? 1 : 0);
1066 constexpr
unsigned int z_index = (dim >= 3 ? 2 : 0);
1077 for (
unsigned int d = 0;
d < dim; ++
d)
1078 result[
d] = std::clamp(result[
d], 0.0, 1.0);
1082 result[x_index] = std::clamp(result[x_index], 0.0, 1.0);
1084 std::clamp(result[y_index], 0.0, 1.0 - result[x_index]);
1087 result[x_index] = std::clamp(result[x_index], 0.0, 1.0);
1089 std::clamp(result[y_index], 0.0, 1.0 - result[x_index]);
1091 std::clamp(result[z_index],
1093 1.0 - result[x_index] - result[y_index]);
1097 result[x_index] = std::clamp(result[x_index], 0.0, 1.0);
1099 std::clamp(result[y_index], 0.0, 1.0 - result[x_index]);
1100 result[z_index] = std::clamp(result[z_index], 0.0, 1.0);
1104 result[x_index] = std::clamp(result[x_index], -1.0, 1.0);
1105 result[y_index] = std::clamp(result[y_index], -1.0, 1.0);
1108 const auto x_abs = std::abs(result[x_index]);
1109 const auto y_abs = std::abs(result[y_index]);
1112 result[z_index] = std::clamp(result[z_index], 0.0, 1.0 - x_abs);
1114 result[z_index] = std::clamp(result[z_index], 0.0, 1.0 - y_abs);
1167 "The reference cell kind just read does not correspond to one of the "
1168 "valid choices. There must be an error."));
1180 #include "reference_cell.inst"
OutputOperator< VectorType > & operator<<(OutputOperator< VectorType > &out, unsigned int step)
numbers::NumberTraits< Number >::real_type distance_square(const Point< dim, Number > &p) const
std::istream & operator>>(std::istream &in, Point< dim, Number > &p)
ArrayView< const unsigned int > faces_for_given_vertex(const unsigned int vertex_index) const
unsigned int n_vertices() const
Quadrature< dim > get_gauss_type_quadrature(const unsigned n_points_1d) const
bool is_hyper_cube() const
unsigned int vtk_linear_type() 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()
const Quadrature< dim > & get_nodal_type_quadrature() const
static constexpr ndarray< unsigned int, 2, 2 > line_vertex_permutations
bool contains_point(const Point< dim > &p, const double tolerance=0) const
unsigned int gmsh_element_type() 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
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
unsigned int vtk_quadratic_type() const
unsigned int vtk_lagrange_type() 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
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
Point< dim > closest_point(const Point< dim > &p) const
static constexpr ndarray< unsigned int, 6, 3 > triangle_vertex_permutations
constexpr DEAL_II_HOST SymmetricTensor< 2, dim, Number > invert(const SymmetricTensor< 2, dim, Number > &t)
constexpr numbers::NumberTraits< Number >::real_type norm_square() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > epsilon(const Tensor< 2, dim, Number > &Grad_u)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
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 Vertex
constexpr const ReferenceCell Line
static const unsigned int invalid_unsigned_int
typename internal::ndarray::HelperArray< T, Ns... >::type ndarray