diff --git a/conda/meta.yaml b/conda/meta.yaml index 49fa83b90..1fb53ea24 100644 --- a/conda/meta.yaml +++ b/conda/meta.yaml @@ -22,8 +22,8 @@ requirements: - mpmath - plasmaboundaries >=0.1.8 - plotly - - brep_part_finder >=0.4.4 # [not win] - - brep_to_h5m >=0.3.1 # [not win] + - brep_part_finder >=0.5.0 # [not win] + - brep_to_h5m >=0.4.0 # [not win] - moab * nompi_tempest_* # - jupyter-cadquery not available on conda diff --git a/paramak/reactor.py b/paramak/reactor.py index 69b4b0502..fbe64f4dd 100644 --- a/paramak/reactor.py +++ b/paramak/reactor.py @@ -265,8 +265,6 @@ def export_dagmc_h5m( for shape in shapes_to_convert: tags.append(shape.name) - print(tags) - output_filename = export_solids_to_dagmc_h5m( solids=[shape.solid for shape in shapes_to_convert], filename=filename, diff --git a/paramak/utils.py b/paramak/utils.py index 6f335618e..8f469ec2d 100644 --- a/paramak/utils.py +++ b/paramak/utils.py @@ -96,6 +96,9 @@ def export_solids_to_dagmc_h5m( bounding_box_atol: float = 0.000001, tags: List[str] = None, ): + if verbose: + print("solids", solids, "\n") + print("tags", tags, "\n") if len(tags) != len(solids): msg = ( "When specifying tags then there must be one tag for " @@ -104,42 +107,43 @@ def export_solids_to_dagmc_h5m( ) raise ValueError(msg) + compound_expanded_tags = [] + # solids could contain compounds + for tag, solid in zip(tags, solids): + # before accessing the .val() check it exists + if hasattr(solid, "val"): + # if it is a compound then we need more material tags + if isinstance(solid.val(), cq.occ_impl.shapes.Compound): + additional_tags = [tag] * len(solid.val().Solids()) + compound_expanded_tags = compound_expanded_tags + additional_tags + else: + compound_expanded_tags.append(tag) + # if it is a compound then we need more material tags + elif isinstance(solid, cq.occ_impl.shapes.Compound): + additional_tags = [tag] * len(solid.Solids()) + compound_expanded_tags = compound_expanded_tags + additional_tags + else: + compound_expanded_tags.append(tag) + + if verbose: + print("compound_expanded_tags", compound_expanded_tags, "\n") + # a local import is used here as these packages need Moab to work - from brep_to_h5m import brep_to_h5m + from brep_to_h5m import mesh_brep, mesh_to_h5m_in_memory_method import brep_part_finder as bpf # saves the reactor as a Brep file with merged surfaces brep_shape = export_solids_to_brep_object(solids=solids) - brep_file_part_properties = bpf.get_brep_part_properties_from_shape(brep_shape) + brep_file_part_properties = bpf.get_part_properties_from_shapes(brep_shape) if verbose: - print("brep_file_part_properties", brep_file_part_properties) + print("brep_file_part_properties", brep_file_part_properties, "\n") - shape_properties = {} - for counter, solid in enumerate(solids): - sub_solid_descriptions = [] + shape_properties = bpf.get_part_properties_from_shapes(solids) + # for counter, solid in enumerate(solids): - # checks if the solid is a cq.Compound or not - if isinstance(solid, cq.occ_impl.shapes.Compound): - iterable_solids = solid.Solids() - else: - iterable_solids = solid.val().Solids() - - for sub_solid in iterable_solids: - part_bb = sub_solid.BoundingBox() - part_center = sub_solid.Center() - sub_solid_description = { - "volume": sub_solid.Volume(), - "center": (part_center.x, part_center.y, part_center.z), - "bounding_box": ( - (part_bb.xmin, part_bb.ymin, part_bb.zmin), - (part_bb.xmax, part_bb.ymax, part_bb.zmax), - ), - } - sub_solid_descriptions.append(sub_solid_description) - - shape_properties[tags[counter]] = sub_solid_descriptions + # shape_properties[counter] = bpf.get_part_properties_from_shape(solid) if verbose: print("shape_properties", shape_properties) @@ -148,34 +152,50 @@ def export_solids_to_dagmc_h5m( # using the volume, center, bounding box that we know about when creating the # CAD geometry in the first place - key_and_part_id = bpf.get_dict_of_part_ids( + brep_and_shape_part_ids = bpf.get_matching_part_ids( brep_part_properties=brep_file_part_properties, shape_properties=shape_properties, volume_atol=volume_atol, center_atol=center_atol, bounding_box_atol=bounding_box_atol, ) + if verbose: + print(f"brep_and_shape_part_ids={brep_and_shape_part_ids}") + + material_tags_in_brep_order = [] + for (brep_id, shape_id) in brep_and_shape_part_ids: + material_tags_in_brep_order.append(compound_expanded_tags[shape_id - 1]) if verbose: - print(f"key_and_part_id={key_and_part_id}") + print(f"material_tags_in_brep_order={material_tags_in_brep_order}") # gmsh requires an actual brep file to load tmp_brep_filename = mkstemp(suffix=".brep", prefix="paramak_")[1] brep_shape.exportBrep(tmp_brep_filename) - brep_to_h5m( + gmsh, volumes = mesh_brep( brep_filename=tmp_brep_filename, - volumes_with_tags=key_and_part_id, - h5m_filename=filename, min_mesh_size=min_mesh_size, max_mesh_size=max_mesh_size, - delete_intermediate_stl_files=True, ) - # temporary brep is deleted using os.remove - remove(tmp_brep_filename) + if verbose: + gmsh_filename = mkstemp(suffix=".msh", prefix="paramak_")[1] + print(f"written gmsh file to {gmsh_filename}") + gmsh.write(gmsh_filename) + + h5m_filename = mesh_to_h5m_in_memory_method( + volumes=volumes, + material_tags=material_tags_in_brep_order, + h5m_filename=filename, + ) + + if not verbose: + print(f"written brep file to {tmp_brep_filename}") + # temporary brep is deleted using os.remove + remove(tmp_brep_filename) - return filename + return h5m_filename def get_bounding_box(solid) -> Tuple[Tuple[float, float, float], Tuple[float, float, float]]: diff --git a/setup.cfg b/setup.cfg index c7e0d6e02..7ef783c5a 100644 --- a/setup.cfg +++ b/setup.cfg @@ -35,8 +35,8 @@ install_requires= plasmaboundaries >= 0.1.8 jupyter-client < 7 jupyter-cadquery >= 3.2.0 - brep_part_finder >= 0.4.4 - brep_to_h5m >= 0.3.1 + brep_part_finder >= 0.5.0 + brep_to_h5m >= 0.4.0 setuptools_scm [options.extras_require] diff --git a/tests/tests_h5m/test_reactor_export_h5m.py b/tests/tests_h5m/test_reactor_export_h5m.py index 60d5fb240..8407acf7e 100644 --- a/tests/tests_h5m/test_reactor_export_h5m.py +++ b/tests/tests_h5m/test_reactor_export_h5m.py @@ -97,7 +97,7 @@ def test_dagmc_h5m_export_with_graveyard(self): named in the resulting h5m file, includes the optional graveyard""" self.test_reactor_3.rotation_angle = 180 - self.test_reactor_3.export_dagmc_h5m("dagmc_reactor.h5m", include_graveyard={"size": 250}) + self.test_reactor_3.export_dagmc_h5m("dagmc_reactor.h5m", include_graveyard={"size": 250}, verbose=True) vols = di.get_volumes_from_h5m("dagmc_reactor.h5m") assert vols == [1, 2, 3, 4]