DOLFIN-X
DOLFIN-X C++ interface
cells.h
1 // Copyright (C) 2019 Jorgen S. Dokken
2 //
3 // This file is part of DOLFINX (https://www.fenicsproject.org)
4 //
5 // SPDX-License-Identifier: LGPL-3.0-or-later
6 
7 #pragma once
8 
9 #include <Eigen/Dense>
10 #include <cstdint>
11 #include <dolfinx/mesh/cell_types.h>
12 #include <vector>
13 
17 {
18 /*
19  The FIAT ordering is used for the geometry nodes, and is shown below
20  for a range of cell types.
21 
22  Triangle: Triangle6: Triangle10:
23  v
24  ^
25  |
26  2 2 2
27  |`\ |`\ | \
28  | `\ | `\ 6 4
29  | `\ 4 `3 | \
30  | `\ | `\ 5 9 3
31  | `\ | `\ | \
32  0----------1 --> u 0-----5----1 0---7---8---1
33 
34  Quadrilateral: Quadrilateral9: Quadrilateral16:
35  v
36  ^
37  |
38  1-----------3 1-----7-----4 1---9--13---5
39  | | | | | |
40  | | | | 3 11 15 7
41  | | 2 8 5 | |
42  | | | | 2 10 14 6
43  | | | | | |
44  0-----------2 --> u 0-----6-----3 0---8--12---4
45 
46  Tetrahedron: Tetrahedron10: Tetrahedron20
47  v
48  /
49  2 2 2
50  ,/|`\ ,/|`\ ,/|`\
51  ,/ | `\ ,/ | `\ 13 | `9
52  ,/ '. `\ ,8 '. `6 ,/ 4 `\
53  ,/ | `\ ,/ 4 `\ 12 19 | `8
54  ,/ | `\ ,/ | `\ ,/ | `\
55  0-----------'.--------1 -> u 0--------9--'.--------1 0-----14----'.--15----1
56  `\. | ,/ `\. | ,/ `\. 17 | 16 ,/
57  `\. | ,/ `\. | ,5 10. 18 5 ,6
58  `\. '. ,/ `7. '. ,/ `\. '. 7
59  `\. |/ `\. |/ 11. |/
60  `3 `3 `3
61  `\.
62  w
63 
64  Hexahedron: Hexahedron27:
65  v
66  2----------6 3----21----12
67  |\ ^ |\ |\ |\
68  | \ | | \ | 5 23 | 14
69  | \ | | \ 6 \ 24 15 \
70  | 3------+---7 | 4----22+---13
71  | | +-- |-- | -> u | 8 | 26 | 17|
72  0---+---\--4 | 0---+18----9 |
73  \ | \ \ | \ 7 25\ 16
74  \ | \ \ | 2 | 20 11|
75  \| w \| \| \|
76  1----------5 1----19----10
77 */
78 
87 std::vector<std::uint8_t> perm_vtk(mesh::CellType type, int num_nodes);
88 
97 std::vector<std::uint8_t> perm_gmsh(mesh::CellType type, int num_nodes);
98 
104 std::vector<std::uint8_t> transpose(const std::vector<std::uint8_t>& map);
105 
113 Eigen::Array<std::int64_t, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>
115  const Eigen::Ref<const Eigen::Array<
116  std::int64_t, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>>& cells,
117  const std::vector<std::uint8_t>& p);
118 
119 } // namespace dolfinx::io::cells
dolfinx::io::cells::perm_gmsh
std::vector< std::uint8_t > perm_gmsh(mesh::CellType type, int num_nodes)
Permutation array to map from Gmsh to DOLFINX node ordering.
Definition: cells.cpp:315
dolfinx::mesh::CellType
CellType
Cell type identifier.
Definition: cell_types.h:21
dolfinx::io::cells::perm_vtk
std::vector< std::uint8_t > perm_vtk(mesh::CellType type, int num_nodes)
Permutation array to map from VTK to DOLFINX node ordering.
Definition: cells.cpp:283
dolfinx::io::cells
Functions for the re-ordering of input mesh topology to the DOLFINX ordering, and transpose orderings...
Definition: cells.h:17
dolfinx::io::cells::compute_permutation
Eigen::Array< std::int64_t, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > compute_permutation(const Eigen::Ref< const Eigen::Array< std::int64_t, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &cells, const std::vector< std::uint8_t > &p)
Permute cell topology by applying a permutation array for each cell.
Definition: cells.cpp:357
dolfinx::io::cells::transpose
std::vector< std::uint8_t > transpose(const std::vector< std::uint8_t > &map)
Compute the transpose of a re-ordering map.
Definition: cells.cpp:348