DOLFINx
DOLFINx 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 <cstdint>
10#include <dolfinx/mesh/cell_types.h>
11#include <vector>
12#include <xtensor/xtensor.hpp>
13
17{
18/*
19 The basix 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 2-----------3 2-----7-----3 2--10--11---3
39 | | | | | |
40 | | | | 7 14 15 9
41 | | 5 8 6 | |
42 | | | | 6 12 13 8
43 | | | | | |
44 0-----------1 --> u 0-----4-----1 0---4---5---1
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 w
66 6----------7 6----19----7
67 /| ^ v /| /| /|
68 / | | / / | 17 | 25 18|
69 / | | / / | / 14 24 / 15
70 4----------5 | 4----16----5 |
71 | | +--|---|--> u |22 | 26 | 23|
72 | 2------+---3 | 2----13+---3
73 | / | / 10 / 21 12 /
74 | / | / | 9 20 | 11
75 |/ |/ |/ |/
76 0----------1 0-----8----1
77
78*/
79
88std::vector<std::uint8_t> perm_vtk(mesh::CellType type, int num_nodes);
89
98std::vector<std::uint8_t> perm_gmsh(mesh::CellType type, int num_nodes);
99
105std::vector<std::uint8_t> transpose(const std::vector<std::uint8_t>& map);
106
114xt::xtensor<std::int64_t, 2>
115compute_permutation(const xt::xtensor<std::int64_t, 2>& cells,
116 const std::vector<std::uint8_t>& p);
117
122std::int8_t get_vtk_cell_type(mesh::CellType cell, int dim);
123
124} // namespace dolfinx::io::cells
Functions for the re-ordering of input mesh topology to the DOLFINx ordering, and transpose orderings...
Definition: cells.h:17
std::vector< std::uint8_t > transpose(const std::vector< std::uint8_t > &map)
Compute the transpose of a re-ordering map.
Definition: cells.cpp:359
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:291
xt::xtensor< std::int64_t, 2 > compute_permutation(const xt::xtensor< std::int64_t, 2 > &cells, const std::vector< std::uint8_t > &p)
Permute cell topology by applying a permutation array for each cell.
Definition: cells.cpp:368
std::int8_t get_vtk_cell_type(mesh::CellType cell, int dim)
Get VTK cell identifier.
Definition: cells.cpp:383
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:326
CellType
Cell type identifier.
Definition: cell_types.h:22