diff --git a/festim/exports/surface_quantity.py b/festim/exports/surface_quantity.py index 2594ed6dc..3553be9b1 100644 --- a/festim/exports/surface_quantity.py +++ b/festim/exports/surface_quantity.py @@ -12,7 +12,7 @@ class SurfaceQuantity: Attributes: field (festim.Species): species for which the surface flux is computed - surface (festim.SurfaceSubdomain1D): surface subdomain + surface (festim.SurfaceSubdomain): surface subdomain filename (str): name of the file to which the surface flux is exported t (list): list of time values data (list): list of values of the surface quantity @@ -47,10 +47,8 @@ def surface(self): @surface.setter def surface(self, value): - if not isinstance(value, (int, F.SurfaceSubdomain1D)) or isinstance( - value, bool - ): - raise TypeError("surface should be an int or F.SurfaceSubdomain1D") + if not isinstance(value, (int, F.SurfaceSubdomain)) or isinstance(value, bool): + raise TypeError("surface should be an int or F.SurfaceSubdomain") self._surface = value @property diff --git a/festim/heat_transfer_problem.py b/festim/heat_transfer_problem.py index f60d9900f..2ba68ce0c 100644 --- a/festim/heat_transfer_problem.py +++ b/festim/heat_transfer_problem.py @@ -112,26 +112,15 @@ def define_meshtags_and_measures(self): self.facet_meshtags = self.mesh.define_surface_meshtags() self.volume_meshtags = self.mesh.define_volume_meshtags() - elif isinstance(self.mesh, F.Mesh1D): + elif ( + isinstance(self.mesh, F.Mesh) + and self.facet_meshtags is None + and self.volume_meshtags is None + ): self.facet_meshtags, self.volume_meshtags = self.mesh.define_meshtags( surface_subdomains=self.surface_subdomains, volume_subdomains=self.volume_subdomains, ) - elif isinstance(self.mesh, F.Mesh): - # FIXME # refer to issue #647 - # create empty mesh tags for now - facet_indices = np.array([], dtype=np.int32) - facet_tags = np.array([], dtype=np.int32) - self.facet_meshtags = meshtags( - self.mesh.mesh, self.mesh.fdim, facet_indices, facet_tags - ) - - num_cells = self.mesh.mesh.topology.index_map(self.mesh.vdim).size_local - mesh_cell_indices = np.arange(num_cells, dtype=np.int32) - tags_volumes = np.full(num_cells, 1, dtype=np.int32) - self.volume_meshtags = meshtags( - self.mesh.mesh, self.mesh.vdim, mesh_cell_indices, tags_volumes - ) # check volume ids are unique vol_ids = [vol.id for vol in self.volume_subdomains] diff --git a/festim/hydrogen_transport_problem.py b/festim/hydrogen_transport_problem.py index bd5514b08..7207cac54 100644 --- a/festim/hydrogen_transport_problem.py +++ b/festim/hydrogen_transport_problem.py @@ -9,8 +9,6 @@ import tqdm.auto as tqdm import festim as F -from dolfinx.mesh import meshtags - class HydrogenTransportProblem: """ @@ -474,30 +472,16 @@ def define_meshtags_and_measures(self): self.facet_meshtags = self.mesh.define_surface_meshtags() self.volume_meshtags = self.mesh.define_volume_meshtags() - elif isinstance(self.mesh, F.Mesh1D): + elif ( + isinstance(self.mesh, F.Mesh) + and self.facet_meshtags is None + and self.volume_meshtags is None + ): self.facet_meshtags, self.volume_meshtags = self.mesh.define_meshtags( surface_subdomains=self.surface_subdomains, volume_subdomains=self.volume_subdomains, ) - elif isinstance(self.mesh, F.Mesh): - if not self.facet_meshtags: - # create empty facet_meshtags - facet_indices = np.array([], dtype=np.int32) - facet_tags = np.array([], dtype=np.int32) - self.facet_meshtags = meshtags( - self.mesh.mesh, self.mesh.fdim, facet_indices, facet_tags - ) - - if not self.volume_meshtags: - # create meshtags with all cells tagged as 1 - num_cells = self.mesh.mesh.topology.index_map(self.mesh.vdim).size_local - mesh_cell_indices = np.arange(num_cells, dtype=np.int32) - tags_volumes = np.full(num_cells, 1, dtype=np.int32) - self.volume_meshtags = meshtags( - self.mesh.mesh, self.mesh.vdim, mesh_cell_indices, tags_volumes - ) - # check volume ids are unique vol_ids = [vol.id for vol in self.volume_subdomains] if len(vol_ids) != len(np.unique(vol_ids)): diff --git a/festim/mesh/mesh.py b/festim/mesh/mesh.py index 819ac8e8c..83537330a 100644 --- a/festim/mesh/mesh.py +++ b/festim/mesh/mesh.py @@ -1,5 +1,7 @@ import ufl import dolfinx +import numpy as np +from dolfinx.mesh import meshtags class Mesh: @@ -52,3 +54,42 @@ def fdim(self): @property def n(self): return ufl.FacetNormal(self.mesh) + + def define_meshtags(self, surface_subdomains, volume_subdomains): + """Defines the facet and volume meshtags of the mesh + + Args: + surface_subdomains (list of festim.SufaceSubdomains): the surface subdomains of the model + volume_subdomains (list of festim.VolumeSubdomains): the volume subdomains of the model + + Returns: + dolfinx.mesh.MeshTags: the facet meshtags + dolfinx.mesh.MeshTags: the volume meshtags + """ + # find all cells in domain and mark them as 0 + num_cells = self.mesh.topology.index_map(self.vdim).size_local + mesh_cell_indices = np.arange(num_cells, dtype=np.int32) + tags_volumes = np.full(num_cells, 0, dtype=np.int32) + + # find all facets in domain and mark them as 0 + num_facets = self.mesh.topology.index_map(self.fdim).size_local + mesh_facet_indices = np.arange(num_facets, dtype=np.int32) + tags_facets = np.full(num_facets, 0, dtype=np.int32) + + for surf in surface_subdomains: + # find all facets in subdomain and mark them as surf.id + entities = surf.locate_boundary_facet_indices(self.mesh) + tags_facets[entities] = surf.id + + for vol in volume_subdomains: + # find all cells in subdomain and mark them as vol.id + entities = vol.locate_subdomain_entities(self.mesh, self.vdim) + tags_volumes[entities] = vol.id + + # define mesh tags + facet_meshtags = meshtags(self.mesh, self.fdim, mesh_facet_indices, tags_facets) + volume_meshtags = meshtags( + self.mesh, self.vdim, mesh_cell_indices, tags_volumes + ) + + return facet_meshtags, volume_meshtags diff --git a/festim/mesh/mesh_1d.py b/festim/mesh/mesh_1d.py index 706d29eec..d9c5ae030 100644 --- a/festim/mesh/mesh_1d.py +++ b/festim/mesh/mesh_1d.py @@ -4,7 +4,6 @@ import ufl import numpy as np import festim as F -from dolfinx.mesh import meshtags class Mesh1D(F.Mesh): @@ -74,44 +73,6 @@ def check_borders(self, volume_subdomains): raise ValueError("borders dont match domain borders") def define_meshtags(self, surface_subdomains, volume_subdomains): - """Defines the facet and volume meshtags of the mesh - - Args: - surface_subdomains (list of festim.SufaceSubdomains): the surface subdomains of the model - volume_subdomains (list of festim.VolumeSubdomains): the volume subdomains of the model - - Returns: - dolfinx.mesh.MeshTags: the facet meshtags - dolfinx.mesh.MeshTags: the volume meshtags - """ - facet_indices, tags_facets = [], [] - - # find all cells in domain and mark them as 0 - num_cells = self.mesh.topology.index_map(self.vdim).size_local - mesh_cell_indices = np.arange(num_cells, dtype=np.int32) - tags_volumes = np.full(num_cells, 0, dtype=np.int32) - - for surf in surface_subdomains: - facet_index = surf.locate_boundary_facet_indices(self.mesh, self.fdim) - facet_indices.append(facet_index) - tags_facets.append(surf.id) - - for vol in volume_subdomains: - # find all cells in subdomain and mark them as sub_dom.id - entities = vol.locate_subdomain_entities(self.mesh, self.vdim) - tags_volumes[entities] = vol.id - # check if all borders are defined self.check_borders(volume_subdomains) - - # dofs and tags need to be in np.in32 format for meshtags - facet_indices = np.array(facet_indices, dtype=np.int32) - tags_facets = np.array(tags_facets, dtype=np.int32) - - # define mesh tags - facet_meshtags = meshtags(self.mesh, self.fdim, facet_indices, tags_facets) - volume_meshtags = meshtags( - self.mesh, self.vdim, mesh_cell_indices, tags_volumes - ) - - return facet_meshtags, volume_meshtags + return super().define_meshtags(surface_subdomains, volume_subdomains) diff --git a/festim/subdomain/surface_subdomain.py b/festim/subdomain/surface_subdomain.py index 2ff373fbf..578545804 100644 --- a/festim/subdomain/surface_subdomain.py +++ b/festim/subdomain/surface_subdomain.py @@ -1,3 +1,7 @@ +import dolfinx.mesh +import numpy as np + + class SurfaceSubdomain: """ Surface subdomain class @@ -9,6 +13,24 @@ class SurfaceSubdomain: def __init__(self, id): self.id = id + def locate_boundary_facet_indices(self, mesh): + """Locates the dof of the surface subdomain within the function space + and return the index of the dof + + Args: + mesh (dolfinx.mesh.Mesh): the mesh of the simulation + + Returns: + index (np.array): the first value in the list of surface facet + indices of the subdomain + """ + fdim = mesh.topology.dim - 1 + # By default, all entities are included + indices = dolfinx.mesh.locate_entities_boundary( + mesh, fdim, lambda x: np.full(x.shape[1], True, dtype=bool) + ) + return indices + def find_surface_from_id(id: int, surfaces: list): """Returns the correct surface subdomain object from a list of surface ids diff --git a/festim/subdomain/surface_subdomain_1d.py b/festim/subdomain/surface_subdomain_1d.py index 23469b6ed..7393b32fa 100644 --- a/festim/subdomain/surface_subdomain_1d.py +++ b/festim/subdomain/surface_subdomain_1d.py @@ -23,19 +23,20 @@ def __init__(self, id, x) -> None: super().__init__(id) self.x = x - def locate_boundary_facet_indices(self, mesh, fdim): + def locate_boundary_facet_indices(self, mesh): """Locates the dof of the surface subdomain within the function space and return the index of the dof Args: mesh (dolfinx.mesh.Mesh): the mesh of the simulation - fdim (int): the dimension of the model facets Returns: index (np.array): the first value in the list of surface facet indices of the subdomain """ + assert mesh.geometry.dim == 1, "This method is only for 1D meshes" + fdim = 0 indices = dolfinx.mesh.locate_entities_boundary( mesh, fdim, lambda x: np.isclose(x[0], self.x) ) - return indices[0] + return indices diff --git a/festim/subdomain/volume_subdomain.py b/festim/subdomain/volume_subdomain.py index 6a7b28b76..3c476ceb8 100644 --- a/festim/subdomain/volume_subdomain.py +++ b/festim/subdomain/volume_subdomain.py @@ -1,3 +1,8 @@ +from dolfinx.mesh import locate_entities + +import numpy as np + + class VolumeSubdomain: """ Volume subdomain class @@ -10,6 +15,22 @@ def __init__(self, id, material): self.id = id self.material = material + def locate_subdomain_entities(self, mesh, vdim): + """Locates all cells in subdomain borders within domain + + Args: + mesh (dolfinx.mesh.Mesh): the mesh of the model + vdim (int): the dimension of the volumes of the mesh, + for 1D this is always 1 + + Returns: + entities (np.array): the entities of the subdomain + """ + # By default, all entities are included + # return array like x full of True + entities = locate_entities(mesh, vdim, lambda x: np.full(x.shape[1], True)) + return entities + def find_volume_from_id(id: int, volumes: list): """Returns the correct volume subdomain object from a list of volume ids diff --git a/test/system_tests/test_1_mobile_1_trap.py b/test/system_tests/test_1_mobile_1_trap.py index 3345f0f0d..b2b7eb68e 100644 --- a/test/system_tests/test_1_mobile_1_trap.py +++ b/test/system_tests/test_1_mobile_1_trap.py @@ -3,7 +3,7 @@ from dolfinx import fem import ufl from .tools import error_L2 -from dolfinx.mesh import meshtags, create_unit_square, create_unit_cube, locate_entities +from dolfinx.mesh import create_unit_square, create_unit_cube, locate_entities from mpi4py import MPI @@ -182,8 +182,8 @@ def u_exact_alt(mod): assert L2_error < 5e-4 -def test_1_mobile_1_trap_MMS_2D(): - """Tests that a steady simulation can be run in a 2D domain with +def test_1_mobile_1_trap_MMS_3D(): + """Tests that a steady simulation can be run in a 3D domain with 1 mobile and 1 trapped species""" def u_exact(mod): @@ -230,38 +230,23 @@ def v_exact(mod): my_model = F.HydrogenTransportProblem() my_model.mesh = F.Mesh(mesh=test_mesh_3d) - # create facet meshtags - boundaries = [ - (1, lambda x: np.isclose(x[0], 0)), - (2, lambda x: np.isclose(x[0], 1)), - ] - facet_indices, facet_markers = [], [] - fdim = test_mesh_3d.topology.dim - 1 - for marker, locator in boundaries: - facets = locate_entities(test_mesh_3d, fdim, locator) - facet_indices.append(facets) - facet_markers.append(np.full_like(facets, marker)) - facet_indices = np.hstack(facet_indices).astype(np.int32) - facet_markers = np.hstack(facet_markers).astype(np.int32) - sorted_facets = np.argsort(facet_indices) - my_facet_meshtags = meshtags( - test_mesh_3d, fdim, facet_indices[sorted_facets], facet_markers[sorted_facets] - ) + my_mat = F.Material(name="mat", D_0=D_0, E_D=E_D) + vol = F.VolumeSubdomain(id=1, material=my_mat) - # create volume meshtags - vdim = test_mesh_3d.topology.dim - num_cells = test_mesh_3d.topology.index_map(vdim).size_local - mesh_cell_indices = np.arange(num_cells, dtype=np.int32) - tags_volumes = np.full(num_cells, 1, dtype=np.int32) - my_volume_meshtags = meshtags(test_mesh_3d, vdim, mesh_cell_indices, tags_volumes) + class LefSurface(F.SurfaceSubdomain): + def locate_boundary_facet_indices(self, mesh): + fdim = mesh.topology.dim - 1 + indices = locate_entities(mesh, fdim, lambda x: np.isclose(x[0], 0)) + return indices - my_model.facet_meshtags = my_facet_meshtags - my_model.volume_meshtags = my_volume_meshtags + class RightSurface(F.SurfaceSubdomain): + def locate_boundary_facet_indices(self, mesh): + fdim = mesh.topology.dim - 1 + indices = locate_entities(mesh, fdim, lambda x: np.isclose(x[0], 1)) + return indices - my_mat = F.Material(name="mat", D_0=D_0, E_D=E_D) - vol = F.VolumeSubdomain(id=1, material=my_mat) - left = F.SurfaceSubdomain(id=1) - right = F.SurfaceSubdomain(id=2) + left = LefSurface(id=1) + right = RightSurface(id=2) my_model.subdomains = [vol, left, right] @@ -311,8 +296,8 @@ def v_exact(mod): assert L2_error_trapped < 9e-03 -def test_1_mobile_1_trap_MMS_3D(): - """Tests that a steady simulation can be run in a 3D domain with +def test_1_mobile_1_trap_MMS_2D(): + """Tests that a steady simulation can be run in a 2D domain with 1 mobile and 1 trapped species""" def u_exact(mod): @@ -359,38 +344,22 @@ def v_exact(mod): my_model = F.HydrogenTransportProblem() my_model.mesh = F.Mesh(mesh=test_mesh_2d) - # create facet meshtags - boundaries = [ - (1, lambda x: np.isclose(x[0], 0)), - (2, lambda x: np.isclose(x[0], 1)), - ] - facet_indices, facet_markers = [], [] - fdim = test_mesh_2d.topology.dim - 1 - for marker, locator in boundaries: - facets = locate_entities(test_mesh_2d, fdim, locator) - facet_indices.append(facets) - facet_markers.append(np.full_like(facets, marker)) - facet_indices = np.hstack(facet_indices).astype(np.int32) - facet_markers = np.hstack(facet_markers).astype(np.int32) - sorted_facets = np.argsort(facet_indices) - my_facet_meshtags = meshtags( - test_mesh_2d, fdim, facet_indices[sorted_facets], facet_markers[sorted_facets] - ) - - # create volume meshtags - vdim = test_mesh_2d.topology.dim - num_cells = test_mesh_2d.topology.index_map(vdim).size_local - mesh_cell_indices = np.arange(num_cells, dtype=np.int32) - tags_volumes = np.full(num_cells, 1, dtype=np.int32) - my_volume_meshtags = meshtags(test_mesh_2d, vdim, mesh_cell_indices, tags_volumes) + class LefSurface(F.SurfaceSubdomain): + def locate_boundary_facet_indices(self, mesh): + fdim = mesh.topology.dim - 1 + indices = locate_entities(mesh, fdim, lambda x: np.isclose(x[0], 0)) + return indices - my_model.facet_meshtags = my_facet_meshtags - my_model.volume_meshtags = my_volume_meshtags + class RightSurface(F.SurfaceSubdomain): + def locate_boundary_facet_indices(self, mesh): + fdim = mesh.topology.dim - 1 + indices = locate_entities(mesh, fdim, lambda x: np.isclose(x[0], 1)) + return indices my_mat = F.Material(name="mat", D_0=D_0, E_D=E_D) vol = F.VolumeSubdomain(id=1, material=my_mat) - left = F.SurfaceSubdomain(id=1) - right = F.SurfaceSubdomain(id=2) + left = LefSurface(id=1) + right = RightSurface(id=2) my_model.subdomains = [vol, left, right] diff --git a/test/system_tests/test_1_mobile_species.py b/test/system_tests/test_1_mobile_species.py index a324c7043..e5b9f3bfb 100644 --- a/test/system_tests/test_1_mobile_species.py +++ b/test/system_tests/test_1_mobile_species.py @@ -3,7 +3,7 @@ from dolfinx import fem import ufl from .tools import error_L2 -from dolfinx.mesh import meshtags, create_unit_square, create_unit_cube, locate_entities +from dolfinx.mesh import create_unit_square, create_unit_cube, locate_entities from mpi4py import MPI test_mesh_1d = F.Mesh1D(np.linspace(0, 1, 10000)) @@ -159,39 +159,22 @@ def u_exact(mod): my_model = F.HydrogenTransportProblem() my_model.mesh = F.Mesh(mesh=test_mesh_2d) - # create facet meshtags - boundaries = [ - (1, lambda x: np.isclose(x[0], 0)), - (2, lambda x: np.isclose(x[0], 1)), - ] - facet_indices, facet_markers = [], [] - fdim = test_mesh_2d.topology.dim - 1 - for marker, locator in boundaries: - facets = locate_entities(test_mesh_2d, fdim, locator) - facet_indices.append(facets) - facet_markers.append(np.full_like(facets, marker)) - facet_indices = np.hstack(facet_indices).astype(np.int32) - facet_markers = np.hstack(facet_markers).astype(np.int32) - sorted_facets = np.argsort(facet_indices) - my_facet_meshtags = meshtags( - test_mesh_2d, fdim, facet_indices[sorted_facets], facet_markers[sorted_facets] - ) - - # create volume meshtags - vdim = test_mesh_2d.topology.dim - num_cells = test_mesh_2d.topology.index_map(vdim).size_local - mesh_cell_indices = np.arange(num_cells, dtype=np.int32) - tags_volumes = np.full(num_cells, 1, dtype=np.int32) - my_volume_meshtags = meshtags(test_mesh_2d, vdim, mesh_cell_indices, tags_volumes) - - my_model.facet_meshtags = my_facet_meshtags - my_model.volume_meshtags = my_volume_meshtags + class LefSurface(F.SurfaceSubdomain): + def locate_boundary_facet_indices(self, mesh): + fdim = mesh.topology.dim - 1 + indices = locate_entities(mesh, fdim, lambda x: np.isclose(x[0], 0)) + return indices + + class RightSurface(F.SurfaceSubdomain): + def locate_boundary_facet_indices(self, mesh): + fdim = mesh.topology.dim - 1 + indices = locate_entities(mesh, fdim, lambda x: np.isclose(x[0], 1)) + return indices my_mat = F.Material(name="mat", D_0=D_0, E_D=E_D) vol = F.VolumeSubdomain(id=1, material=my_mat) - left = F.SurfaceSubdomain(id=1) - right = F.SurfaceSubdomain(id=2) - + left = LefSurface(id=1) + right = RightSurface(id=2) my_model.subdomains = [vol, left, right] H = F.Species("H") @@ -242,38 +225,22 @@ def u_exact(mod): my_model = F.HydrogenTransportProblem() my_model.mesh = F.Mesh(mesh=test_mesh_3d) - # create facet meshtags - boundaries = [ - (1, lambda x: np.isclose(x[0], 0)), - (2, lambda x: np.isclose(x[0], 1)), - ] - facet_indices, facet_markers = [], [] - fdim = test_mesh_3d.topology.dim - 1 - for marker, locator in boundaries: - facets = locate_entities(test_mesh_3d, fdim, locator) - facet_indices.append(facets) - facet_markers.append(np.full_like(facets, marker)) - facet_indices = np.hstack(facet_indices).astype(np.int32) - facet_markers = np.hstack(facet_markers).astype(np.int32) - sorted_facets = np.argsort(facet_indices) - my_facet_meshtags = meshtags( - test_mesh_3d, fdim, facet_indices[sorted_facets], facet_markers[sorted_facets] - ) - - # create volume meshtags - vdim = test_mesh_3d.topology.dim - num_cells = test_mesh_3d.topology.index_map(vdim).size_local - mesh_cell_indices = np.arange(num_cells, dtype=np.int32) - tags_volumes = np.full(num_cells, 1, dtype=np.int32) - my_volume_meshtags = meshtags(test_mesh_3d, vdim, mesh_cell_indices, tags_volumes) - - my_model.facet_meshtags = my_facet_meshtags - my_model.volume_meshtags = my_volume_meshtags + class LefSurface(F.SurfaceSubdomain): + def locate_boundary_facet_indices(self, mesh): + fdim = mesh.topology.dim - 1 + indices = locate_entities(mesh, fdim, lambda x: np.isclose(x[0], 0)) + return indices + + class RightSurface(F.SurfaceSubdomain): + def locate_boundary_facet_indices(self, mesh): + fdim = mesh.topology.dim - 1 + indices = locate_entities(mesh, fdim, lambda x: np.isclose(x[0], 1)) + return indices my_mat = F.Material(name="mat", D_0=D_0, E_D=E_D) vol = F.VolumeSubdomain(id=1, material=my_mat) - left = F.SurfaceSubdomain(id=1) - right = F.SurfaceSubdomain(id=2) + left = LefSurface(id=1) + right = RightSurface(id=2) my_model.subdomains = [vol, left, right] diff --git a/test/test_subdomains.py b/test/test_subdomains.py index 411a0eeed..2da15c72b 100644 --- a/test/test_subdomains.py +++ b/test/test_subdomains.py @@ -1,6 +1,8 @@ import numpy as np import festim as F import pytest +import dolfinx.mesh +from mpi4py import MPI def test_different_surface_ids(): @@ -180,3 +182,22 @@ def test_volume_subdomain_properties(): assert len(my_model.surface_subdomains) == 3 for subdomain in my_model.surface_subdomains: assert isinstance(subdomain, F.SurfaceSubdomain) + + +mesh1d = dolfinx.mesh.create_unit_interval(MPI.COMM_WORLD, nx=1) +mesh2d = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, nx=1, ny=1) +mesh3d = dolfinx.mesh.create_unit_cube(MPI.COMM_WORLD, nx=1, ny=1, nz=1) + + +@pytest.mark.parametrize("mesh,number_facets", [(mesh1d, 2), (mesh2d, 4), (mesh3d, 12)]) +def test_surface_subdomain_default(mesh, number_facets): + """Test that the default surface subdomain locates all boundary facets + + Args: + mesh (dolfinx.mesh.Mesh): a mesh + number_facets (int): the expected number of boundary facets + """ + my_subdomain = F.SurfaceSubdomain(id=1) + + indices = my_subdomain.locate_boundary_facet_indices(mesh) + assert len(indices) == number_facets diff --git a/test/test_surface_quantity.py b/test/test_surface_quantity.py index f264de8b6..a7234884e 100644 --- a/test/test_surface_quantity.py +++ b/test/test_surface_quantity.py @@ -18,7 +18,7 @@ def test_surface_flux_export_compute(): # define mesh ds measure facet_indices = np.array( - dummy_surface.locate_boundary_facet_indices(my_mesh.mesh, 0), + dummy_surface.locate_boundary_facet_indices(my_mesh.mesh), dtype=np.int32, ) tags_facets = np.array( @@ -152,10 +152,10 @@ def test_writer(tmp_path, value): def test_surface_setter_raises_TypeError(): """Test that a TypeError is raised when the surface is not a - F.SurfaceSubdomain1D""" + F.SurfaceSubdomain""" with pytest.raises( - TypeError, match="surface should be an int or F.SurfaceSubdomain1D" + TypeError, match="surface should be an int or F.SurfaceSubdomain" ): F.SurfaceQuantity( field=F.Species("H"),