From 9469765b13602c61ac84b2e570e62d6c79ef99b9 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Tue, 10 Sep 2024 14:23:36 +0100 Subject: [PATCH 01/35] Make mesh cells all be cell complexes --- firedrake/mesh.py | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 906547ae72..a9fca2cecf 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -1216,7 +1216,6 @@ def _ufl_cell(self): # TODO: this needs to be updated for mixed-cell meshes. nfacets = self._comm.allreduce(nfacets, op=MPI.MAX) - # Note that the geometric dimension of the cell is not set here # despite it being a property of a UFL cell. It will default to # equal the topological dimension. @@ -1224,7 +1223,7 @@ def _ufl_cell(self): # represent a mesh topology (as here) have geometric dimension # equal their topological dimension. This is reflected in the # corresponding UFL mesh. - return ufl.Cell(_cells[tdim][nfacets]) + return constructCellComplex(_cells[tdim][nfacets]) @utils.cached_property def _ufl_mesh(self): @@ -1955,7 +1954,7 @@ def _mark_entity_classes(self): @utils.cached_property def _ufl_cell(self): - return ufl.Cell(_cells[0][0]) + return constructCellComplex(_cells[0][0]) @utils.cached_property def _ufl_mesh(self): From 7aed12d9ec9c33736d756787b4b3bcc5b144b2ea Mon Sep 17 00:00:00 2001 From: India Marsden Date: Mon, 25 Nov 2024 12:25:54 +0000 Subject: [PATCH 02/35] DROP BEFORE MERGE --- .github/workflows/build.yml | 3 +++ 1 file changed, 3 insertions(+) diff --git a/.github/workflows/build.yml b/.github/workflows/build.yml index 909db6990a..3ceca9c802 100644 --- a/.github/workflows/build.yml +++ b/.github/workflows/build.yml @@ -83,6 +83,9 @@ jobs: --install defcon \ --install gadopt \ --install asQ \ + --package-branch tsfc indiamai/new_def_integration \ + --package-branch fiat indiamai/union_different_coeff_size \ + --package-branch FInAT indiamai/new_def_integration \ || (cat firedrake-install.log && /bin/false) - name: Install test dependencies run: | From a457a5679f8ee118da2d0476e389d131d17f8444 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Tue, 26 Nov 2024 11:01:52 +0000 Subject: [PATCH 03/35] First try at adding redefining_fe to packages --- .github/workflows/build.yml | 1 + requirements-git.txt | 1 + 2 files changed, 2 insertions(+) diff --git a/.github/workflows/build.yml b/.github/workflows/build.yml index 3ceca9c802..478d4a95bd 100644 --- a/.github/workflows/build.yml +++ b/.github/workflows/build.yml @@ -86,6 +86,7 @@ jobs: --package-branch tsfc indiamai/new_def_integration \ --package-branch fiat indiamai/union_different_coeff_size \ --package-branch FInAT indiamai/new_def_integration \ + --package-branch redefining_fe indiamai/dual_sets \ || (cat firedrake-install.log && /bin/false) - name: Install test dependencies run: | diff --git a/requirements-git.txt b/requirements-git.txt index c037220ed1..a6369a694d 100644 --- a/requirements-git.txt +++ b/requirements-git.txt @@ -6,3 +6,4 @@ git+https://github.com/dolfin-adjoint/pyadjoint.git#egg=pyadjoint-ad git+https://github.com/firedrakeproject/loopy.git@main#egg=loopy git+https://github.com/firedrakeproject/pytest-mpi.git@main#egg=pytest-mpi git+https://github.com/firedrakeproject/petsc.git@firedrake#egg=petsc +git+https://github.com/indiamai/redefining_fe.git#egg=redefining_fe \ No newline at end of file From 05086228e65e0015a576f6cd1b8c16599d158e55 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Mon, 2 Dec 2024 13:32:03 +0000 Subject: [PATCH 04/35] Fix import accidently deleted in merge --- firedrake/mesh.py | 1 + 1 file changed, 1 insertion(+) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index a9fca2cecf..06492165ac 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -37,6 +37,7 @@ from firedrake.adjoint_utils import MeshGeometryMixin from pyadjoint import stop_annotating import gem +from redefining_fe import constructCellComplex try: import netgen From 94a96dcb192f4159cbce45279ef176c58034885b Mon Sep 17 00:00:00 2001 From: India Marsden Date: Mon, 9 Dec 2024 10:55:45 +0000 Subject: [PATCH 05/35] Add changes from tsfc --- tsfc/finatinterface.py | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/tsfc/finatinterface.py b/tsfc/finatinterface.py index b7e3d0ad72..bcc1ff159e 100644 --- a/tsfc/finatinterface.py +++ b/tsfc/finatinterface.py @@ -25,6 +25,7 @@ import FIAT import finat import finat.ufl +import redefining_fe import ufl __all__ = ("as_fiat_cell", "create_base_element", @@ -111,6 +112,8 @@ def as_fiat_cell(cell): :arg cell: the :class:`ufl.Cell` to convert.""" if not isinstance(cell, ufl.AbstractCell): raise ValueError("Expecting a UFL Cell") + if isinstance(cell, redefining_fe.cells.CellComplexToUFL): + return cell.to_fiat() return FIAT.ufc_cell(cell) @@ -220,6 +223,9 @@ def convert_finiteelement(element, **kwargs): return lmbda(cell, element.degree(), **finat_kwargs), set() +@convert.register(redefining_fe.triples.IndiaTripleUFL) +def convert_india_def(element, **kwargs): + return finat.fiat_elements.IndiaDefElement(element.triple), set() # Element modifiers and compound element types @convert.register(finat.ufl.BrokenElement) From f0813214062c2695620c58b684825c3f72679c78 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Mon, 9 Dec 2024 10:57:53 +0000 Subject: [PATCH 06/35] DROP BEFORE MERGE remove finat and tsfc branches --- .github/workflows/build.yml | 2 -- 1 file changed, 2 deletions(-) diff --git a/.github/workflows/build.yml b/.github/workflows/build.yml index 478d4a95bd..5b1afe7455 100644 --- a/.github/workflows/build.yml +++ b/.github/workflows/build.yml @@ -83,9 +83,7 @@ jobs: --install defcon \ --install gadopt \ --install asQ \ - --package-branch tsfc indiamai/new_def_integration \ --package-branch fiat indiamai/union_different_coeff_size \ - --package-branch FInAT indiamai/new_def_integration \ --package-branch redefining_fe indiamai/dual_sets \ || (cat firedrake-install.log && /bin/false) - name: Install test dependencies From e755dbffc9aa921470659b578c4b146746bd3bb6 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Mon, 9 Dec 2024 10:58:31 +0000 Subject: [PATCH 07/35] Lint --- tsfc/finatinterface.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/tsfc/finatinterface.py b/tsfc/finatinterface.py index bcc1ff159e..9ff7e5d9d6 100644 --- a/tsfc/finatinterface.py +++ b/tsfc/finatinterface.py @@ -223,10 +223,12 @@ def convert_finiteelement(element, **kwargs): return lmbda(cell, element.degree(), **finat_kwargs), set() + @convert.register(redefining_fe.triples.IndiaTripleUFL) def convert_india_def(element, **kwargs): return finat.fiat_elements.IndiaDefElement(element.triple), set() + # Element modifiers and compound element types @convert.register(finat.ufl.BrokenElement) def convert_brokenelement(element, **kwargs): From 59bca7388603d209b49587257341a1baa15d9986 Mon Sep 17 00:00:00 2001 From: India Marsden <37078108+indiamai@users.noreply.github.com> Date: Sat, 11 Jan 2025 12:09:08 +0000 Subject: [PATCH 08/35] Rename redefining fe to fuse (#3960) Rename package --- .github/workflows/build.yml | 4 ++-- firedrake/mesh.py | 2 +- pyproject.toml | 1 + requirements-git.txt | 2 +- 4 files changed, 5 insertions(+), 4 deletions(-) diff --git a/.github/workflows/build.yml b/.github/workflows/build.yml index f4ce3da9e0..fac37eef50 100644 --- a/.github/workflows/build.yml +++ b/.github/workflows/build.yml @@ -84,8 +84,8 @@ jobs: --install defcon \ --install gadopt \ --install asQ \ - --package-branch fiat indiamai/union_different_coeff_size \ - --package-branch redefining_fe indiamai/dual_sets \ + --package-branch fiat indiamai/integrate_fuse \ + --package-branch fuse main \ || (cat firedrake-install.log && /bin/false) - name: Install test dependencies run: | diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 547b7b2e87..65009f719e 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -37,7 +37,7 @@ from firedrake.adjoint_utils import MeshGeometryMixin from pyadjoint import stop_annotating import gem -from redefining_fe import constructCellComplex +from fuse import constructCellComplex try: import netgen diff --git a/pyproject.toml b/pyproject.toml index bb934e4182..892483e37e 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -39,6 +39,7 @@ dependencies = [ "pyadjoint-ad @ git+https://github.com/dolfin-adjoint/pyadjoint.git", "loopy @ git+https://github.com/firedrakeproject/loopy.git@main", "libsupermesh @ git+https://github.com/firedrakeproject/libsupermesh.git", + "fuse-element @ git+https://github.com/indiamai/fuse.git", ] classifiers = [ "Development Status :: 5 - Production/Stable", diff --git a/requirements-git.txt b/requirements-git.txt index ed6fb71bba..990ec151fb 100644 --- a/requirements-git.txt +++ b/requirements-git.txt @@ -5,4 +5,4 @@ git+https://github.com/firedrakeproject/loopy.git@main#egg=loopy git+https://github.com/firedrakeproject/pytest-mpi.git@main#egg=pytest-mpi git+https://github.com/firedrakeproject/petsc.git@firedrake#egg=petsc git+https://github.com/firedrakeproject/libsupermesh.git#egg=libsupermesh -git+https://github.com/indiamai/redefining_fe.git#egg=redefining_fe \ No newline at end of file +git+https://github.com/indiamai/fuse.git#egg=fuse \ No newline at end of file From 6aa685c2475b4b7a541c44e7e206b105443e0230 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Mon, 13 Jan 2025 16:10:01 +0000 Subject: [PATCH 09/35] Add original code back in as comment --- firedrake/mesh.py | 1 + 1 file changed, 1 insertion(+) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 65009f719e..15464ce2a2 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -1224,6 +1224,7 @@ def _ufl_cell(self): # represent a mesh topology (as here) have geometric dimension # equal their topological dimension. This is reflected in the # corresponding UFL mesh. + # return ufl.Cell(_cells[tdim][nfacets]) return constructCellComplex(_cells[tdim][nfacets]) @utils.cached_property From 361d81a4507b341e013bba0709a4a8e37233bb82 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Fri, 17 Jan 2025 13:26:47 +0000 Subject: [PATCH 10/35] Clean up mesh.py --- firedrake/mesh.py | 14 +++++++++++--- 1 file changed, 11 insertions(+), 3 deletions(-) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 15464ce2a2..a6b5156b68 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -16,6 +16,7 @@ import rtree from textwrap import dedent from pathlib import Path +from firedrake.petsc import PETSc from pyop2 import op2 from pyop2.mpi import ( @@ -1263,6 +1264,7 @@ def cell_closure(self): cell = self.ufl_cell() assert tdim == cell.topological_dimension() + if self.submesh_parent is not None: return dmcommon.submesh_create_cell_closure_cell_submesh(plex, self.submesh_parent.topology_dm, @@ -1274,10 +1276,16 @@ def cell_closure(self): entity_per_cell = np.zeros(len(topology), dtype=IntType) for d, ents in topology.items(): entity_per_cell[d] = len(ents) - return dmcommon.closure_ordering(plex, vertex_numbering, cell_numbering, entity_per_cell) - + # elif hasattr(cell, "to_fiat"): + # TODO + # topology = cell.to_fiat().topology + # entity_per_cell = np.zeros(len(topology), dtype=IntType) + # for d, ents in topology.items(): + # entity_per_cell[d] = len(ents) + # return dmcommon.closure_ordering(plex, vertex_numbering, + # cell_numbering, entity_per_cell) elif cell.cellname() == "quadrilateral": from firedrake_citations import Citations Citations().register("Homolya2016") @@ -1308,6 +1316,7 @@ def cell_closure(self): @utils.cached_property def entity_orientations(self): + plex = self.topology_dm return dmcommon.entity_orientations(self, self.cell_closure) @PETSc.Log.EventDecorator() @@ -1333,7 +1342,6 @@ def merge_ids(x, y, datatype): op.Free() else: unique_markers = None - local_facet_number, facet_cell = \ dmcommon.facet_numbering(dm, kind, facets, self._cell_numbering, From ff9764e4d5939f67cd31a9ecc72fe47dc522d459 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Fri, 17 Jan 2025 13:32:38 +0000 Subject: [PATCH 11/35] lint --- firedrake/mesh.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index a6b5156b68..8581a198ec 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -16,7 +16,6 @@ import rtree from textwrap import dedent from pathlib import Path -from firedrake.petsc import PETSc from pyop2 import op2 from pyop2.mpi import ( @@ -1264,7 +1263,7 @@ def cell_closure(self): cell = self.ufl_cell() assert tdim == cell.topological_dimension() - + if self.submesh_parent is not None: return dmcommon.submesh_create_cell_closure_cell_submesh(plex, self.submesh_parent.topology_dm, @@ -1316,7 +1315,6 @@ def cell_closure(self): @utils.cached_property def entity_orientations(self): - plex = self.topology_dm return dmcommon.entity_orientations(self, self.cell_closure) @PETSc.Log.EventDecorator() From 6710bd110c43512305238aaf1708b951dd32e44c Mon Sep 17 00:00:00 2001 From: India Marsden Date: Mon, 20 Jan 2025 15:36:42 +0000 Subject: [PATCH 12/35] temp testing file --- tests/firedrake/regression/test_fuse.py | 21 +++++++++++++++++++++ 1 file changed, 21 insertions(+) create mode 100644 tests/firedrake/regression/test_fuse.py diff --git a/tests/firedrake/regression/test_fuse.py b/tests/firedrake/regression/test_fuse.py new file mode 100644 index 0000000000..8009bbbff3 --- /dev/null +++ b/tests/firedrake/regression/test_fuse.py @@ -0,0 +1,21 @@ +from test_helmholtz import helmholtz +from test_poisson_strong_bcs import run_test +import pytest +import numpy as np + + + +@pytest.mark.parametrize(['params', 'degree', 'quadrilateral'], + [(p, d, q) + for p in [{}, {'snes_type': 'ksponly', 'ksp_type': 'preonly', 'pc_type': 'lu'}] + for d in (1, 2, 3) + for q in [False, True]]) +def test_poisson_analytic(params, degree, quadrilateral): + assert (run_test(2, degree, parameters=params, quadrilateral=False) < 1.e-9) + +def test_helmholtz(): + diff = np.array([helmholtz(i)[0] for i in range(3, 6)]) + print("l2 error norms:", diff) + conv = np.log2(diff[:-1] / diff[1:]) + print("convergence order:", conv) + assert (np.array(conv) > 2.8).all() \ No newline at end of file From 82965b98b1a77af31b66afc5730e2efe0a011d39 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Wed, 22 Jan 2025 10:58:12 +0000 Subject: [PATCH 13/35] Add small subset of tests for easy temporary fuse testing --- tests/firedrake/regression/test_fuse.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/tests/firedrake/regression/test_fuse.py b/tests/firedrake/regression/test_fuse.py index 8009bbbff3..1958fca8d3 100644 --- a/tests/firedrake/regression/test_fuse.py +++ b/tests/firedrake/regression/test_fuse.py @@ -13,9 +13,12 @@ def test_poisson_analytic(params, degree, quadrilateral): assert (run_test(2, degree, parameters=params, quadrilateral=False) < 1.e-9) -def test_helmholtz(): - diff = np.array([helmholtz(i)[0] for i in range(3, 6)]) +@pytest.mark.parametrize(['conv_num', 'degree'], + [(p, d) + for p,d in zip([1.8, 2.8, 3.8],[1, 2, 3])]) +def test_helmholtz(conv_num, degree): + diff = np.array([helmholtz(i, degree=degree)[0] for i in range(3, 6)]) print("l2 error norms:", diff) conv = np.log2(diff[:-1] / diff[1:]) print("convergence order:", conv) - assert (np.array(conv) > 2.8).all() \ No newline at end of file + assert (np.array(conv) > conv_num).all() \ No newline at end of file From 3e13b4fd579d5c817db2382b7ddbdf1078d59952 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Tue, 10 Sep 2024 14:23:36 +0100 Subject: [PATCH 14/35] Make mesh cells all be cell complexes --- firedrake/mesh.py | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 0be58d6d95..be9724a7b8 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -1216,7 +1216,6 @@ def _ufl_cell(self): # TODO: this needs to be updated for mixed-cell meshes. nfacets = self._comm.allreduce(nfacets, op=MPI.MAX) - # Note that the geometric dimension of the cell is not set here # despite it being a property of a UFL cell. It will default to # equal the topological dimension. @@ -1224,7 +1223,7 @@ def _ufl_cell(self): # represent a mesh topology (as here) have geometric dimension # equal their topological dimension. This is reflected in the # corresponding UFL mesh. - return ufl.Cell(_cells[tdim][nfacets]) + return constructCellComplex(_cells[tdim][nfacets]) @utils.cached_property def _ufl_mesh(self): @@ -1955,7 +1954,7 @@ def _mark_entity_classes(self): @utils.cached_property def _ufl_cell(self): - return ufl.Cell(_cells[0][0]) + return constructCellComplex(_cells[0][0]) @utils.cached_property def _ufl_mesh(self): From 79e8f03a4dd8cd69b652db4f1b6f9e3b6ad186af Mon Sep 17 00:00:00 2001 From: India Marsden Date: Mon, 2 Dec 2024 13:32:03 +0000 Subject: [PATCH 15/35] Fix import accidently deleted in merge --- firedrake/mesh.py | 1 + 1 file changed, 1 insertion(+) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index be9724a7b8..3819287595 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -37,6 +37,7 @@ from firedrake.adjoint_utils import MeshGeometryMixin from pyadjoint import stop_annotating import gem +from redefining_fe import constructCellComplex try: import netgen From 779efdd525f5d8db83f0bcf3cef0cc61c4d3e70e Mon Sep 17 00:00:00 2001 From: India Marsden <37078108+indiamai@users.noreply.github.com> Date: Sat, 11 Jan 2025 12:09:08 +0000 Subject: [PATCH 16/35] Rename redefining fe to fuse (#3960) Rename package --- firedrake/mesh.py | 2 +- pyproject.toml | 1 + requirements-git.txt | 1 + 3 files changed, 3 insertions(+), 1 deletion(-) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 3819287595..179cec153b 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -37,7 +37,7 @@ from firedrake.adjoint_utils import MeshGeometryMixin from pyadjoint import stop_annotating import gem -from redefining_fe import constructCellComplex +from fuse import constructCellComplex try: import netgen diff --git a/pyproject.toml b/pyproject.toml index 93aab13cf9..2359c43e60 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -38,6 +38,7 @@ dependencies = [ "fenics-ufl @ git+https://github.com/firedrakeproject/ufl.git", "firedrake-fiat @ git+https://github.com/firedrakeproject/fiat.git", "islpy>=2025.1.5; sys_platform == 'darwin'", + "fuse-element @ git+https://github.com/indiamai/fuse.git", ] classifiers = [ "Development Status :: 5 - Production/Stable", diff --git a/requirements-git.txt b/requirements-git.txt index 4cf888c821..791ce0cdf7 100644 --- a/requirements-git.txt +++ b/requirements-git.txt @@ -4,3 +4,4 @@ git+https://github.com/dolfin-adjoint/pyadjoint.git#egg=pyadjoint-ad git+https://github.com/firedrakeproject/loopy.git@main#egg=loopy git+https://github.com/firedrakeproject/petsc.git@firedrake#egg=petsc git+https://github.com/firedrakeproject/libsupermesh.git#egg=libsupermesh +git+https://github.com/indiamai/fuse.git#egg=fuse \ No newline at end of file From 9be37906b47fc518d7422b88faff2f09a9fe0cdb Mon Sep 17 00:00:00 2001 From: India Marsden Date: Mon, 13 Jan 2025 16:10:01 +0000 Subject: [PATCH 17/35] Add original code back in as comment --- firedrake/mesh.py | 1 + 1 file changed, 1 insertion(+) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 179cec153b..91b0011edc 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -1224,6 +1224,7 @@ def _ufl_cell(self): # represent a mesh topology (as here) have geometric dimension # equal their topological dimension. This is reflected in the # corresponding UFL mesh. + # return ufl.Cell(_cells[tdim][nfacets]) return constructCellComplex(_cells[tdim][nfacets]) @utils.cached_property From 1be38761fbbe7fb8ee392e845765d855dbad207a Mon Sep 17 00:00:00 2001 From: India Marsden Date: Fri, 17 Jan 2025 13:26:47 +0000 Subject: [PATCH 18/35] Clean up mesh.py --- firedrake/mesh.py | 14 +++++++++++--- 1 file changed, 11 insertions(+), 3 deletions(-) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 91b0011edc..205f6a2a7d 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -16,6 +16,7 @@ import rtree from textwrap import dedent from pathlib import Path +from firedrake.petsc import PETSc from pyop2 import op2 from pyop2.mpi import ( @@ -1263,6 +1264,7 @@ def cell_closure(self): cell = self.ufl_cell() assert tdim == cell.topological_dimension() + if self.submesh_parent is not None: return dmcommon.submesh_create_cell_closure_cell_submesh(plex, self.submesh_parent.topology_dm, @@ -1274,10 +1276,16 @@ def cell_closure(self): entity_per_cell = np.zeros(len(topology), dtype=IntType) for d, ents in topology.items(): entity_per_cell[d] = len(ents) - return dmcommon.closure_ordering(plex, vertex_numbering, cell_numbering, entity_per_cell) - + # elif hasattr(cell, "to_fiat"): + # TODO + # topology = cell.to_fiat().topology + # entity_per_cell = np.zeros(len(topology), dtype=IntType) + # for d, ents in topology.items(): + # entity_per_cell[d] = len(ents) + # return dmcommon.closure_ordering(plex, vertex_numbering, + # cell_numbering, entity_per_cell) elif cell.cellname() == "quadrilateral": from firedrake_citations import Citations Citations().register("Homolya2016") @@ -1308,6 +1316,7 @@ def cell_closure(self): @utils.cached_property def entity_orientations(self): + plex = self.topology_dm return dmcommon.entity_orientations(self, self.cell_closure) @PETSc.Log.EventDecorator() @@ -1333,7 +1342,6 @@ def merge_ids(x, y, datatype): op.Free() else: unique_markers = None - local_facet_number, facet_cell = \ dmcommon.facet_numbering(dm, kind, facets, self._cell_numbering, From 60cff449923c95382344a20949ac7259afc6488c Mon Sep 17 00:00:00 2001 From: India Marsden Date: Fri, 17 Jan 2025 13:32:38 +0000 Subject: [PATCH 19/35] lint --- firedrake/mesh.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 205f6a2a7d..72c6c0284c 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -16,7 +16,6 @@ import rtree from textwrap import dedent from pathlib import Path -from firedrake.petsc import PETSc from pyop2 import op2 from pyop2.mpi import ( @@ -1264,7 +1263,7 @@ def cell_closure(self): cell = self.ufl_cell() assert tdim == cell.topological_dimension() - + if self.submesh_parent is not None: return dmcommon.submesh_create_cell_closure_cell_submesh(plex, self.submesh_parent.topology_dm, @@ -1316,7 +1315,6 @@ def cell_closure(self): @utils.cached_property def entity_orientations(self): - plex = self.topology_dm return dmcommon.entity_orientations(self, self.cell_closure) @PETSc.Log.EventDecorator() From 8642f0dc56135986cbec44483688d708be468782 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Mon, 20 Jan 2025 15:36:42 +0000 Subject: [PATCH 20/35] temp testing file --- tests/firedrake/regression/test_fuse.py | 21 +++++++++++++++++++++ 1 file changed, 21 insertions(+) create mode 100644 tests/firedrake/regression/test_fuse.py diff --git a/tests/firedrake/regression/test_fuse.py b/tests/firedrake/regression/test_fuse.py new file mode 100644 index 0000000000..8009bbbff3 --- /dev/null +++ b/tests/firedrake/regression/test_fuse.py @@ -0,0 +1,21 @@ +from test_helmholtz import helmholtz +from test_poisson_strong_bcs import run_test +import pytest +import numpy as np + + + +@pytest.mark.parametrize(['params', 'degree', 'quadrilateral'], + [(p, d, q) + for p in [{}, {'snes_type': 'ksponly', 'ksp_type': 'preonly', 'pc_type': 'lu'}] + for d in (1, 2, 3) + for q in [False, True]]) +def test_poisson_analytic(params, degree, quadrilateral): + assert (run_test(2, degree, parameters=params, quadrilateral=False) < 1.e-9) + +def test_helmholtz(): + diff = np.array([helmholtz(i)[0] for i in range(3, 6)]) + print("l2 error norms:", diff) + conv = np.log2(diff[:-1] / diff[1:]) + print("convergence order:", conv) + assert (np.array(conv) > 2.8).all() \ No newline at end of file From c5890b1bdee79c26bef1fdd53c019bc931d78a05 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Wed, 22 Jan 2025 10:58:12 +0000 Subject: [PATCH 21/35] Add small subset of tests for easy temporary fuse testing --- tests/firedrake/regression/test_fuse.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/tests/firedrake/regression/test_fuse.py b/tests/firedrake/regression/test_fuse.py index 8009bbbff3..1958fca8d3 100644 --- a/tests/firedrake/regression/test_fuse.py +++ b/tests/firedrake/regression/test_fuse.py @@ -13,9 +13,12 @@ def test_poisson_analytic(params, degree, quadrilateral): assert (run_test(2, degree, parameters=params, quadrilateral=False) < 1.e-9) -def test_helmholtz(): - diff = np.array([helmholtz(i)[0] for i in range(3, 6)]) +@pytest.mark.parametrize(['conv_num', 'degree'], + [(p, d) + for p,d in zip([1.8, 2.8, 3.8],[1, 2, 3])]) +def test_helmholtz(conv_num, degree): + diff = np.array([helmholtz(i, degree=degree)[0] for i in range(3, 6)]) print("l2 error norms:", diff) conv = np.log2(diff[:-1] / diff[1:]) print("convergence order:", conv) - assert (np.array(conv) > 2.8).all() \ No newline at end of file + assert (np.array(conv) > conv_num).all() \ No newline at end of file From 3dcf135670ab15196cd69c93f98da370e3407a85 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Tue, 29 Apr 2025 10:34:29 +0100 Subject: [PATCH 22/35] Change repo links --- pyproject.toml | 2 +- requirements-git.txt | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/pyproject.toml b/pyproject.toml index 2359c43e60..226c687618 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -38,7 +38,7 @@ dependencies = [ "fenics-ufl @ git+https://github.com/firedrakeproject/ufl.git", "firedrake-fiat @ git+https://github.com/firedrakeproject/fiat.git", "islpy>=2025.1.5; sys_platform == 'darwin'", - "fuse-element @ git+https://github.com/indiamai/fuse.git", + "fuse-element @ git+https://github.com/firedrakeproject/fuse.git", ] classifiers = [ "Development Status :: 5 - Production/Stable", diff --git a/requirements-git.txt b/requirements-git.txt index 791ce0cdf7..2533340062 100644 --- a/requirements-git.txt +++ b/requirements-git.txt @@ -4,4 +4,4 @@ git+https://github.com/dolfin-adjoint/pyadjoint.git#egg=pyadjoint-ad git+https://github.com/firedrakeproject/loopy.git@main#egg=loopy git+https://github.com/firedrakeproject/petsc.git@firedrake#egg=petsc git+https://github.com/firedrakeproject/libsupermesh.git#egg=libsupermesh -git+https://github.com/indiamai/fuse.git#egg=fuse \ No newline at end of file +git+https://github.com/firedrakeproject/fuse.git#egg=fuse \ No newline at end of file From 00ef926d42cd6fdc7e58f969b49d1e0f3bcbedd7 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Tue, 29 Apr 2025 10:41:22 +0100 Subject: [PATCH 23/35] lint --- tests/firedrake/regression/test_fuse.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/tests/firedrake/regression/test_fuse.py b/tests/firedrake/regression/test_fuse.py index 1958fca8d3..cbae8cdcd4 100644 --- a/tests/firedrake/regression/test_fuse.py +++ b/tests/firedrake/regression/test_fuse.py @@ -4,7 +4,6 @@ import numpy as np - @pytest.mark.parametrize(['params', 'degree', 'quadrilateral'], [(p, d, q) for p in [{}, {'snes_type': 'ksponly', 'ksp_type': 'preonly', 'pc_type': 'lu'}] @@ -13,12 +12,13 @@ def test_poisson_analytic(params, degree, quadrilateral): assert (run_test(2, degree, parameters=params, quadrilateral=False) < 1.e-9) + @pytest.mark.parametrize(['conv_num', 'degree'], [(p, d) - for p,d in zip([1.8, 2.8, 3.8],[1, 2, 3])]) + for p, d in zip([1.8, 2.8, 3.8], [1, 2, 3])]) def test_helmholtz(conv_num, degree): diff = np.array([helmholtz(i, degree=degree)[0] for i in range(3, 6)]) print("l2 error norms:", diff) conv = np.log2(diff[:-1] / diff[1:]) print("convergence order:", conv) - assert (np.array(conv) > conv_num).all() \ No newline at end of file + assert (np.array(conv) > conv_num).all() From d26e3ad5c6cea84444645c35a98f54ee39c71efd Mon Sep 17 00:00:00 2001 From: India Marsden <37078108+indiamai@users.noreply.github.com> Date: Sat, 11 Jan 2025 12:09:08 +0000 Subject: [PATCH 24/35] Rename redefining fe to fuse (#3960) Rename package --- requirements-git.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/requirements-git.txt b/requirements-git.txt index 09703174da..2533340062 100644 --- a/requirements-git.txt +++ b/requirements-git.txt @@ -4,4 +4,4 @@ git+https://github.com/dolfin-adjoint/pyadjoint.git#egg=pyadjoint-ad git+https://github.com/firedrakeproject/loopy.git@main#egg=loopy git+https://github.com/firedrakeproject/petsc.git@firedrake#egg=petsc git+https://github.com/firedrakeproject/libsupermesh.git#egg=libsupermesh -git+https://github.com/firedrakeproject/fuse.git#egg=fuse +git+https://github.com/firedrakeproject/fuse.git#egg=fuse \ No newline at end of file From cdc0b11f7fd4ed16089c446df2ac5299a0d8acb1 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Tue, 4 Feb 2025 16:10:53 +0000 Subject: [PATCH 25/35] Brute forcing a fuse quad in --- firedrake/cython/dmcommon.pyx | 42 +++++++++++++++++++++++++++++---- firedrake/mesh.py | 44 +++++++++++++++++++++++------------ 2 files changed, 67 insertions(+), 19 deletions(-) diff --git a/firedrake/cython/dmcommon.pyx b/firedrake/cython/dmcommon.pyx index b67c0c21ec..7ea1cc3f8d 100644 --- a/firedrake/cython/dmcommon.pyx +++ b/firedrake/cython/dmcommon.pyx @@ -363,10 +363,20 @@ cdef inline PetscInt _reorder_plex_cone(PETSc.DM dm, # 0 2 # | | # +---1---+ + # FUSE +---1---+ + # | | + # 0 2 + # | | + # +---3---+ + # plex_cone_new[0] = plex_cone_old[0] + # plex_cone_new[1] = plex_cone_old[2] + # plex_cone_new[2] = plex_cone_old[1] + # plex_cone_new[3] = plex_cone_old[3] + print("SWITCHING") plex_cone_new[0] = plex_cone_old[0] - plex_cone_new[1] = plex_cone_old[2] - plex_cone_new[2] = plex_cone_old[1] - plex_cone_new[3] = plex_cone_old[3] + plex_cone_new[1] = plex_cone_old[3] + plex_cone_new[2] = plex_cone_old[2] + plex_cone_new[3] = plex_cone_old[1] elif dm.getCellType(p) == PETSc.DM.PolytopeType.HEXAHEDRON: # UFCHexahedron: +-------+ +-------+ # /. | / 5 /| @@ -466,7 +476,22 @@ cdef inline PetscInt _reorder_plex_closure(PETSc.DM dm, # 1 0 3 # | | # 6---2---7 - raise NotImplementedError(f"Not implemented for {dm.getCellType(p)}") + # + # FUSE 1---5---3 + # | | + # 4 8 6 + # | | + # 0---7---2 + fiat_closure[0] = plex_closure[2 * 6] + fiat_closure[1] = plex_closure[2 * 5] + fiat_closure[2] = plex_closure[2 * 7] + fiat_closure[3] = plex_closure[2 * 8] + fiat_closure[4] = plex_closure[2 * 1] + fiat_closure[5] = plex_closure[2 * 4] + fiat_closure[6] = plex_closure[2 * 3] + fiat_closure[7] = plex_closure[2 * 2] + fiat_closure[8] = plex_closure[2 * 0] + # raise NotImplementedError(f"Not implemented for {dm.getCellType(p)}") elif dm.getCellType(p) == PETSc.DM.PolytopeType.HEXAHEDRON: # UFCHexahedron: 3--19---7 3--19---7 # 13. | 13 25 15| @@ -953,11 +978,15 @@ cdef inline PetscInt _compute_orientation_simplex(PetscInt *fiat_cone, CHKERR(PetscMalloc1(coneSize, &cone1)) CHKERR(PetscMalloc1(coneSize, &inds)) + print("plex") for k in range(coneSize1): cone1[k] = plex_cone[k] + print(plex_cone[k]) n = 0 + print("fiat") for e in range(coneSize): q = fiat_cone[e] + print(q) for k in range(coneSize1): if q == cone1[k]: inds[n] = k @@ -1011,8 +1040,13 @@ cdef inline PetscInt _compute_orientation_interval_tensor_product(PetscInt *fiat # 0---1 1---0 2---3 3---2 # dim1 = dim + print("plex") for i in range(2 * dim): + print(plex_cone[i]) plex_cone_copy[i] = plex_cone[i] + print("fiat") + for i in range(2 * dim): + print(fiat_cone[i]) for i in range(dim): for j in range(dim1): if plex_cone_copy[2 * j] == fiat_cone[2 * i] or plex_cone_copy[2 * j + 1] == fiat_cone[2 * i]: diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 72c6c0284c..985ca32668 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -830,8 +830,7 @@ def make_cell_node_list(self, global_numbering, entity_dofs, entity_permutations :arg entity_permutations: FInAT element entity permutations :arg offsets: layer offsets for each entity dof (may be None). """ - return dmcommon.get_cell_nodes(self, global_numbering, - entity_dofs, entity_permutations, offsets) + return dmcommon.get_cell_nodes(self, global_numbering, entity_dofs, entity_permutations, offsets) def make_dofs_per_plex_entity(self, entity_dofs): """Returns the number of DoFs per plex entity for each stratum, @@ -1224,6 +1223,12 @@ def _ufl_cell(self): # represent a mesh topology (as here) have geometric dimension # equal their topological dimension. This is reflected in the # corresponding UFL mesh. + cell1 = ufl.Cell(_cells[tdim][nfacets]) + cell2 = constructCellComplex(_cells[tdim][nfacets]) + print(as_fiat_cell(cell1).topology) + print(as_fiat_cell(cell1).vertices) + print(as_fiat_cell(cell2).topology) + print(as_fiat_cell(cell2).vertices) # return ufl.Cell(_cells[tdim][nfacets]) return constructCellComplex(_cells[tdim][nfacets]) @@ -1263,7 +1268,9 @@ def cell_closure(self): cell = self.ufl_cell() assert tdim == cell.topological_dimension() - + # plex.viewFromOptions("-dm_view") + closure, _ = plex.getTransitiveClosure(0) + print("plex", closure) if self.submesh_parent is not None: return dmcommon.submesh_create_cell_closure_cell_submesh(plex, self.submesh_parent.topology_dm, @@ -1271,40 +1278,47 @@ def cell_closure(self): self.submesh_parent._cell_numbering, self.submesh_parent.cell_closure) elif cell.is_simplex(): + print("simplex case") topology = FIAT.ufc_cell(cell).get_topology() entity_per_cell = np.zeros(len(topology), dtype=IntType) for d, ents in topology.items(): entity_per_cell[d] = len(ents) return dmcommon.closure_ordering(plex, vertex_numbering, cell_numbering, entity_per_cell) - # elif hasattr(cell, "to_fiat"): - # TODO - # topology = cell.to_fiat().topology - # entity_per_cell = np.zeros(len(topology), dtype=IntType) - # for d, ents in topology.items(): - # entity_per_cell[d] = len(ents) - # return dmcommon.closure_ordering(plex, vertex_numbering, - # cell_numbering, entity_per_cell) + elif hasattr(cell, "to_fiat"): + # TODO + topology = cell.to_fiat().topology + # entity_per_cell = np.zeros(len(topology), dtype=IntType) + # for d, ents in topology.items(): + # entity_per_cell[d] = len(ents) + # return dmcommon.closure_ordering(plex, vertex_numbering, + # cell_numbering, entity_per_cell) + topology = FIAT.ufc_cell(cell).get_topology() + closureSize = sum([len(ents) for _, ents in topology.items()]) + return dmcommon.create_cell_closure(plex, cell_numbering, closureSize) elif cell.cellname() == "quadrilateral": + print("quadcase") from firedrake_citations import Citations Citations().register("Homolya2016") Citations().register("McRae2016") # Quadrilateral mesh cell_ranks = dmcommon.get_cell_remote_ranks(plex) - + print(plex.getTransitiveClosure(0)) facet_orientations = dmcommon.quadrilateral_facet_orientations( plex, vertex_numbering, cell_ranks) - + print(facet_orientations) cell_orientations = dmcommon.orientations_facet2cell( plex, vertex_numbering, cell_ranks, facet_orientations, cell_numbering) - + print(cell_orientations) dmcommon.exchange_cell_orientations(plex, cell_numbering, cell_orientations) - return dmcommon.quadrilateral_closure_ordering( + res = dmcommon.quadrilateral_closure_ordering( plex, vertex_numbering, cell_numbering, cell_orientations) + print(res) + return res elif cell.cellname() == "hexahedron": # TODO: Should change and use create_cell_closure() for all cell types. topology = FIAT.ufc_cell(cell).get_topology() From f97141359a101b90d7776a2a99ef771d2ec61d07 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Thu, 20 Feb 2025 17:46:17 +0000 Subject: [PATCH 26/35] move more things to generic as cell, further work on qudas --- firedrake/cython/dmcommon.pyx | 21 ++++++++------------- firedrake/functionspace.py | 3 ++- firedrake/mesh.py | 30 ++++++++++++------------------ firedrake/tsfc_interface.py | 2 +- 4 files changed, 23 insertions(+), 33 deletions(-) diff --git a/firedrake/cython/dmcommon.pyx b/firedrake/cython/dmcommon.pyx index 7ea1cc3f8d..ef6dfdde48 100644 --- a/firedrake/cython/dmcommon.pyx +++ b/firedrake/cython/dmcommon.pyx @@ -368,10 +368,10 @@ cdef inline PetscInt _reorder_plex_cone(PETSc.DM dm, # 0 2 # | | # +---3---+ - # plex_cone_new[0] = plex_cone_old[0] - # plex_cone_new[1] = plex_cone_old[2] - # plex_cone_new[2] = plex_cone_old[1] - # plex_cone_new[3] = plex_cone_old[3] + #plex_cone_new[0] = plex_cone_old[0] + #plex_cone_new[1] = plex_cone_old[2] + #plex_cone_new[2] = plex_cone_old[1] + #plex_cone_new[3] = plex_cone_old[3] print("SWITCHING") plex_cone_new[0] = plex_cone_old[0] plex_cone_new[1] = plex_cone_old[3] @@ -978,15 +978,15 @@ cdef inline PetscInt _compute_orientation_simplex(PetscInt *fiat_cone, CHKERR(PetscMalloc1(coneSize, &cone1)) CHKERR(PetscMalloc1(coneSize, &inds)) - print("plex") + #print("plex") for k in range(coneSize1): cone1[k] = plex_cone[k] - print(plex_cone[k]) + #print(plex_cone[k]) n = 0 - print("fiat") + #print("fiat") for e in range(coneSize): q = fiat_cone[e] - print(q) + #print(q) for k in range(coneSize1): if q == cone1[k]: inds[n] = k @@ -1040,13 +1040,8 @@ cdef inline PetscInt _compute_orientation_interval_tensor_product(PetscInt *fiat # 0---1 1---0 2---3 3---2 # dim1 = dim - print("plex") for i in range(2 * dim): - print(plex_cone[i]) plex_cone_copy[i] = plex_cone[i] - print("fiat") - for i in range(2 * dim): - print(fiat_cone[i]) for i in range(dim): for j in range(dim1): if plex_cone_copy[2 * j] == fiat_cone[2 * i] or plex_cone_copy[2 * j + 1] == fiat_cone[2 * i]: diff --git a/firedrake/functionspace.py b/firedrake/functionspace.py index 509f7a2863..826134be0d 100644 --- a/firedrake/functionspace.py +++ b/firedrake/functionspace.py @@ -5,6 +5,7 @@ backwards-compatibility, argument checking, and dispatch. """ import ufl +from ufl.cell import as_cell import finat.ufl from pyop2.utils import flatten @@ -62,7 +63,7 @@ def make_scalar_element(mesh, family, degree, vfamily, vdegree, variant): degree=degree, variant=variant) # If second element was passed in, use it lb = finat.ufl.FiniteElement(vfamily, - cell=ufl.interval, + cell=as_cell("interval"), degree=vdegree, variant=variant) # Now make the TensorProductElement return finat.ufl.TensorProductElement(la, lb) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 985ca32668..f15334b9ec 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -47,6 +47,7 @@ # Only for docstring import mpi4py # noqa: F401 from finat.element_factory import as_fiat_cell +from ufl.cell import as_cell __all__ = [ @@ -1223,14 +1224,7 @@ def _ufl_cell(self): # represent a mesh topology (as here) have geometric dimension # equal their topological dimension. This is reflected in the # corresponding UFL mesh. - cell1 = ufl.Cell(_cells[tdim][nfacets]) - cell2 = constructCellComplex(_cells[tdim][nfacets]) - print(as_fiat_cell(cell1).topology) - print(as_fiat_cell(cell1).vertices) - print(as_fiat_cell(cell2).topology) - print(as_fiat_cell(cell2).vertices) - # return ufl.Cell(_cells[tdim][nfacets]) - return constructCellComplex(_cells[tdim][nfacets]) + return as_cell(_cells[tdim][nfacets]) @utils.cached_property def _ufl_mesh(self): @@ -1270,7 +1264,7 @@ def cell_closure(self): assert tdim == cell.topological_dimension() # plex.viewFromOptions("-dm_view") closure, _ = plex.getTransitiveClosure(0) - print("plex", closure) + # print("plex", closure) if self.submesh_parent is not None: return dmcommon.submesh_create_cell_closure_cell_submesh(plex, self.submesh_parent.topology_dm, @@ -1278,7 +1272,7 @@ def cell_closure(self): self.submesh_parent._cell_numbering, self.submesh_parent.cell_closure) elif cell.is_simplex(): - print("simplex case") + # print("simplex case") topology = FIAT.ufc_cell(cell).get_topology() entity_per_cell = np.zeros(len(topology), dtype=IntType) for d, ents in topology.items(): @@ -1297,27 +1291,27 @@ def cell_closure(self): closureSize = sum([len(ents) for _, ents in topology.items()]) return dmcommon.create_cell_closure(plex, cell_numbering, closureSize) elif cell.cellname() == "quadrilateral": - print("quadcase") + # print("quadcase") from firedrake_citations import Citations Citations().register("Homolya2016") Citations().register("McRae2016") # Quadrilateral mesh cell_ranks = dmcommon.get_cell_remote_ranks(plex) - print(plex.getTransitiveClosure(0)) + # print(plex.getTransitiveClosure(0)) facet_orientations = dmcommon.quadrilateral_facet_orientations( plex, vertex_numbering, cell_ranks) - print(facet_orientations) + # print(facet_orientations) cell_orientations = dmcommon.orientations_facet2cell( plex, vertex_numbering, cell_ranks, facet_orientations, cell_numbering) - print(cell_orientations) + # print(cell_orientations) dmcommon.exchange_cell_orientations(plex, cell_numbering, cell_orientations) res = dmcommon.quadrilateral_closure_ordering( plex, vertex_numbering, cell_numbering, cell_orientations) - print(res) + # print(res) return res elif cell.cellname() == "hexahedron": # TODO: Should change and use create_cell_closure() for all cell types. @@ -1758,7 +1752,7 @@ def __init__(self, mesh, layers, periodic=False, name=None): @utils.cached_property def _ufl_cell(self): - return ufl.TensorProductCell(self._base_mesh.ufl_cell(), ufl.interval) + return ufl.TensorProductCell(self._base_mesh.ufl_cell(), as_cell("interval")) @utils.cached_property def _ufl_mesh(self): @@ -3256,9 +3250,9 @@ def ExtrudedMesh(mesh, layers, layer_height=None, extrusion_type='uniform', peri if extrusion_type == 'radial_hedgehog': helement = helement.reconstruct(family="DG", variant="equispaced") if periodic: - velement = finat.ufl.FiniteElement("DP", ufl.interval, 1, variant="equispaced") + velement = finat.ufl.FiniteElement("DP", as_cell("interval"), 1, variant="equispaced") else: - velement = finat.ufl.FiniteElement("Lagrange", ufl.interval, 1) + velement = finat.ufl.FiniteElement("Lagrange", as_cell("interval"), 1) element = finat.ufl.TensorProductElement(helement, velement) if gdim is None: diff --git a/firedrake/tsfc_interface.py b/firedrake/tsfc_interface.py index d6c663f5de..2c2a7ddfba 100644 --- a/firedrake/tsfc_interface.py +++ b/firedrake/tsfc_interface.py @@ -190,7 +190,6 @@ def compile_form(form, name, parameters=None, split=True, interface=None, diagon _ = parameters parameters = default_parameters["form_compiler"].copy() parameters.update(_) - kernels = [] numbering = form.terminal_numbering() if split: @@ -217,6 +216,7 @@ def compile_form(form, name, parameters=None, split=True, interface=None, diagon numbering[c] for c in extract_firedrake_constants(f) ) prefix = name + "".join(map(str, (i for i in idx if i is not None))) + breakpoint() tsfc_kernel = TSFCKernel( f, prefix, From ac9a9710203db6f1620c4d4c89c3b281635830f8 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Wed, 26 Feb 2025 11:25:17 +0000 Subject: [PATCH 27/35] clean up --- firedrake/cython/dmcommon.pyx | 16 ++++++++-------- firedrake/mesh.py | 29 +++++++++++------------------ firedrake/tsfc_interface.py | 1 - tsfc/driver.py | 1 - 4 files changed, 19 insertions(+), 28 deletions(-) diff --git a/firedrake/cython/dmcommon.pyx b/firedrake/cython/dmcommon.pyx index ef6dfdde48..23bf9157c9 100644 --- a/firedrake/cython/dmcommon.pyx +++ b/firedrake/cython/dmcommon.pyx @@ -368,15 +368,15 @@ cdef inline PetscInt _reorder_plex_cone(PETSc.DM dm, # 0 2 # | | # +---3---+ - #plex_cone_new[0] = plex_cone_old[0] - #plex_cone_new[1] = plex_cone_old[2] - #plex_cone_new[2] = plex_cone_old[1] - #plex_cone_new[3] = plex_cone_old[3] - print("SWITCHING") plex_cone_new[0] = plex_cone_old[0] - plex_cone_new[1] = plex_cone_old[3] - plex_cone_new[2] = plex_cone_old[2] - plex_cone_new[3] = plex_cone_old[1] + plex_cone_new[1] = plex_cone_old[2] + plex_cone_new[2] = plex_cone_old[1] + plex_cone_new[3] = plex_cone_old[3] + # print("SWITCHING") + #plex_cone_new[0] = plex_cone_old[0] + #plex_cone_new[1] = plex_cone_old[3] + #plex_cone_new[2] = plex_cone_old[2] + #plex_cone_new[3] = plex_cone_old[1] elif dm.getCellType(p) == PETSc.DM.PolytopeType.HEXAHEDRON: # UFCHexahedron: +-------+ +-------+ # /. | / 5 /| diff --git a/firedrake/mesh.py b/firedrake/mesh.py index f15334b9ec..2dcc3b51eb 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -1264,7 +1264,6 @@ def cell_closure(self): assert tdim == cell.topological_dimension() # plex.viewFromOptions("-dm_view") closure, _ = plex.getTransitiveClosure(0) - # print("plex", closure) if self.submesh_parent is not None: return dmcommon.submesh_create_cell_closure_cell_submesh(plex, self.submesh_parent.topology_dm, @@ -1272,46 +1271,40 @@ def cell_closure(self): self.submesh_parent._cell_numbering, self.submesh_parent.cell_closure) elif cell.is_simplex(): - # print("simplex case") topology = FIAT.ufc_cell(cell).get_topology() entity_per_cell = np.zeros(len(topology), dtype=IntType) for d, ents in topology.items(): entity_per_cell[d] = len(ents) return dmcommon.closure_ordering(plex, vertex_numbering, cell_numbering, entity_per_cell) - elif hasattr(cell, "to_fiat"): - # TODO - topology = cell.to_fiat().topology - # entity_per_cell = np.zeros(len(topology), dtype=IntType) - # for d, ents in topology.items(): - # entity_per_cell[d] = len(ents) - # return dmcommon.closure_ordering(plex, vertex_numbering, - # cell_numbering, entity_per_cell) - topology = FIAT.ufc_cell(cell).get_topology() - closureSize = sum([len(ents) for _, ents in topology.items()]) - return dmcommon.create_cell_closure(plex, cell_numbering, closureSize) + # elif hasattr(cell, "to_fiat"): + # # TODO + # topology = cell.to_fiat().topology + # # entity_per_cell = np.zeros(len(topology), dtype=IntType) + # # for d, ents in topology.items(): + # # entity_per_cell[d] = len(ents) + # # return dmcommon.closure_ordering(plex, vertex_numbering, + # # cell_numbering, entity_per_cell) + # topology = FIAT.ufc_cell(cell).get_topology() + # closureSize = sum([len(ents) for _, ents in topology.items()]) + # return dmcommon.create_cell_closure(plex, cell_numbering, closureSize) elif cell.cellname() == "quadrilateral": - # print("quadcase") from firedrake_citations import Citations Citations().register("Homolya2016") Citations().register("McRae2016") # Quadrilateral mesh cell_ranks = dmcommon.get_cell_remote_ranks(plex) - # print(plex.getTransitiveClosure(0)) facet_orientations = dmcommon.quadrilateral_facet_orientations( plex, vertex_numbering, cell_ranks) - # print(facet_orientations) cell_orientations = dmcommon.orientations_facet2cell( plex, vertex_numbering, cell_ranks, facet_orientations, cell_numbering) - # print(cell_orientations) dmcommon.exchange_cell_orientations(plex, cell_numbering, cell_orientations) res = dmcommon.quadrilateral_closure_ordering( plex, vertex_numbering, cell_numbering, cell_orientations) - # print(res) return res elif cell.cellname() == "hexahedron": # TODO: Should change and use create_cell_closure() for all cell types. diff --git a/firedrake/tsfc_interface.py b/firedrake/tsfc_interface.py index 2c2a7ddfba..2758e53e2c 100644 --- a/firedrake/tsfc_interface.py +++ b/firedrake/tsfc_interface.py @@ -216,7 +216,6 @@ def compile_form(form, name, parameters=None, split=True, interface=None, diagon numbering[c] for c in extract_firedrake_constants(f) ) prefix = name + "".join(map(str, (i for i in idx if i is not None))) - breakpoint() tsfc_kernel = TSFCKernel( f, prefix, diff --git a/tsfc/driver.py b/tsfc/driver.py index 32a65f0542..2f9d4587a8 100644 --- a/tsfc/driver.py +++ b/tsfc/driver.py @@ -57,7 +57,6 @@ def compile_form(form, prefix="form", parameters=None, interface=None, diagonal= :returns: list of kernels """ cpu_time = time.time() - assert isinstance(form, Form) GREEN = "\033[1;37;32m%s\033[0m" From f10b66aed1f20325af368b82e8d0052f34868edf Mon Sep 17 00:00:00 2001 From: India Marsden Date: Thu, 27 Feb 2025 13:42:05 +0000 Subject: [PATCH 28/35] Fix small issue in interpolation --- firedrake/interpolation.py | 4 +++- tests/firedrake/regression/test_fuse.py | 3 ++- 2 files changed, 5 insertions(+), 2 deletions(-) diff --git a/firedrake/interpolation.py b/firedrake/interpolation.py index 05d3438318..f0cc8eb1d5 100644 --- a/firedrake/interpolation.py +++ b/firedrake/interpolation.py @@ -1281,7 +1281,9 @@ def rebuild_dg(element, expr, rt_var_name): # weights.append(w) # assert len(weights) == num_points # but for now we just fix the values to what we know works: - if element.degree != 0 or not isinstance(element.cell, FIAT.reference_element.Point): + + # not isinstance(element.cell, FIAT.reference_element.Point) + if element.degree != 0 or element.cell.get_dimension() != 0: raise NotImplementedError("Cross mesh interpolation only implemented for P0DG on vertex cells.") num_points = 1 weights = [1.]*num_points diff --git a/tests/firedrake/regression/test_fuse.py b/tests/firedrake/regression/test_fuse.py index cbae8cdcd4..45ada80ea1 100644 --- a/tests/firedrake/regression/test_fuse.py +++ b/tests/firedrake/regression/test_fuse.py @@ -15,7 +15,8 @@ def test_poisson_analytic(params, degree, quadrilateral): @pytest.mark.parametrize(['conv_num', 'degree'], [(p, d) - for p, d in zip([1.8, 2.8, 3.8], [1, 2, 3])]) + for p, d in zip([1.8, 2.8], [1, 2])]) + # 3.8, 3 def test_helmholtz(conv_num, degree): diff = np.array([helmholtz(i, degree=degree)[0] for i in range(3, 6)]) print("l2 error norms:", diff) From 8254f33eff9f888a39c6fd39a65da0a6d917376b Mon Sep 17 00:00:00 2001 From: India Marsden Date: Tue, 6 May 2025 12:07:03 +0100 Subject: [PATCH 29/35] change quad orderin based on given cell, add cg3 to fuse tests --- firedrake/cython/dmcommon.pyx | 24 ++++++++++++++---------- firedrake/mesh.py | 23 +++++++++++++++-------- tests/firedrake/regression/test_fuse.py | 7 ++++--- 3 files changed, 33 insertions(+), 21 deletions(-) diff --git a/firedrake/cython/dmcommon.pyx b/firedrake/cython/dmcommon.pyx index 23bf9157c9..4b8728a338 100644 --- a/firedrake/cython/dmcommon.pyx +++ b/firedrake/cython/dmcommon.pyx @@ -368,15 +368,18 @@ cdef inline PetscInt _reorder_plex_cone(PETSc.DM dm, # 0 2 # | | # +---3---+ - plex_cone_new[0] = plex_cone_old[0] - plex_cone_new[1] = plex_cone_old[2] - plex_cone_new[2] = plex_cone_old[1] - plex_cone_new[3] = plex_cone_old[3] - # print("SWITCHING") - #plex_cone_new[0] = plex_cone_old[0] - #plex_cone_new[1] = plex_cone_old[3] - #plex_cone_new[2] = plex_cone_old[2] - #plex_cone_new[3] = plex_cone_old[1] + if "fuse" in dm.name: + # FUSE rules + plex_cone_new[0] = plex_cone_old[0] + plex_cone_new[1] = plex_cone_old[3] + plex_cone_new[2] = plex_cone_old[2] + plex_cone_new[3] = plex_cone_old[1] + else: + # UFC rules + plex_cone_new[0] = plex_cone_old[0] + plex_cone_new[1] = plex_cone_old[2] + plex_cone_new[2] = plex_cone_old[1] + plex_cone_new[3] = plex_cone_old[3] elif dm.getCellType(p) == PETSc.DM.PolytopeType.HEXAHEDRON: # UFCHexahedron: +-------+ +-------+ # /. | / 5 /| @@ -1180,7 +1183,7 @@ def entity_orientations(mesh, entity_cone_map[i] = entity_cone_list[i] for i in range(len(entity_cone_list_offset)): entity_cone_map_offset[i] = entity_cone_list_offset[i] - # + dm = mesh.topology_dm dim = dm.getDimension() numCells = cell_closure.shape[0] @@ -1517,6 +1520,7 @@ def get_cell_nodes(mesh, perm_offset += ceil_ndofs[i] * num_orientations_c[i] else: # FInAT element must eventually add entity_permutations() method + if extruded_periodic_1_layer: for j in range(ceil_ndofs[i]): cell_nodes[cell, flat_index[k]] = off + j % offset[flat_index[k]] diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 2dcc3b51eb..dabef608e4 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -335,6 +335,7 @@ def local_facet_orientation_dat(self): local_facet_start = offsets[-3] local_facet_end = offsets[-2] map_from_cell_to_facet_orientations = self.mesh.entity_orientations[:, local_facet_start:local_facet_end] + # Make output data; # this is a map from an exterior/interior facet to the corresponding local facet orientation/orientations. # Halo data are required by design, but not actually used. @@ -1264,21 +1265,19 @@ def cell_closure(self): assert tdim == cell.topological_dimension() # plex.viewFromOptions("-dm_view") closure, _ = plex.getTransitiveClosure(0) + if hasattr(cell, "to_fiat"): + plex.setName('firedrake_default_topology_fuse') + if self.submesh_parent is not None: return dmcommon.submesh_create_cell_closure_cell_submesh(plex, self.submesh_parent.topology_dm, cell_numbering, self.submesh_parent._cell_numbering, self.submesh_parent.cell_closure) - elif cell.is_simplex(): - topology = FIAT.ufc_cell(cell).get_topology() - entity_per_cell = np.zeros(len(topology), dtype=IntType) - for d, ents in topology.items(): - entity_per_cell[d] = len(ents) - return dmcommon.closure_ordering(plex, vertex_numbering, - cell_numbering, entity_per_cell) # elif hasattr(cell, "to_fiat"): - # # TODO + # breakpoint() + # plex.setName('firedrake_default_topology_fuse') + # # TODO find better way of branching here # topology = cell.to_fiat().topology # # entity_per_cell = np.zeros(len(topology), dtype=IntType) # # for d, ents in topology.items(): @@ -1287,7 +1286,15 @@ def cell_closure(self): # # cell_numbering, entity_per_cell) # topology = FIAT.ufc_cell(cell).get_topology() # closureSize = sum([len(ents) for _, ents in topology.items()]) + # breakpoint() # return dmcommon.create_cell_closure(plex, cell_numbering, closureSize) + elif cell.is_simplex(): + topology = FIAT.ufc_cell(cell).get_topology() + entity_per_cell = np.zeros(len(topology), dtype=IntType) + for d, ents in topology.items(): + entity_per_cell[d] = len(ents) + return dmcommon.closure_ordering(plex, vertex_numbering, + cell_numbering, entity_per_cell) elif cell.cellname() == "quadrilateral": from firedrake_citations import Citations Citations().register("Homolya2016") diff --git a/tests/firedrake/regression/test_fuse.py b/tests/firedrake/regression/test_fuse.py index 45ada80ea1..34726410df 100644 --- a/tests/firedrake/regression/test_fuse.py +++ b/tests/firedrake/regression/test_fuse.py @@ -1,5 +1,6 @@ from test_helmholtz import helmholtz from test_poisson_strong_bcs import run_test +from fuse.cells import firedrake_triangle import pytest import numpy as np @@ -15,9 +16,9 @@ def test_poisson_analytic(params, degree, quadrilateral): @pytest.mark.parametrize(['conv_num', 'degree'], [(p, d) - for p, d in zip([1.8, 2.8], [1, 2])]) - # 3.8, 3 -def test_helmholtz(conv_num, degree): + for p, d in zip([1.8, 2.8, 3.8], [1, 2, 3])]) +def test_helmholtz(mocker, conv_num, degree): + # mocker.patch('firedrake.mesh.as_cell', return_value=firedrake_triangle().to_ufl("triangle")) diff = np.array([helmholtz(i, degree=degree)[0] for i in range(3, 6)]) print("l2 error norms:", diff) conv = np.log2(diff[:-1] / diff[1:]) From 0f155d50995665ef35a59f6569183a6d7a8aba96 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Tue, 6 May 2025 13:12:37 +0100 Subject: [PATCH 30/35] DROP BEFORE MERGE --- .github/workflows/core.yml | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/.github/workflows/core.yml b/.github/workflows/core.yml index 3435fdc56c..e75c31d32c 100644 --- a/.github/workflows/core.yml +++ b/.github/workflows/core.yml @@ -189,6 +189,10 @@ jobs: --no-binary h5py \ --extra-index-url https://download.pytorch.org/whl/cpu \ './firedrake-repo[ci,docs]' + + pip install -I "fenics-ufl @ git+https://github.com/firedrakeproject/ufl.git@indiamai/integrate-fuse" + pip install -I "firedrake-fiat @ git+https://github.com/firedrakeproject/fiat.git@indiamai/integrate_fuse" + pip install -I "fuse-element @ git+https://github.com/firedrakeproject/fuse.git@indiamai/compare-topologies" firedrake-clean pip list From 1636fdacae4871a947d9ad396fbdaec3f0d63b70 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Wed, 7 May 2025 13:58:13 +0100 Subject: [PATCH 31/35] DROP BEFORE MERGE make makefile python3 --- Makefile | 22 +++++++++++----------- 1 file changed, 11 insertions(+), 11 deletions(-) diff --git a/Makefile b/Makefile index ef491f3efe..c0fdd96241 100644 --- a/Makefile +++ b/Makefile @@ -4,7 +4,7 @@ all: modules .PHONY: modules modules: @echo " Building extension modules" - @python setup.py build_ext --inplace > build.log 2>&1 || cat build.log + @python3 setup.py build_ext --inplace > build.log 2>&1 || cat build.log .PHONY: lint lint: srclint actionlint dockerlint @@ -20,18 +20,18 @@ endif .PHONY: srclint srclint: @echo " Linting firedrake" - @python -m flake8 $(FLAKE8_FORMAT) firedrake + @python3 -m flake8 $(FLAKE8_FORMAT) firedrake @echo " Linting firedrake scripts" - @python -m flake8 $(FLAKE8_FORMAT) firedrake/scripts --filename=* - @python -m flake8 $(FLAKE8_FORMAT) scripts --filename=* + @python3 -m flake8 $(FLAKE8_FORMAT) firedrake/scripts --filename=* + @python3 -m flake8 $(FLAKE8_FORMAT) scripts --filename=* @echo " Linting firedrake tests" - @python -m flake8 $(FLAKE8_FORMAT) tests + @python3 -m flake8 $(FLAKE8_FORMAT) tests @echo " Linting PyOP2" - @python -m flake8 $(FLAKE8_FORMAT) pyop2 + @python3 -m flake8 $(FLAKE8_FORMAT) pyop2 @echo " Linting PyOP2 scripts" - @python -m flake8 $(FLAKE8_FORMAT) pyop2/scripts --filename=* + @python3 -m flake8 $(FLAKE8_FORMAT) pyop2/scripts --filename=* @echo " Linting TSFC" - @python -m flake8 $(FLAKE8_FORMAT) tsfc + @python3 -m flake8 $(FLAKE8_FORMAT) tsfc .PHONY: actionlint actionlint: @@ -58,7 +58,7 @@ dockerlint: .PHONY: clean clean: @echo " Cleaning extension modules" - @python setup.py clean > /dev/null 2>&1 + @python3 setup.py clean > /dev/null 2>&1 @echo " RM firedrake/cython/dmplex.*.so" -@rm -f firedrake/cython/dmplex.so > /dev/null 2>&1 @echo " RM firedrake/cython/dmplex.c" @@ -117,9 +117,9 @@ check: .PHONY: durations durations: @echo " Generate timings to optimise pytest-split" - python -m pytest --store-durations -m "parallel[1] or not parallel" tests/ + python3 -m pytest --store-durations -m "parallel[1] or not parallel" tests/ # use ':' to ensure that only rank 0 writes to the durations file for nprocs in 2 3 4 6 7 8; do \ - mpiexec -n 1 python -m pytest --store-durations -m parallel[$${nprocs}] tests/ : \ + mpiexec -n 1 python3 -m pytest --store-durations -m parallel[$${nprocs}] tests/ : \ -n $$(( $${nprocs} - 1 )) pytest -m parallel[$${nprocs}] -q tests/ ; \ done From 63d74137d0f536b93495d5938ed679f802160af3 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Fri, 16 May 2025 16:30:39 +0100 Subject: [PATCH 32/35] remove fuse special casing --- firedrake/cython/dmcommon.pyx | 26 ++++++++++++++------------ 1 file changed, 14 insertions(+), 12 deletions(-) diff --git a/firedrake/cython/dmcommon.pyx b/firedrake/cython/dmcommon.pyx index 4b8728a338..35bd01d8ba 100644 --- a/firedrake/cython/dmcommon.pyx +++ b/firedrake/cython/dmcommon.pyx @@ -368,18 +368,20 @@ cdef inline PetscInt _reorder_plex_cone(PETSc.DM dm, # 0 2 # | | # +---3---+ - if "fuse" in dm.name: - # FUSE rules - plex_cone_new[0] = plex_cone_old[0] - plex_cone_new[1] = plex_cone_old[3] - plex_cone_new[2] = plex_cone_old[2] - plex_cone_new[3] = plex_cone_old[1] - else: - # UFC rules - plex_cone_new[0] = plex_cone_old[0] - plex_cone_new[1] = plex_cone_old[2] - plex_cone_new[2] = plex_cone_old[1] - plex_cone_new[3] = plex_cone_old[3] + #if "fuse" in dm.name: + # FUSE rules + print("FUSE") + #plex_cone_new[0] = plex_cone_old[0] + #plex_cone_new[1] = plex_cone_old[3] + #plex_cone_new[2] = plex_cone_old[2] + #plex_cone_new[3] = plex_cone_old[1] + #else: + # UFC rules + #print("UFC") + plex_cone_new[0] = plex_cone_old[0] + plex_cone_new[1] = plex_cone_old[2] + plex_cone_new[2] = plex_cone_old[1] + plex_cone_new[3] = plex_cone_old[3] elif dm.getCellType(p) == PETSc.DM.PolytopeType.HEXAHEDRON: # UFCHexahedron: +-------+ +-------+ # /. | / 5 /| From cdec59f149dcc849859606be28f85cf823860790 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Thu, 29 May 2025 14:52:06 +0100 Subject: [PATCH 33/35] Getting Tetrahedra closer --- firedrake/cython/dmcommon.pyx | 42 +++++++++++++++++-- firedrake/mesh.py | 20 +++------ tests/firedrake/regression/test_fuse.py | 17 ++++++-- .../regression/test_steady_advection_3D.py | 7 +++- 4 files changed, 65 insertions(+), 21 deletions(-) diff --git a/firedrake/cython/dmcommon.pyx b/firedrake/cython/dmcommon.pyx index 35bd01d8ba..a04f39c988 100644 --- a/firedrake/cython/dmcommon.pyx +++ b/firedrake/cython/dmcommon.pyx @@ -468,7 +468,39 @@ cdef inline PetscInt _reorder_plex_closure(PETSc.DM dm, # 8 9 # \ / # 14 - raise NotImplementedError(f"Not implemented for {dm.getCellType(p)}") + # fuse tet + fiat_closure[0] = plex_closure[2 * 13] + fiat_closure[1] = plex_closure[2 * 11] + fiat_closure[2] = plex_closure[2 * 14] + fiat_closure[3] = plex_closure[2 * 12] + fiat_closure[4] = plex_closure[2 * 7] + fiat_closure[5] = plex_closure[2 * 8] + fiat_closure[6] = plex_closure[2 * 10] + fiat_closure[7] = plex_closure[2 * 6] + fiat_closure[8] = plex_closure[2 * 5] + fiat_closure[9] = plex_closure[2 * 9] + fiat_closure[10] = plex_closure[2 * 4] + fiat_closure[11] = plex_closure[2 * 1] + fiat_closure[12] = plex_closure[2 * 3] + fiat_closure[13] = plex_closure[2 * 2] + fiat_closure[14] = plex_closure[2 * 0] + # ufc tet fuse + #fiat_closure[0] = plex_closure[2 * 14] + #fiat_closure[1] = plex_closure[2 * 13] + #fiat_closure[2] = plex_closure[2 * 12] + #fiat_closure[3] = plex_closure[2 * 11] + #fiat_closure[4] = plex_closure[2 * 5] + #fiat_closure[5] = plex_closure[2 * 7] + #fiat_closure[6] = plex_closure[2 * 6] + #fiat_closure[7] = plex_closure[2 * 8] + #fiat_closure[8] = plex_closure[2 * 9] + #fiat_closure[9] = plex_closure[2 * 10] + #fiat_closure[10] = plex_closure[2 * 2] + #fiat_closure[11] = plex_closure[2 * 4] + #fiat_closure[12] = plex_closure[2 * 1] + #fiat_closure[13] = plex_closure[2 * 3] + #fiat_closure[14] = plex_closure[2 * 0] + # raise NotImplementedError(f"Not implemented for {dm.getCellType(p)}") elif dm.getCellType(p) == PETSc.DM.PolytopeType.QUADRILATERAL: # UFCQuadrilateral: 1---7---3 # | | @@ -1098,21 +1130,25 @@ cdef inline PetscInt _compute_orientation(PETSc.DM dm, const PetscInt *cone = NULL p = cell_closure[cell, e] + #print("p", p) if dm.getCellType(p) == PETSc.DM.PolytopeType.POINT: return 0 CHKERR(DMPlexGetConeSize(dm.dm, p, &coneSize)) CHKERR(DMPlexGetCone(dm.dm, p, &cone)) + #print(coneSize) if (entity_cone_map_offset[e + 1] - entity_cone_map_offset[e]) != coneSize: raise RuntimeError("FIAT entity cone size != plex point cone size") offset = entity_cone_map_offset[e] for i in range(coneSize): fiat_cone[i] = cell_closure[cell, entity_cone_map[offset + i]] if dm.getCellType(p) == PETSc.DM.PolytopeType.SEGMENT or \ - dm.getCellType(p) == PETSc.DM.PolytopeType.TRIANGLE or \ - dm.getCellType(p) == PETSc.DM.PolytopeType.TETRAHEDRON: + dm.getCellType(p) == PETSc.DM.PolytopeType.TRIANGLE: # UFCInterval <- PETSc.DM.PolytopeType.SEGMENT # UFCTriangle <- PETSc.DM.PolytopeType.TRIANGLE + return _compute_orientation_simplex(fiat_cone, cone, coneSize) + elif dm.getCellType(p) == PETSc.DM.PolytopeType.TETRAHEDRON: # UFCTetrahedron <- PETSc.DM.PolytopeType.TETRAHEDRON + # _reorder_plex_cone(dm, p, cone, plex_cone) return _compute_orientation_simplex(fiat_cone, cone, coneSize) elif dm.getCellType(p) == PETSc.DM.PolytopeType.QUADRILATERAL: # UFCQuadrilateral <- PETSc.DM.PolytopeType.QUADRILATERAL diff --git a/firedrake/mesh.py b/firedrake/mesh.py index dabef608e4..6265281ed4 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -1274,20 +1274,12 @@ def cell_closure(self): cell_numbering, self.submesh_parent._cell_numbering, self.submesh_parent.cell_closure) - # elif hasattr(cell, "to_fiat"): - # breakpoint() - # plex.setName('firedrake_default_topology_fuse') - # # TODO find better way of branching here - # topology = cell.to_fiat().topology - # # entity_per_cell = np.zeros(len(topology), dtype=IntType) - # # for d, ents in topology.items(): - # # entity_per_cell[d] = len(ents) - # # return dmcommon.closure_ordering(plex, vertex_numbering, - # # cell_numbering, entity_per_cell) - # topology = FIAT.ufc_cell(cell).get_topology() - # closureSize = sum([len(ents) for _, ents in topology.items()]) - # breakpoint() - # return dmcommon.create_cell_closure(plex, cell_numbering, closureSize) + elif hasattr(cell, "to_fiat") and cell.cellname() == "tetrahedron": + plex.setName('firedrake_default_topology_fuse') + # TODO find better way of branching here + topology = cell.to_fiat().topology + closureSize = sum([len(ents) for _, ents in topology.items()]) + return dmcommon.create_cell_closure(plex, cell_numbering, closureSize) elif cell.is_simplex(): topology = FIAT.ufc_cell(cell).get_topology() entity_per_cell = np.zeros(len(topology), dtype=IntType) diff --git a/tests/firedrake/regression/test_fuse.py b/tests/firedrake/regression/test_fuse.py index 34726410df..2aac60d398 100644 --- a/tests/firedrake/regression/test_fuse.py +++ b/tests/firedrake/regression/test_fuse.py @@ -1,6 +1,8 @@ from test_helmholtz import helmholtz from test_poisson_strong_bcs import run_test -from fuse.cells import firedrake_triangle +from test_steady_advection_3D import run_near_to_far +from fuse.cells import ufc_triangle +from firedrake import * import pytest import numpy as np @@ -11,16 +13,25 @@ for d in (1, 2, 3) for q in [False, True]]) def test_poisson_analytic(params, degree, quadrilateral): - assert (run_test(2, degree, parameters=params, quadrilateral=False) < 1.e-9) + assert (run_test(2, degree, parameters=params, quadrilateral=quadrilateral) < 1.e-9) @pytest.mark.parametrize(['conv_num', 'degree'], [(p, d) for p, d in zip([1.8, 2.8, 3.8], [1, 2, 3])]) def test_helmholtz(mocker, conv_num, degree): - # mocker.patch('firedrake.mesh.as_cell', return_value=firedrake_triangle().to_ufl("triangle")) + # mocker.patch('firedrake.mesh.as_cell', return_value=ufc_triangle().to_ufl("triangle")) diff = np.array([helmholtz(i, degree=degree)[0] for i in range(3, 6)]) print("l2 error norms:", diff) conv = np.log2(diff[:-1] / diff[1:]) print("convergence order:", conv) assert (np.array(conv) > conv_num).all() + + +def test_3D_advection(): + mesh = UnitTetrahedronMesh() + dg0 = FunctionSpace(mesh, "DG", 0) + dg1 = FunctionSpace(mesh, "DG", 1) + rt1 = FunctionSpace(mesh, "RT", 1) + + run_near_to_far(mesh, dg0, rt1) \ No newline at end of file diff --git a/tests/firedrake/regression/test_steady_advection_3D.py b/tests/firedrake/regression/test_steady_advection_3D.py index e90d15336e..8de1013532 100644 --- a/tests/firedrake/regression/test_steady_advection_3D.py +++ b/tests/firedrake/regression/test_steady_advection_3D.py @@ -10,7 +10,8 @@ @pytest.fixture(scope='module') def mesh(): - return UnitCubeMesh(3, 3, 3) + return UnitTetrahedronMesh() + # return UnitCubeMesh(3, 3, 3) @pytest.fixture(scope='module') @@ -50,6 +51,7 @@ def run_near_to_far(mesh, DG0, W): L = -inflow * inner(dot(u0, n), phi) * ds(3) # inflow at near wall + out = Function(DG0) solve(a == L, out) @@ -85,6 +87,9 @@ def run_up_to_down(mesh, DG1, W): a3 = inner(un * D, phi) * ds(5) # outflow at lower wall a = a1 + a2 + a3 + W = VectorFunctionSpace(mesh, DG1.ufl_element()) + X = assemble(interpolate(mesh.coordinates, W)) + print(X.dat.data) L = -inflow * inner(dot(u0, n), phi) * ds(6) # inflow at upper wall out = Function(DG1) From 8702d2ad4ea2765f48218a09d0cb5702178e41d7 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Fri, 30 May 2025 14:55:19 +0100 Subject: [PATCH 34/35] add simple 3d helmholtz test to fuse tests --- tests/firedrake/regression/test_fuse.py | 15 +++++++++------ .../regression/test_steady_advection_3D.py | 3 +-- 2 files changed, 10 insertions(+), 8 deletions(-) diff --git a/tests/firedrake/regression/test_fuse.py b/tests/firedrake/regression/test_fuse.py index 2aac60d398..ee92d6827c 100644 --- a/tests/firedrake/regression/test_fuse.py +++ b/tests/firedrake/regression/test_fuse.py @@ -28,10 +28,13 @@ def test_helmholtz(mocker, conv_num, degree): assert (np.array(conv) > conv_num).all() -def test_3D_advection(): - mesh = UnitTetrahedronMesh() - dg0 = FunctionSpace(mesh, "DG", 0) - dg1 = FunctionSpace(mesh, "DG", 1) - rt1 = FunctionSpace(mesh, "RT", 1) +@pytest.mark.parametrize(['conv_num', 'degree'], + [(p, d) + for p, d in zip([ 2.8, 3.8], [2, 3])]) +def test_helmholtz_3d(mocker, conv_num, degree): + diff = np.array([helmholtz(i, degree=degree, mesh=UnitCubeMesh(2 ** i, 2 ** i, 2 ** i))[0] for i in range(2, 4)]) + print("l2 error norms:", diff) + conv = np.log2(diff[:-1] / diff[1:]) + print("convergence order:", conv) + assert (np.array(conv) > conv_num).all() - run_near_to_far(mesh, dg0, rt1) \ No newline at end of file diff --git a/tests/firedrake/regression/test_steady_advection_3D.py b/tests/firedrake/regression/test_steady_advection_3D.py index 8de1013532..f146ffdeb8 100644 --- a/tests/firedrake/regression/test_steady_advection_3D.py +++ b/tests/firedrake/regression/test_steady_advection_3D.py @@ -10,8 +10,7 @@ @pytest.fixture(scope='module') def mesh(): - return UnitTetrahedronMesh() - # return UnitCubeMesh(3, 3, 3) + return UnitCubeMesh(3, 3, 3) @pytest.fixture(scope='module') From 7f55d3cca82b7a8adce6a2ccdb25e3c5b5e15458 Mon Sep 17 00:00:00 2001 From: India Marsden Date: Fri, 31 Oct 2025 10:43:24 +0000 Subject: [PATCH 35/35] use new as_cell --- firedrake/mesh.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/firedrake/mesh.py b/firedrake/mesh.py index 04c3c1efe5..0d60862196 100644 --- a/firedrake/mesh.py +++ b/firedrake/mesh.py @@ -49,7 +49,7 @@ # Only for docstring import mpi4py # noqa: F401 from finat.element_factory import as_fiat_cell -from ufl.cell import as_cell +from finat.ufl import as_cell __all__ = [