Skip to content

Renumber simplex reference cells #1028

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Draft
mscroggs wants to merge 9 commits into
base: main
Choose a base branch
from
Draft

Renumber simplex reference cells #1028

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
9 commits
Select commit Hold shift + click to select a range
2793f55
add test that asserts reference cell numbering follows DefElement con…
mscroggs Jun 17, 2026
cfc46f0
update convention for triangle cell
mscroggs Jun 17, 2026
5bbb919
update convention for tets
mscroggs Jun 17, 2026
bbf7199
updating tests
mscroggs Jun 17, 2026
cb9ea4b
update more tests
mscroggs Jun 17, 2026
86691a3
fix jacobians
mscroggs Jun 18, 2026
d38f6d6
update readme images
mscroggs Jun 18, 2026
e43240e
update demo
mscroggs Jun 18, 2026
78e546c
ffc branch
mscroggs Jun 18, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion .github/workflows/fenicsx-refs.env
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
ufl_repository=FEniCS/ufl
ufl_ref=main
ffcx_repository=FEniCS/ffcx
ffcx_ref=main
ffcx_ref=mscroggs/reorder-simplices
dolfinx_repository=FEniCS/dolfinx
dolfinx_ref=main
32 changes: 16 additions & 16 deletions cpp/basix/cell.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -63,7 +63,7 @@ std::vector<std::vector<std::vector<int>>> cell::topology(cell::type celltype)
// Vertices
t[0] = {{0}, {1}, {2}};
// Edges
t[1] = {{1, 2}, {0, 2}, {0, 1}};
t[1] = {{0, 1}, {0, 2}, {1, 2}};
// Cell
t[2] = {{0, 1, 2}};
return t;
Expand All @@ -85,9 +85,9 @@ std::vector<std::vector<std::vector<int>>> cell::topology(cell::type celltype)
// Vertices
t[0] = {{0}, {1}, {2}, {3}};
// Edges
t[1] = {{2, 3}, {1, 3}, {1, 2}, {0, 3}, {0, 2}, {0, 1}};
t[1] = {{0, 1}, {0, 2}, {0, 3}, {1, 2}, {1, 3}, {2, 3}};
// Faces
t[2] = {{1, 2, 3}, {0, 2, 3}, {0, 1, 3}, {0, 1, 2}};
t[2] = {{0, 1, 2}, {0, 1, 3}, {0, 2, 3}, {1, 2, 3}};
// Cell
t[3] = {{0, 1, 2, 3}};
return t;
Expand Down Expand Up @@ -159,9 +159,9 @@ cell::sub_entity_connectivity(cell::type celltype)
{
std::vector<std::vector<std::vector<std::vector<int>>>> t(3);
// Vertices
t[0] = {{{0}, {1, 2}, {0}}, {{1}, {0, 2}, {0}}, {{2}, {0, 1}, {0}}};
t[0] = {{{0}, {0, 1}, {0}}, {{1}, {0, 2}, {0}}, {{2}, {1, 2}, {0}}};
// Edges
t[1] = {{{1, 2}, {0}, {0}}, {{0, 2}, {1}, {0}}, {{0, 1}, {2}, {0}}};
t[1] = {{{0, 1}, {0}, {0}}, {{0, 2}, {1}, {0}}, {{1, 2}, {2}, {0}}};
// Face
t[2] = {{{0, 1, 2}, {0, 1, 2}, {0}}};
return t;
Expand All @@ -187,21 +187,21 @@ cell::sub_entity_connectivity(cell::type celltype)
{
std::vector<std::vector<std::vector<std::vector<int>>>> t(4);
// Vertices
t[0] = {{{0}, {3, 4, 5}, {1, 2, 3}, {0}},
{{1}, {1, 2, 5}, {0, 2, 3}, {0}},
{{2}, {0, 2, 4}, {0, 1, 3}, {0}},
{{3}, {0, 1, 3}, {0, 1, 2}, {0}}};
t[0] = {{{0}, {0, 1, 2}, {0, 1, 2}, {0}},
{{1}, {0, 3, 4}, {0, 1, 3}, {0}},
{{2}, {1, 3, 5}, {0, 2, 3}, {0}},
{{3}, {2, 4, 5}, {1, 2, 3}, {0}}};
// Edges
t[1] = {
{{2, 3}, {0}, {0, 1}, {0}}, {{1, 3}, {1}, {0, 2}, {0}},
{{1, 2}, {2}, {0, 3}, {0}}, {{0, 3}, {3}, {1, 2}, {0}},
{{0, 2}, {4}, {1, 3}, {0}}, {{0, 1}, {5}, {2, 3}, {0}},
{{0, 1}, {0}, {0, 1}, {0}}, {{0, 2}, {1}, {0, 2}, {0}},
{{0, 3}, {2}, {1, 2}, {0}}, {{1, 2}, {3}, {0, 3}, {0}},
{{1, 3}, {4}, {1, 3}, {0}}, {{2, 3}, {5}, {2, 3}, {0}},
};
// Faces
t[2] = {{{1, 2, 3}, {0, 1, 2}, {0}, {0}},
{{0, 2, 3}, {0, 3, 4}, {1}, {0}},
{{0, 1, 3}, {1, 3, 5}, {2}, {0}},
{{0, 1, 2}, {2, 4, 5}, {3}, {0}}};
t[2] = {{{0, 1, 2}, {0, 1, 3}, {0}, {0}},
{{0, 1, 3}, {0, 2, 4}, {1}, {0}},
{{0, 2, 3}, {1, 2, 5}, {2}, {0}},
{{1, 2, 3}, {3, 4, 5}, {3}, {0}}};
// Volume
t[3] = {{{0, 1, 2, 3}, {0, 1, 2, 3, 4, 5}, {0, 1, 2, 3}, {0}}};
return t;
Expand Down
36 changes: 20 additions & 16 deletions cpp/basix/dof-transformations.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -85,13 +85,15 @@ mapinfo_t<T> get_mapinfo(cell::type cell_type)
{
mapinfo_t<T> mapinfo;
auto& data = mapinfo.try_emplace(cell::type::interval).first->second;
auto map = [](auto pt) -> std::array<T, 3> { return {pt[1], pt[0], 0}; };
auto map = [](auto pt) -> std::array<T, 3> {
return {1 - pt[0], pt[1], 0};
};
mdex::mdarray<T, md::extents<std::size_t, 2, 2>> J(
md::extents<std::size_t, 2, 2>{}, {0., 1., 1., 0.});
md::extents<std::size_t, 2, 2>{}, {-1., 0., 0., 1.});

T detJ = -1;
mdex::mdarray<T, md::extents<std::size_t, 2, 2>> K(
md::extents<std::size_t, 2, 2>{}, {0., 1., 1., 0.});
md::extents<std::size_t, 2, 2>{}, {-1., 0., 0., 1.});
data.push_back(std::tuple(map, J, detJ, K));
return mapinfo;
}
Expand All @@ -115,45 +117,48 @@ mapinfo_t<T> get_mapinfo(cell::type cell_type)
mapinfo_t<T> mapinfo;
{
auto& data = mapinfo.try_emplace(cell::type::interval).first->second;
auto map
= [](auto pt) -> std::array<T, 3> { return {pt[0], pt[2], pt[1]}; };
auto map = [](auto pt) -> std::array<T, 3> {
return {1 - pt[0], pt[1], pt[2]};
};
mdex::mdarray<T, md::extents<std::size_t, 3, 3>> J(
md::extents<std::size_t, 3, 3>{},
{1., 0., 0., 0., 0., 1., 0., 1., 0.});
{-1., 0., 0., 0., 1., 0., 0., 0., 1.});

