From 18e029b19dca98dc51d781a74ee93b0741e448b5 Mon Sep 17 00:00:00 2001 From: Henrik Finsberg Date: Tue, 24 Sep 2024 09:46:46 +0200 Subject: [PATCH] Let entities be a list of tuples and add tests for on_boundary = True and False --- examples/mark_entities.py | 18 +++------- src/scifem/mesh.py | 28 +++++++--------- tests/test_mesh.py | 69 +++++++++++++++++++++++++++++++++++---- 3 files changed, 78 insertions(+), 37 deletions(-) diff --git a/examples/mark_entities.py b/examples/mark_entities.py index 5a6fb51..2d1dfd1 100644 --- a/examples/mark_entities.py +++ b/examples/mark_entities.py @@ -38,15 +38,12 @@ def inner(x): from scifem import create_meshtags -cell_tag = create_meshtags( - mesh, - mesh.topology.dim, - [{"tag": 1, "locator": left}, {"tag": 3, "locator": right}, {"tag": 7, "locator": inner}], -) +cell_tag = create_meshtags(mesh, mesh.topology.dim, [(1, left), (3, right), (7, inner)]) # Next we can plot these marked entities import pyvista + pyvista.start_xvfb() vtk_grid = dolfinx.plot.vtk_mesh(mesh, cell_tag.dim, cell_tag.indices) grid = pyvista.UnstructuredGrid(*vtk_grid) @@ -72,9 +69,7 @@ def top(x): return x[1] > 0.9 -facet_tags = create_meshtags( - mesh, mesh.topology.dim - 1, [{"tag": 2, "locator": top}, {"tag": 7, "locator": circle}] -) +facet_tags = create_meshtags(mesh, mesh.topology.dim - 1, [(2, top), (7, circle)]) facet_grid = dolfinx.plot.vtk_mesh(mesh, facet_tags.dim, facet_tags.indices) @@ -91,12 +86,7 @@ def top(x): # We can also exclude interior facets by adding `on_boundary: True` (by default this is set to False). boundary_facet_tags = create_meshtags( - mesh, - mesh.topology.dim - 1, - [ - {"tag": 2, "locator": top, "on_boundary": True}, - {"tag": 7, "locator": circle, "on_boundary": False}, - ], + mesh, mesh.topology.dim - 1, [(2, top, True), (7, circle, False)] ) diff --git a/src/scifem/mesh.py b/src/scifem/mesh.py index 3d61186..eb245b4 100644 --- a/src/scifem/mesh.py +++ b/src/scifem/mesh.py @@ -6,15 +6,11 @@ __all__ = ["create_meshtags"] -# Workaround to make id and locator required but on_boundary optional -# see https://peps.python.org/pep-0655/ -class _TaggedEntities(typing.TypedDict): - tag: int - locator: typing.Callable[[npt.NDArray[np.floating]], npt.NDArray[np.bool_]] - - -class TaggedEntities(_TaggedEntities, total=False): - on_boundary: bool +# (tag, locator, on_boundary) where on_boundary is optional +TaggedEntities = ( + tuple[int, typing.Callable[[npt.NDArray[np.floating]], npt.NDArray[np.bool_]]] + | tuple[int, typing.Callable[[npt.NDArray[np.floating]], npt.NDArray[np.bool_]], bool] +) def create_meshtags( @@ -27,12 +23,12 @@ def create_meshtags( Args: domain: A ``dolfinx.mesh.Mesh`` object dim: Dimension of the entities to mark - entities_list: A list of dictionaries with the following keys + entities_list: A list of tuples with the following elements: - - ``tag``: The tag to assign to the entities - - ``locator``: A function that takes a point and returns a boolean array + - ``index 0``: The tag to assign to the entities + - ``index 1``: A function that takes a point and returns a boolean array indicating whether the point is inside the entity - - ``on_boundary``: Optional, if True, the entities will be marked on the boundary + - ``index 2``: Optional, if True, the entities will be marked on the boundary Returns: A ``dolfinx.mesh.MeshTags`` object with the corresponding entities marked. @@ -60,9 +56,9 @@ def create_meshtags( # Concatenate and sort the arrays based on indices for tagged_entity in entities_list: - on_boundary = tagged_entity.get("on_boundary", False) - entities = locate_entities(on_boundary)(domain, dim, tagged_entity["locator"]) - markers[entities] = tagged_entity["tag"] + on_boundary = False if len(tagged_entity) == 2 else tagged_entity[2] + entities = locate_entities(on_boundary)(domain, dim, tagged_entity[1]) + markers[entities] = tagged_entity[0] facets = np.flatnonzero(markers != -1).astype(np.int32) return dolfinx.mesh.meshtags(domain, dim, facets, markers[facets]) diff --git a/tests/test_mesh.py b/tests/test_mesh.py index e5d1158..7006757 100644 --- a/tests/test_mesh.py +++ b/tests/test_mesh.py @@ -9,18 +9,18 @@ @pytest.mark.parametrize( "entities_list, set_values", [ - ([{"tag": 1, "locator": lambda x: x[0] <= 1.0}], {1}), - ([{"tag": 2, "locator": lambda x: x[0] <= 1.0}], {2}), + ([(1, lambda x: x[0] <= 1.0)], {1}), + ([(2, lambda x: x[0] <= 1.0)], {2}), ( [ - {"tag": 1, "locator": lambda x: x[0] <= 1.0}, - {"tag": 2, "locator": lambda x: x[0] >= 0.0}, + (1, lambda x: x[0] <= 1.0), + (2, lambda x: x[0] >= 0.0), ], {2}, ), ], ) -def test_create_meshtags_celltags_all(entities_list, set_values): +def test_create_celltags_celltags_all(entities_list, set_values): mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 3, 3) cell_tag = scifem.create_meshtags(mesh, mesh.topology.dim, entities_list=entities_list) @@ -32,14 +32,69 @@ def test_create_meshtags_celltags_all(entities_list, set_values): assert set(cell_tag.values) == set_values +@pytest.mark.parametrize( + "entities_list, set_values", + [ + ([(1, lambda x: x[0] <= 1.0, False)], {1}), + ([(2, lambda x: x[0] <= 1.0, False)], {2}), + ( + [ + (1, lambda x: x[0] <= 1.0, False), + (2, lambda x: x[0] >= 0.0, False), + ], + {2}, + ), + ], +) +def test_create_facet_tags_all_on_boundary_False(entities_list, set_values): + mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 3, 3) + facet_tag = scifem.create_meshtags(mesh, mesh.topology.dim - 1, entities_list=entities_list) + im = mesh.topology.index_map(mesh.topology.dim - 1) + + assert facet_tag.dim == mesh.topology.dim - 1 + assert facet_tag.indices.shape[0] == im.size_local + im.num_ghosts + assert facet_tag.values.shape[0] == im.size_local + im.num_ghosts + assert facet_tag.values.dtype == np.int32 + assert set(facet_tag.values) == set_values + + +@pytest.mark.parametrize( + "entities_list, set_values", + [ + ([(1, lambda x: x[0] <= 1.0, True)], {1}), + ([(2, lambda x: x[0] <= 1.0, True)], {2}), + ( + [ + (1, lambda x: x[0] <= 1.0, True), + (2, lambda x: x[0] >= 0.0, True), + ], + {2}, + ), + ], +) +def test_create_facet_tags_all_on_boundary_True(entities_list, set_values): + mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 3, 3) + facet_tag = scifem.create_meshtags(mesh, mesh.topology.dim - 1, entities_list=entities_list) + + mesh.topology.create_entities(1) + mesh.topology.create_connectivity(mesh.topology.dim - 1, mesh.topology.dim) + facet_indices = dolfinx.mesh.exterior_facet_indices(mesh.topology) + assert facet_tag.dim == mesh.topology.dim - 1 + assert facet_tag.indices.shape[0] == len(facet_indices) + assert facet_tag.values.shape[0] == len(facet_indices) + assert facet_tag.values.dtype == np.int32 + assert set(facet_tag.values) == set_values + + @pytest.mark.parametrize( "entities_list", [ - ([{"tag": 1, "locator": lambda x: x[0] < 0.0}]), + ([(1, lambda x: x[0] < 0.0)]), + ([(1, lambda x: x[0] < 0.0)]), ([]), ], ) -def test_create_meshtags_celltags_empty(entities_list): +def test_create_celltags_empty(entities_list): mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 3, 3) cell_tag = scifem.create_meshtags(mesh, mesh.topology.dim, entities_list=entities_list)