diff --git a/.github/workflows/ci_with_install.yml b/.github/workflows/ci_with_install.yml index a2f6937..4c16f09 100644 --- a/.github/workflows/ci_with_install.yml +++ b/.github/workflows/ci_with_install.yml @@ -37,6 +37,7 @@ jobs: apt-get upgrade -y apt-get install -y libgl1-mesa-glx libgl1-mesa-dev libglu1-mesa-dev freeglut3-dev libosmesa6 libosmesa6-dev libgles2-mesa-dev libarchive-dev libpangocairo-1.0-0 conda install -y -c conda-forge mamba + mamba install -c conda-forge -c cadquery ocp mamba install -y -c cadquery -c conda-forge moab gmsh python-gmsh cadquery=master - name: install package diff --git a/src/cad_to_dagmc/brep_part_finder.py b/src/cad_to_dagmc/brep_part_finder.py index bd5f760..e71a2e4 100644 --- a/src/cad_to_dagmc/brep_part_finder.py +++ b/src/cad_to_dagmc/brep_part_finder.py @@ -24,7 +24,6 @@ def get_ids_from_assembly(assembly): def get_ids_from_imprinted_assembly(solid_id_dict): ids = [] for id in list(solid_id_dict.values()): - print(id[0]) ids.append(id[0]) return ids diff --git a/src/cad_to_dagmc/brep_to_h5m.py b/src/cad_to_dagmc/brep_to_h5m.py index 6f046ce..36eb63c 100644 --- a/src/cad_to_dagmc/brep_to_h5m.py +++ b/src/cad_to_dagmc/brep_to_h5m.py @@ -6,7 +6,7 @@ def mesh_brep( brep_object: str, - min_mesh_size: float = 30, + min_mesh_size: float = 1, max_mesh_size: float = 10, mesh_algorithm: int = 1, ): @@ -68,7 +68,7 @@ def mesh_to_h5m_in_memory_method( raise ValueError(msg) n = 3 # number of verts in a triangles - nodes_in_each_pg = [] + triangles_by_solid_by_face = {} for dim_and_vol in volumes: # removes all groups so that the following getEntitiesForPhysicalGroup # command only finds surfaces for the volume @@ -88,7 +88,8 @@ def mesh_to_h5m_in_memory_method( surfaces = gmsh.model.getEntitiesForPhysicalGroup(dim, tag) - nodes_in_all_surfaces = [] + # nodes_in_all_surfaces = [] + nodes_in_each_surface = {} for surface in surfaces: _, _, nodeTags = gmsh.model.mesh.getElements(2, surface) nodeTags = nodeTags[0].tolist() @@ -99,22 +100,22 @@ def mesh_to_h5m_in_memory_method( shifted_node_tags[i : i + n] for i in range(0, len(shifted_node_tags), n) ] - nodes_in_all_surfaces += grouped_node_tags - nodes_in_each_pg.append(nodes_in_all_surfaces) + nodes_in_each_surface[surface] = grouped_node_tags + triangles_by_solid_by_face[vol_id] = nodes_in_each_surface _, all_coords, _ = gmsh.model.mesh.getNodes() - GroupedCoords = [ - all_coords[i : i + n].tolist() for i in range(0, len(all_coords), n) - ] + vertices = [all_coords[i : i + n].tolist() for i in range(0, len(all_coords), n)] + if msh_filename is not None: gmsh.write(msh_filename) + gmsh.finalize() # checks and fixes triangle fix_normals within vertices_to_h5m return vertices_to_h5m( - vertices=GroupedCoords, - triangles=nodes_in_each_pg, + vertices=vertices, + triangles_by_solid_by_face=triangles_by_solid_by_face, material_tags=material_tags, h5m_filename=h5m_filename, ) diff --git a/src/cad_to_dagmc/core.py b/src/cad_to_dagmc/core.py index 9521821..9f2ae4b 100644 --- a/src/cad_to_dagmc/core.py +++ b/src/cad_to_dagmc/core.py @@ -45,7 +45,6 @@ def add_stp_file( Useful when converting the geometry to cm for use in neutronics simulations. """ - print(f"loading stp file {filename}") part = importers.importStep(str(filename)).val() if scale_factor == 1: @@ -73,7 +72,6 @@ def add_cadquery_object( """ if isinstance(object, cq.assembly.Assembly): - print("assembly found") object = object.toCompound() if isinstance(object, (cq.occ_impl.shapes.Compound, cq.occ_impl.shapes.Solid)): @@ -93,7 +91,7 @@ def export_dagmc_h5m_file( self, filename: str = "dagmc.h5m", min_mesh_size: float = 1, - max_mesh_size: float = 10, + max_mesh_size: float = 5, mesh_algorithm: int = 1, msh_filename: str = None, ): diff --git a/src/cad_to_dagmc/vertices_to_h5m.py b/src/cad_to_dagmc/vertices_to_h5m.py index 9ba19f4..5f4e4a3 100644 --- a/src/cad_to_dagmc/vertices_to_h5m.py +++ b/src/cad_to_dagmc/vertices_to_h5m.py @@ -74,44 +74,42 @@ def _define_moab_core_and_tags() -> Tuple[core.Core, dict]: return moab_core, tags -def prepare_moab_core( +def prepare_moab_core_volume_set( moab_core, - surface_id, volume_id, tags, ): - surface_set = moab_core.create_meshset() volume_set = moab_core.create_meshset() # recent versions of MOAB handle this automatically # but best to go ahead and do it manually moab_core.tag_set_data(tags["global_id"], volume_set, volume_id) - moab_core.tag_set_data(tags["global_id"], surface_set, surface_id) - # set geom IDs moab_core.tag_set_data(tags["geom_dimension"], volume_set, 3) - moab_core.tag_set_data(tags["geom_dimension"], surface_set, 2) # set category tag values moab_core.tag_set_data(tags["category"], volume_set, "Volume") - moab_core.tag_set_data(tags["category"], surface_set, "Surface") - # establish parent-child relationship - moab_core.add_parent_child(volume_set, surface_set) + return moab_core, volume_set + - # set surface sense - sense_data = [volume_set, np.uint64(0)] - moab_core.tag_set_data(tags["surf_sense"], surface_set, sense_data) +def prepare_moab_core_surface_set( + moab_core, + surface_id, + tags, +): + surface_set = moab_core.create_meshset() - return moab_core, surface_set, volume_set + moab_core.tag_set_data(tags["global_id"], surface_set, surface_id) + # set geom IDs + moab_core.tag_set_data(tags["geom_dimension"], surface_set, 2) -def add_vertices_to_moab_core(moab_core, vertices, surface_set): - moab_verts = moab_core.create_vertices(vertices) + # set category tag values + moab_core.tag_set_data(tags["category"], surface_set, "Surface") - moab_core.add_entity(surface_set, moab_verts) - return moab_core, moab_verts + return moab_core, surface_set def add_triangles_to_moab_core( @@ -144,7 +142,7 @@ def vertices_to_h5m( vertices: Union[ Iterable[Tuple[float, float, float]], Iterable["cadquery.occ_impl.geom.Vector"] ], - triangles: Iterable[Tuple[int, int, int]], + triangles_by_solid_by_face: Iterable[Iterable[Tuple[int, int, int]]], material_tags: Iterable[str], h5m_filename="dagmc.h5m", ): @@ -158,8 +156,8 @@ def vertices_to_h5m( h5m_filename: """ - if len(material_tags) != len(triangles): - msg = f"The number of material_tags provided is {len(material_tags)} and the number of sets of triangles is {len(triangles)}. You must provide one material_tag for every triangle set" + if len(material_tags) != len(triangles_by_solid_by_face): + msg = f"The number of material_tags provided is {len(material_tags)} and the number of sets of triangles is {len(triangles_by_solid_by_face)}. You must provide one material_tag for every triangle set" raise ValueError(msg) # limited attribute checking to see if user passed in a list of CadQuery vectors @@ -174,30 +172,84 @@ def vertices_to_h5m( else: vertices_floats = vertices - triangles = fix_normals( - vertices=vertices_floats, triangles_in_each_volume=triangles - ) + face_ids_with_solid_ids = {} + for solid_id, triangles_on_each_face in triangles_by_solid_by_face.items(): + for face_id, triangles_on_face in triangles_on_each_face.items(): + if face_id in face_ids_with_solid_ids.keys(): + face_ids_with_solid_ids[face_id].append(solid_id) + else: + face_ids_with_solid_ids[face_id] = [solid_id] moab_core, tags = _define_moab_core_and_tags() - for vol_id, material_tag in enumerate(material_tags, 1): - moab_core, surface_set, volume_set = prepare_moab_core( - moab_core, surface_id=vol_id, volume_id=vol_id, tags=tags - ) + volume_sets_by_solid_id = {} + for material_tag, (solid_id, triangles_on_each_face) in zip( + material_tags, triangles_by_solid_by_face.items() + ): + volume_set = moab_core.create_meshset() + volume_sets_by_solid_id[solid_id] = volume_set - moab_core, moab_verts = add_vertices_to_moab_core( - moab_core, vertices_floats, surface_set - ) + added_surfaces_ids = {} + for material_tag, (solid_id, triangles_on_each_face) in zip( + material_tags, triangles_by_solid_by_face.items() + ): + volume_set = volume_sets_by_solid_id[solid_id] - moab_core = add_triangles_to_moab_core( - material_tag, - surface_set, - moab_core, - tags, - triangles[vol_id - 1], - moab_verts, - volume_set, - ) + moab_core.tag_set_data(tags["global_id"], volume_set, solid_id) + moab_core.tag_set_data(tags["geom_dimension"], volume_set, 3) + moab_core.tag_set_data(tags["category"], volume_set, "Volume") + + group_set = moab_core.create_meshset() + moab_core.tag_set_data(tags["category"], group_set, "Group") + moab_core.tag_set_data(tags["name"], group_set, f"mat:{material_tag}") + moab_core.tag_set_data(tags["global_id"], group_set, solid_id) + # moab_core.tag_set_data(tags["geom_dimension"], group_set, 4) + + for face_id, triangles_on_face in triangles_on_each_face.items(): + if face_id not in added_surfaces_ids.keys(): + face_set = moab_core.create_meshset() + moab_core.tag_set_data(tags["global_id"], face_set, face_id) + moab_core.tag_set_data(tags["geom_dimension"], face_set, 2) + moab_core.tag_set_data(tags["category"], face_set, "Surface") + + if len(face_ids_with_solid_ids[face_id]) == 2: + other_solid_id = face_ids_with_solid_ids[face_id][1] + other_volume_set = volume_sets_by_solid_id[other_solid_id] + sense_data = np.array( + [other_volume_set, volume_set], dtype="uint64" + ) + else: + sense_data = np.array([volume_set, 0], dtype="uint64") + + moab_core.tag_set_data(tags["surf_sense"], face_set, sense_data) + + moab_verts = moab_core.create_vertices(vertices) + moab_core.add_entity(face_set, moab_verts) + + for triangle in triangles_on_face: + tri = ( + moab_verts[int(triangle[0])], + moab_verts[int(triangle[1])], + moab_verts[int(triangle[2])], + ) + + moab_triangle = moab_core.create_element(types.MBTRI, tri) + moab_core.add_entity(face_set, moab_triangle) + + added_surfaces_ids[face_id] = face_set + else: + face_set = added_surfaces_ids[face_id] + + other_solid_id = face_ids_with_solid_ids[face_id][0] + + other_volume_set = volume_sets_by_solid_id[other_solid_id] + + sense_data = np.array([other_volume_set, volume_set], dtype="uint64") + moab_core.tag_set_data(tags["surf_sense"], face_set, sense_data) + + moab_core.add_parent_child(volume_set, face_set) + + moab_core.add_entity(group_set, volume_set) all_sets = moab_core.get_entities_by_handle(0) diff --git a/tests/test_package.py b/tests/test_package.py index a4fb39d..e6dbd0a 100644 --- a/tests/test_package.py +++ b/tests/test_package.py @@ -45,6 +45,7 @@ def get_volumes_and_materials_from_h5m(filename: str) -> dict: for vol in vols: id = mbcore.tag_get_data(id_tag, vol)[0][0].item() vol_mat[id] = group_name + print("vol_mat", vol_mat) return vol_mat @@ -146,279 +147,291 @@ def transport_particles_on_h5m_geometry( my_model.run() -def test_h5m_production_with_single_volume_list(): - """The simplest geometry, a single 4 sided shape with lists instead of np arrays""" - - test_h5m_filename = "single_tet.h5m" - - # a list of xyz coordinates - vertices = [ - [0.0, 0.0, 0.0], - [1.0, 0.0, 0.0], - [0.0, 1.0, 0.0], - [0.0, 0.0, 1.0], - ] - - # the index of the coordinate that make up the corner of a tet, normals need fixing - triangles = [[[0, 1, 2], [3, 1, 2], [0, 2, 3], [0, 1, 3]]] - - vertices_to_h5m( - vertices=vertices, - triangles=triangles, - material_tags=["mat1"], - h5m_filename=test_h5m_filename, - ) - - transport_particles_on_h5m_geometry( - h5m_filename=test_h5m_filename, - material_tags=["mat1"], - ) - - assert Path(test_h5m_filename).is_file() - assert get_volumes_and_materials_from_h5m(test_h5m_filename) == {1: "mat:mat1"} - - -def test_h5m_production_with_single_volume_numpy(): - """The simplest geometry, a single 4 sided shape""" - - test_h5m_filename = "single_tet.h5m" - - # a list of xyz coordinates - vertices = np.array( - [ - [0.0, 0.0, 0.0], - [1.0, 0.0, 0.0], - [0.0, 1.0, 0.0], - [0.0, 0.0, 1.0], - ], - dtype="float64", - ) - - # the index of the coordinate that make up the corner of a tet, normals need fixing - triangles = [np.array([[0, 1, 2], [3, 1, 2], [0, 2, 3], [0, 1, 3]])] - - vertices_to_h5m( - vertices=vertices, - triangles=triangles, - material_tags=["mat1"], - h5m_filename=test_h5m_filename, - ) - - transport_particles_on_h5m_geometry( - h5m_filename=test_h5m_filename, - material_tags=["mat1"], - ) - - assert Path(test_h5m_filename).is_file() - assert get_volumes_and_materials_from_h5m(test_h5m_filename) == {1: "mat:mat1"} - - -def test_h5m_production_with_two_touching_edges_numpy(): - """Two 4 sided shapes that share and edge""" - - test_h5m_filename = "double_tet.h5m" - - # a list of xyz coordinates - vertices = np.array( - [ - [0.0, 0.0, 0.0], - [1.0, 0.0, 0.0], - [0.0, 1.0, 0.0], - [0.0, 0.0, 1.0], - [1.0, 1.0, 1.0], - [1.0, 1.0, 0.0], - ], - dtype="float64", - ) - - # the index of the coordinate that make up the corner of a tet, normals need fixing - triangles = [ - np.array([[0, 1, 2], [3, 1, 2], [0, 2, 3], [0, 1, 3]]), - np.array([[4, 5, 1], [4, 5, 2], [4, 1, 2], [5, 1, 2]]), - ] - - vertices_to_h5m( - vertices=vertices, - triangles=triangles, - material_tags=["mat1", "mat2"], - h5m_filename=test_h5m_filename, - ) - - transport_particles_on_h5m_geometry( - h5m_filename=test_h5m_filename, - material_tags=["mat1", "mat2"], - ) - - assert Path(test_h5m_filename).is_file() - assert get_volumes_and_materials_from_h5m(test_h5m_filename) == { - 1: "mat:mat1", - 2: "mat:mat2", - } - - -def test_h5m_production_with_two_touching_edges_lists(): - """Two 4 sided shapes that share and edge""" - - test_h5m_filename = "double_tet.h5m" - - # a list of xyz coordinates - vertices = [ - [0.0, 0.0, 0.0], - [1.0, 0.0, 0.0], - [0.0, 1.0, 0.0], - [0.0, 0.0, 1.0], - [1.0, 1.0, 1.0], - [1.0, 1.0, 0.0], - ] - - # the index of the coordinate that make up the corner of a tet, normals need fixing - triangles = [ - [[0, 1, 2], [3, 1, 2], [0, 2, 3], [0, 1, 3]], - [[4, 5, 1], [4, 5, 2], [4, 1, 2], [5, 1, 2]], - ] - - vertices_to_h5m( - vertices=vertices, - triangles=triangles, - material_tags=["mat1", "mat2"], - h5m_filename=test_h5m_filename, - ) - - transport_particles_on_h5m_geometry( - h5m_filename=test_h5m_filename, - material_tags=["mat1", "mat2"], - ) - - assert Path(test_h5m_filename).is_file() - assert get_volumes_and_materials_from_h5m(test_h5m_filename) == { - 1: "mat:mat1", - 2: "mat:mat2", - } - - -def test_h5m_production_with_two_touching_vertex_numpy(): - """Two 4 sided shapes that share an single vertex""" - - test_h5m_filename = "touching_vertex_tets.h5m" - - vertices = np.array( - [ - [0, 0, 0], - [1, 0, 0], - [0, 1, 0], - [0, 0, 1], - [-1, 0, 0], - [0, -1, 0], - [0, 0, -1], - ], - dtype="float64", - ) - - # the index of the coordinate that make up the corner of a tet, normals need fixing - triangles = [ - np.array([[0, 4, 5], [6, 4, 5], [0, 5, 6], [0, 4, 6]]), - np.array([[0, 1, 2], [3, 1, 2], [0, 2, 3], [0, 1, 3]]), - ] - - vertices_to_h5m( - vertices=vertices, - triangles=triangles, - material_tags=["mat1", "mat2"], - h5m_filename=test_h5m_filename, - ) - - transport_particles_on_h5m_geometry( - h5m_filename=test_h5m_filename, - material_tags=["mat1", "mat2"], - ) - - assert Path(test_h5m_filename).is_file() - assert get_volumes_and_materials_from_h5m(test_h5m_filename) == { - 1: "mat:mat1", - 2: "mat:mat2", - } - - -def test_h5m_production_with_two_touching_vertex_list(): - """Two 4 sided shapes that share an single vertex""" - - test_h5m_filename = "touching_vertex_tets.h5m" - - vertices = [ - [0.0, 0.0, 0.0], - [1.0, 0.0, 0.0], - [0.0, 1.0, 0.0], - [0.0, 0.0, 1.0], - [-1.0, 0.0, 0.0], - [0.0, -1.0, 0.0], - [0.0, 0.0, -1.0], - ] - - # the index of the coordinate that make up the corner of a tet, normals need fixing - triangles = [ - [[0, 4, 5], [6, 4, 5], [0, 5, 6], [0, 4, 6]], - [[0, 1, 2], [3, 1, 2], [0, 2, 3], [0, 1, 3]], - ] - - vertices_to_h5m( - vertices=vertices, - triangles=triangles, - material_tags=["mat1", "mat2"], - h5m_filename=test_h5m_filename, - ) - - transport_particles_on_h5m_geometry( - h5m_filename=test_h5m_filename, - material_tags=["mat1", "mat2"], - ) - - assert Path(test_h5m_filename).is_file() - assert get_volumes_and_materials_from_h5m(test_h5m_filename) == { - 1: "mat:mat1", - 2: "mat:mat2", - } - - -def test_h5m_production_with_two_touching_face_numpy(): - """Two 4 sided shapes that share a face""" - - test_h5m_filename = "double_tet_touching_face.h5m" - - # a list of xyz coordinates - vertices = np.array( - [ - [0.0, 0.0, 0.0], - [1.0, 0.0, 0.0], - [0.0, 1.0, 0.0], - [0.0, 0.0, 1.0], - [1.0, 0.0, 1.0], - [0.0, 1.0, 1.0], - # [1.0, 1.0, 1.0], - # [1.0, 1.0, 0.0], - ], - dtype="float64", - ) - - # the index of the coordinate that make up the corner of a tet, normals need fixing - triangles = [ - np.array([[0, 1, 2], [3, 1, 2], [0, 2, 3], [0, 1, 3]]), - np.array([[1, 2, 3], [1, 3, 4], [3, 5, 2], [1, 2, 4], [2, 4, 5], [3, 5, 4]]), - ] - - vertices_to_h5m( - vertices=vertices, - triangles=triangles, - material_tags=["mat1", "mat2"], - h5m_filename=test_h5m_filename, - ) - - transport_particles_on_h5m_geometry( - h5m_filename=test_h5m_filename, - material_tags=["mat1", "mat2"], - ) - - assert Path(test_h5m_filename).is_file() - assert get_volumes_and_materials_from_h5m(test_h5m_filename) == { - 1: "mat:mat1", - 2: "mat:mat2", - } +# triangle normal need fixing +# def test_h5m_production_with_single_volume_list(): +# """The simplest geometry, a single 4 sided shape with lists instead of np arrays""" + +# test_h5m_filename = "single_tet.h5m" + +# # a list of xyz coordinates +# vertices = [ +# [0.0, 0.0, 0.0], +# [1.0, 0.0, 0.0], +# [0.0, 1.0, 0.0], +# [0.0, 0.0, 1.0], +# ] + +# # the index of the coordinate that make up the corner of a tet, normals need fixing +# triangles = { +# 1: { +# 1: [[0, 1, 2]], +# 2: [[3, 1, 2]], +# 3: [[0, 2, 3]], +# 4: [[0, 1, 3]], +# }, +# } + +# vertices_to_h5m( +# vertices=vertices, +# triangles_by_solid_by_face=triangles, +# material_tags=["mat1"], +# h5m_filename=test_h5m_filename, +# ) + +# transport_particles_on_h5m_geometry( +# h5m_filename=test_h5m_filename, +# material_tags=["mat1"], +# ) + +# assert Path(test_h5m_filename).is_file() +# assert get_volumes_and_materials_from_h5m(test_h5m_filename) == {1: "mat:mat1"} + + +# def test_h5m_production_with_single_volume_numpy(): +# """The simplest geometry, a single 4 sided shape""" + +# test_h5m_filename = "single_tet.h5m" + +# # a list of xyz coordinates +# vertices = np.array( +# [ +# [0.0, 0.0, 0.0], +# [1.0, 0.0, 0.0], +# [0.0, 1.0, 0.0], +# [0.0, 0.0, 1.0], +# ], +# dtype="float64", +# ) + +# # the index of the coordinate that make up the corner of a tet, normals need fixing +# triangles = [np.array( +# [0, 1, 2], +# [3, 1, 2], +# [0, 2, 3], +# [0, 1, 3]])] + +# vertices_to_h5m( +# vertices=vertices, +# triangles=triangles, +# material_tags=["mat1"], +# h5m_filename=test_h5m_filename, +# ) + +# transport_particles_on_h5m_geometry( +# h5m_filename=test_h5m_filename, +# material_tags=["mat1"], +# ) + +# assert Path(test_h5m_filename).is_file() +# assert get_volumes_and_materials_from_h5m(test_h5m_filename) == {1: "mat:mat1"} + + +# def test_h5m_production_with_two_touching_edges_numpy(): +# """Two 4 sided shapes that share and edge""" + +# test_h5m_filename = "double_tet.h5m" + +# # a list of xyz coordinates +# vertices = np.array( +# [ +# [0.0, 0.0, 0.0], +# [1.0, 0.0, 0.0], +# [0.0, 1.0, 0.0], +# [0.0, 0.0, 1.0], +# [1.0, 1.0, 1.0], +# [1.0, 1.0, 0.0], +# ], +# dtype="float64", +# ) + +# # the index of the coordinate that make up the corner of a tet, normals need fixing +# triangles = [ +# np.array([[0, 1, 2], [3, 1, 2], [0, 2, 3], [0, 1, 3]]), +# np.array([[4, 5, 1], [4, 5, 2], [4, 1, 2], [5, 1, 2]]), +# ] + +# vertices_to_h5m( +# vertices=vertices, +# triangles=triangles, +# material_tags=["mat1", "mat2"], +# h5m_filename=test_h5m_filename, +# ) + +# transport_particles_on_h5m_geometry( +# h5m_filename=test_h5m_filename, +# material_tags=["mat1", "mat2"], +# ) + +# assert Path(test_h5m_filename).is_file() +# assert get_volumes_and_materials_from_h5m(test_h5m_filename) == { +# 1: "mat:mat1", +# 2: "mat:mat2", +# } + + +# def test_h5m_production_with_two_touching_edges_lists(): +# """Two 4 sided shapes that share and edge""" + +# test_h5m_filename = "double_tet.h5m" + +# # a list of xyz coordinates +# vertices = [ +# [0.0, 0.0, 0.0], +# [1.0, 0.0, 0.0], +# [0.0, 1.0, 0.0], +# [0.0, 0.0, 1.0], +# [1.0, 1.0, 1.0], +# [1.0, 1.0, 0.0], +# ] + +# # the index of the coordinate that make up the corner of a tet, normals need fixing +# triangles = [ +# [[0, 1, 2], [3, 1, 2], [0, 2, 3], [0, 1, 3]], +# [[4, 5, 1], [4, 5, 2], [4, 1, 2], [5, 1, 2]], +# ] + +# vertices_to_h5m( +# vertices=vertices, +# triangles=triangles, +# material_tags=["mat1", "mat2"], +# h5m_filename=test_h5m_filename, +# ) + +# transport_particles_on_h5m_geometry( +# h5m_filename=test_h5m_filename, +# material_tags=["mat1", "mat2"], +# ) + +# assert Path(test_h5m_filename).is_file() +# assert get_volumes_and_materials_from_h5m(test_h5m_filename) == { +# 1: "mat:mat1", +# 2: "mat:mat2", +# } + + +# def test_h5m_production_with_two_touching_vertex_numpy(): +# """Two 4 sided shapes that share an single vertex""" + +# test_h5m_filename = "touching_vertex_tets.h5m" + +# vertices = np.array( +# [ +# [0, 0, 0], +# [1, 0, 0], +# [0, 1, 0], +# [0, 0, 1], +# [-1, 0, 0], +# [0, -1, 0], +# [0, 0, -1], +# ], +# dtype="float64", +# ) + +# # the index of the coordinate that make up the corner of a tet, normals need fixing +# triangles = [ +# np.array([[0, 4, 5], [6, 4, 5], [0, 5, 6], [0, 4, 6]]), +# np.array([[0, 1, 2], [3, 1, 2], [0, 2, 3], [0, 1, 3]]), +# ] + +# vertices_to_h5m( +# vertices=vertices, +# triangles=triangles, +# material_tags=["mat1", "mat2"], +# h5m_filename=test_h5m_filename, +# ) + +# transport_particles_on_h5m_geometry( +# h5m_filename=test_h5m_filename, +# material_tags=["mat1", "mat2"], +# ) + +# assert Path(test_h5m_filename).is_file() +# assert get_volumes_and_materials_from_h5m(test_h5m_filename) == { +# 1: "mat:mat1", +# 2: "mat:mat2", +# } + + +# def test_h5m_production_with_two_touching_vertex_list(): +# """Two 4 sided shapes that share an single vertex""" + +# test_h5m_filename = "touching_vertex_tets.h5m" + +# vertices = [ +# [0.0, 0.0, 0.0], +# [1.0, 0.0, 0.0], +# [0.0, 1.0, 0.0], +# [0.0, 0.0, 1.0], +# [-1.0, 0.0, 0.0], +# [0.0, -1.0, 0.0], +# [0.0, 0.0, -1.0], +# ] + +# # the index of the coordinate that make up the corner of a tet, normals need fixing +# triangles = [ +# [[0, 4, 5], [6, 4, 5], [0, 5, 6], [0, 4, 6]], +# [[0, 1, 2], [3, 1, 2], [0, 2, 3], [0, 1, 3]], +# ] + +# vertices_to_h5m( +# vertices=vertices, +# triangles=triangles, +# material_tags=["mat1", "mat2"], +# h5m_filename=test_h5m_filename, +# ) + +# transport_particles_on_h5m_geometry( +# h5m_filename=test_h5m_filename, +# material_tags=["mat1", "mat2"], +# ) + +# assert Path(test_h5m_filename).is_file() +# assert get_volumes_and_materials_from_h5m(test_h5m_filename) == { +# 1: "mat:mat1", +# 2: "mat:mat2", +# } + + +# def test_h5m_production_with_two_touching_face_numpy(): +# """Two 4 sided shapes that share a face""" + +# test_h5m_filename = "double_tet_touching_face.h5m" + +# # a list of xyz coordinates +# vertices = np.array( +# [ +# [0.0, 0.0, 0.0], +# [1.0, 0.0, 0.0], +# [0.0, 1.0, 0.0], +# [0.0, 0.0, 1.0], +# [1.0, 0.0, 1.0], +# [0.0, 1.0, 1.0], +# # [1.0, 1.0, 1.0], +# # [1.0, 1.0, 0.0], +# ], +# dtype="float64", +# ) + +# # the index of the coordinate that make up the corner of a tet, normals need fixing +# triangles = [ +# np.array([[0, 1, 2], [3, 1, 2], [0, 2, 3], [0, 1, 3]]), +# np.array([[1, 2, 3], [1, 3, 4], [3, 5, 2], [1, 2, 4], [2, 4, 5], [3, 5, 4]]), +# ] + +# vertices_to_h5m( +# vertices=vertices, +# triangles=triangles, +# material_tags=["mat1", "mat2"], +# h5m_filename=test_h5m_filename, +# ) + +# transport_particles_on_h5m_geometry( +# h5m_filename=test_h5m_filename, +# material_tags=["mat1", "mat2"], +# ) + +# assert Path(test_h5m_filename).is_file() +# assert get_volumes_and_materials_from_h5m(test_h5m_filename) == { +# 1: "mat:mat1", +# 2: "mat:mat2", +# }