T detJ = -1.0;
mdex::mdarray<T, md::extents<std::size_t, 3, 3>> K(
md::extents<std::size_t, 3, 3>{},
{1., 0., 0., 0., 0., 1., 0., 1., 0.});
{-1., 0., 0., 0., 1., 0., 0., 0., 1.});
data.push_back(std::tuple(map, J, detJ, K));
}

{
auto& data = mapinfo.try_emplace(cell::type::triangle).first->second;
{
auto map
= [](auto pt) -> std::array<T, 3> { return {pt[2], pt[0], pt[1]}; };
auto map = [](auto pt) -> std::array<T, 3> {
return {1 - pt[0] - pt[1], pt[0], pt[2]};
};
mdex::mdarray<T, md::extents<std::size_t, 3, 3>> J(
md::extents<std::size_t, 3, 3>{},
{0., 1., 0., 0., 0., 1., 1., 0., 0.});
{0., 1., 0., -1., -1., 0., 0., 0., 1.});

T detJ = 1.0;
mdex::mdarray<T, md::extents<std::size_t, 3, 3>> K(
md::extents<std::size_t, 3, 3>{},
{0., 0., 1., 1., 0., 0., 0., 1., 0.});
{-1., -1., 0., 1., 0., 0., 0., 0., 1.});
data.push_back(std::tuple(map, J, detJ, K));
}
{
auto map
= [](auto pt) -> std::array<T, 3> { return {pt[0], pt[2], pt[1]}; };
auto map = [](auto pt) -> std::array<T, 3> {
return {pt[1], pt[0], pt[2]};
};
mdex::mdarray<T, md::extents<std::size_t, 3, 3>> J(
md::extents<std::size_t, 3, 3>{},
{1., 0., 0., 0., 0., 1., 0., 1., 0.});
{0., 1., 0., 1., 0., 0., 0., 0., 1.});

T detJ = -1.0;
mdex::mdarray<T, md::extents<std::size_t, 3, 3>> K(
md::extents<std::size_t, 3, 3>{},
{1., 0., 0., 0., 0., 1., 0., 1., 0.});
{0., 1., 0., 1., 0., 0., 0., 0., 1.});
data.push_back(std::tuple(map, J, detJ, K));
}
}
Expand Down Expand Up @@ -433,7 +438,6 @@ std::pair<std::vector<T>, std::array<std::size_t, 2>> compute_transformation(
for (std::size_t k1 = 0; k1 < result.extent(0); ++k1)
tabulated_data(k0, k1, j) = result(k1, k0);
}

// push forward
mdarray_t<T, 3> pushed_data(tabulated_data.extents());
{
Expand Down
24 changes: 12 additions & 12 deletions cpp/basix/finite-element.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -561,15 +561,15 @@ std::vector<int> basix::lex_dof_ordering(element::family family, cell::type cell
{
int n = degree - 1;
for (int i = 0; i < n; ++i)
perm.push_back(3 + 2 * n + i);
perm.push_back(3 + i);
perm.push_back(1);
int dof = 3 + 3 * n;
for (int i = 0; i < n; ++i)
{
perm.push_back(3 + n + i);
for (int j = 0; j < n - 1 - i; ++j)
perm.push_back(dof++);
perm.push_back(3 + i);
perm.push_back(3 + 2 * n + i);
}
perm.push_back(2);
}
Expand Down Expand Up @@ -648,30 +648,30 @@ std::vector<int> basix::lex_dof_ordering(element::family family, cell::type cell
int face3 = 4 + 6 * n + n * (n - 1) * 3 / 2;
int interior = 4 + 6 * n + n * (n - 1) * 2;
for (int i = 0; i < n; ++i)
perm.push_back(4 + 5 * n + i);
perm.push_back(4 + i);
perm.push_back(1);
for (int i = 0; i < n; ++i)
{
perm.push_back(4 + 4 * n + i);
perm.push_back(4 + n + i);
for (int j = 0; j < n - 1 - i; ++j)
perm.push_back(face3++);
perm.push_back(4 + 2 * n + i);
perm.push_back(face0++);
perm.push_back(4 + 3 * n + i);
}
perm.push_back(2);
for (int i = 0; i < n; ++i)
{
perm.push_back(4 + 3 * n + i);
perm.push_back(4 + 2 * n + i);
for (int j = 0; j < n - 1 - i; ++j)
perm.push_back(face2++);
perm.push_back(4 + n + i);
perm.push_back(face1++);
perm.push_back(4 + 4 * n + i);
for (int j = 0; j < n - 1 - i; ++j)
{
perm.push_back(face1++);
perm.push_back(face2++);
for (int k = 0; k < n - 2 - i - j; ++k)
perm.push_back(interior++);
perm.push_back(face0++);
perm.push_back(face3++);
}
perm.push_back(4 + i);
perm.push_back(4 + 5 * n + i);
}
perm.push_back(3);
}
Expand Down
6 changes: 3 additions & 3 deletions demo/python/demo_custom_element.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -235,9 +235,9 @@
x = [[], [], [], []]
for _ in range(3):
x[0].append(np.zeros((0, 2)))
x[1].append(np.array([[1 - p[0], p[0]] for p in pts]))
x[1].append(np.array([[0, p[0]] for p in pts]))
x[1].append(np.array([[p[0], 0] for p in pts]))
x[1].append(np.array([[0, p[0]] for p in pts]))
x[1].append(np.array([[1 - p[0], p[0]] for p in pts]))
x[2].append(np.zeros((0, 2)))

# The interpolation matrices for the edges in this example will be have shape (1, 2, len(pts), 1),
Expand All @@ -248,7 +248,7 @@
M = [[], [], [], []]
for _ in range(3):
M[0].append(np.zeros((0, 2, 0, 1)))
for normal in [[-1, -1], [-1, 0], [0, 1]]:
for normal in [[0, 1], [-1, 0], [-1, -1]]:
mat = np.empty((1, 2, len(wts), 1))
mat[0, 0, :, 0] = normal[0] * wts
mat[0, 1, :, 0] = normal[1] * wts
Expand Down
Binary file modified joss/img/tetrahedron_numbering.png
View file
Open in desktop
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
64 changes: 32 additions & 32 deletions joss/img/tetrahedron_numbering.svg
View file
Open in desktop
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file modified joss/img/triangle_numbering.png
View file
Open in desktop
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
16 changes: 8 additions & 8 deletions joss/img/triangle_numbering.svg
View file
Open in desktop
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
12 changes: 6 additions & 6 deletions test/test_bernstein.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -106,21 +106,21 @@ def get_bernstein_polynomials_entity_order(celltype, degree):
(1 - x - y) ** 2,
x**2,
y**2,
2 * x * y,
2 * y * (1 - x - y),
2 * x * (1 - x - y),
2 * y * (1 - x - y),
2 * x * y,
]
if degree == 3:
return [
(1 - x - y) ** 3,
x**3,
y**3,
3 * x**2 * y,
3 * x * y**2,
3 * y * (1 - x - y) ** 2,
3 * y**2 * (1 - x - y),
3 * x * (1 - x - y) ** 2,
3 * x**2 * (1 - x - y),
3 * y * (1 - x - y) ** 2,
3 * y**2 * (1 - x - y),
3 * x**2 * y,
3 * x * y**2,
6 * x * y * (1 - x - y),
]
raise NotImplementedError()
Expand Down
Loading
Loading