Merge pull request #165 from svalinn/magnet-meshing

parastell/cubit_io.py

@@ -79,14 +79,16 @@ def export_cub5(filename, export_dir=""):
     cubit.cmd(f'save cub5 "{export_path}" overwrite')
 
 
-def export_mesh_cubit(filename, export_dir=""):
+def export_mesh_cubit(filename, export_dir="", delete_upon_export=True):
     """Exports Cubit mesh to H5M file format, first exporting to Exodus format
     via Coreform Cubit and converting to H5M via MOAB.
 
     Arguments:
         filename (str): name of H5M output file, excluding '.h5m' extension.
         export_dir (str): directory to which to export the H5M output file
             (defaults to empty string).
+        delete_upon_export (bool): delete the mesh from the Cubit instance
+            after exporting. Prevents inclusion of mesh in future exports.
     """
     init_cubit()
 
@@ -97,6 +99,11 @@ def export_mesh_cubit(filename, export_dir=""):
     subprocess.run(f"mbconvert {exo_path} {h5m_path}", shell=True)
     Path.unlink(exo_path)
 
+    # Delete any meshes present to prevent inclusion in future Cubit mesh
+    # exports
+    if delete_upon_export:
+        cubit.cmd(f"delete mesh volume all propagate")
+
 
 def export_dagmc_cubit_legacy(
     faceting_tolerance=None,
@@ -151,6 +158,7 @@ def export_dagmc_cubit_native(
     deviation_angle=5.0,
     filename="dagmc",
     export_dir="",
+    delete_upon_export=True,
 ):
     """Exports DAGMC neutronics H5M file of ParaStell components via native
     faceting method for Coreform Cubit.
@@ -165,6 +173,8 @@ def export_dagmc_cubit_native(
             (defaults to 'dagmc').
         export_dir (str): directory to which to export the DAGMC output file
             (defaults to empty string).
+        delete_upon_export (bool): delete the mesh from the Cubit instance
+            after exporting. Prevents inclusion of mesh in future exports.
     """
     init_cubit()
 
@@ -177,3 +187,8 @@ def export_dagmc_cubit_native(
 
     export_path = Path(export_dir) / Path(filename).with_suffix(".h5m")
     cubit.cmd(f'export cf_dagmc "{export_path}" overwrite')
+
+    # Delete any meshes present to prevent inclusion in future Cubit mesh
+    # exports
+    if delete_upon_export:
+        cubit.cmd(f"delete mesh volume all propagate")

parastell/invessel_build.py

@@ -10,6 +10,7 @@
 import pystell.read_vmec as read_vmec
 
 from . import log
+from . import cubit_io
 from .utils import (
     normalize,
     expand_list,
@@ -292,6 +293,12 @@ def merge_layer_surfaces(self):
                 )
                 prev_outer_surface_id = outer_surface_id
 
+    def import_step_cubit(self):
+        """Imports STEP files from in-vessel build into Coreform Cubit."""
+        for name, data in self.radial_build.radial_build.items():
+            vol_id = cubit_io.import_step_cubit(name, self.export_dir)
+            data["vol_id"] = vol_id
+
     def export_step(self, export_dir=""):
         """Export CAD solids as STEP files via CadQuery.
 

parastell/magnet_coils.py

@@ -238,6 +238,18 @@ def build_magnet_coils(self):
 
         self._cut_magnets()
 
+    def import_step_cubit(self):
+        """Import STEP file for magnet set into Coreform Cubit."""
+        first_vol_id = 1
+        if cubit_io.initialized:
+            first_vol_id += cubit.get_last_id("volume")
+
+        last_vol_id = cubit_io.import_step_cubit(
+            self.step_filename, self.export_dir
+        )
+
+        self.volume_ids = list(range(first_vol_id, last_vol_id + 1))
+
     def export_step(self, step_filename="magnet_set", export_dir=""):
         """Export CAD solids as a STEP file via CadQuery.
 
@@ -261,19 +273,28 @@ def export_step(self, step_filename="magnet_set", export_dir=""):
         )
         cq.exporters.export(coil_set, str(export_path))
 
-    def mesh_magnets(self):
-        """Creates tetrahedral mesh of magnet volumes via Coreform Cubit."""
-        self._logger.info("Generating tetrahedral mesh of magnet coils...")
+    def mesh_magnets(self, min_size=20.0, max_size=50.0, max_gradient=1.5):
+        """Creates tetrahedral mesh of magnet volumes via Coreform Cubit.
 
-        last_vol_id = cubit_io.import_step_cubit(
-            self.step_filename, self.export_dir
-        )
+        Arguments:
+            min_size (float): minimum size of mesh elements (defaults to 20.0).
+            max_size (float): maximum size of mesh elements (defaults to 50.0).
+            max_gradient (float): maximum transition in mesh element size
+                (defaults to 1.5).
+        """
+        self._logger.info("Generating tetrahedral mesh of magnet coils...")
 
-        self.volume_ids = range(1, last_vol_id + 1)
+        if not hasattr(self, "volume_ids"):
+            self.import_step_cubit()
 
-        for vol in self.volume_ids:
-            cubit.cmd(f"volume {vol} scheme tetmesh")
-            cubit.cmd(f"mesh volume {vol}")
+        volume_ids_str = " ".join(str(id) for id in self.volume_ids)
+        cubit.cmd(f"volume {volume_ids_str} scheme tetmesh")
+        cubit.cmd(
+            f"volume {volume_ids_str} sizing function type skeleton min_size "
+            f"{min_size} max_size {max_size} max_gradient {max_gradient} "
+            "min_num_layers_3d 1 min_num_layers_2d 1 min_num_layers_1d 1"
+        )
+        cubit.cmd(f"mesh volume {volume_ids_str}")
 
     def export_mesh(self, mesh_filename="magnet_mesh", export_dir=""):
         """Creates tetrahedral mesh of magnet volumes and exports H5M format

parastell/parastell.py

@@ -7,7 +7,6 @@
 import pystell.read_vmec as read_vmec
 
 from . import log
-
 from . import invessel_build as ivb
 from . import magnet_coils as mc
 from . import source_mesh as sm
@@ -233,6 +232,7 @@ def export_magnets(
         export_mesh=False,
         mesh_filename="magnet_mesh",
         export_dir="",
+        **kwargs,
     ):
         """Export magnet components.
 
@@ -246,13 +246,21 @@ def export_magnets(
                 '.h5m' extension (optional, defaults to 'magnet_mesh').
             export_dir (str): directory to which to export output files
                 (optional, defaults to empty string).
+
+        Optional attributes:
+            min_size (float): minimum size of magnet mesh elements (defaults to
+                20.0).
+            max_size (float): maximum size of magnet mesh elements (defaults to
+                50.0).
+            max_gradient (float): maximum transition in magnet mesh element
+                size (defaults to 1.5).
         """
         self.magnet_set.export_step(
             step_filename=step_filename, export_dir=export_dir
         )
 
         if export_mesh:
-            self.magnet_set.mesh_magnets()
+            self.magnet_set.mesh_magnets(**kwargs)
             self.magnet_set.export_mesh(
                 mesh_filename=mesh_filename, export_dir=export_dir
             )
@@ -304,29 +312,6 @@ def export_source_mesh(self, filename="source_mesh", export_dir=""):
         """
         self.source_mesh.export_mesh(filename=filename, export_dir=export_dir)
 
-    def _import_ivb_step(self):
-        """Imports STEP files from in-vessel build into Coreform Cubit.
-        (Internal function not intended to be called externally)
-        """
-        for (
-            name,
-            data,
-        ) in self.invessel_build.radial_build.radial_build.items():
-            vol_id = cubit_io.import_step_cubit(
-                name, self.invessel_build.export_dir
-            )
-            data["vol_id"] = vol_id
-
-    def _import_magnets_step(self):
-        """Import STEP file for magnet set into Coreform Cubit.
-        (Internal function not intended to be called externally)
-        """
-        last_vol_id = cubit_io.import_step_cubit(
-            self.magnet_set.step_filename, self.magnet_set.export_dir
-        )
-
-        self.magnet_set.volume_ids = range(1, last_vol_id + 1)
-
     def _tag_materials_legacy(self):
         """Applies material tags to corresponding CAD volumes for legacy DAGMC
         neutronics model export.
@@ -382,16 +367,17 @@ def build_cubit_model(self, skip_imprint=False, legacy_faceting=True):
             "Building DAGMC neutronics model via Coreform Cubit..."
         )
 
+        # Ensure fresh Cubit instance
         if cubit_io.initialized:
             cubit.cmd("new")
         else:
             cubit_io.init_cubit()
 
         if self.invessel_build:
-            self._import_ivb_step()
+            self.invessel_build.import_step_cubit()
 
         if self.magnet_set:
-            self._import_magnets_step()
+            self.magnet_set.import_step_cubit()
 
         if skip_imprint:
             self.invessel_build.merge_layer_surfaces()

tests/test_magnet_coils.py

@@ -67,6 +67,8 @@ def test_magnet_construction(coil_set):
 
 def test_magnet_exports(coil_set):
 
+    volume_ids_exp = list(range(1, 2))
+
     remove_files()
 
     coil_set.populate_magnet_coils()
@@ -75,6 +77,8 @@ def test_magnet_exports(coil_set):
     assert Path("magnet_set.step").exists()
 
     coil_set.mesh_magnets()
+    assert coil_set.volume_ids == volume_ids_exp
+
     coil_set.export_mesh()
     assert Path("magnet_mesh.h5m").exists()
 

tests/test_parastell.py

@@ -114,9 +114,27 @@ def test_parastell(stellarator):
 
     assert Path(filename_exp).with_suffix(".h5m").exists()
 
+    chamber_volume_id_exp = 1
+    component_volume_id_exp = 2
+    magnet_volume_ids_exp = list(range(3, 4))
     filename_exp = "dagmc"
 
     stellarator.build_cubit_model()
+
+    assert (
+        stellarator.invessel_build.radial_build.radial_build["chamber"][
+            "vol_id"
+        ]
+        == chamber_volume_id_exp
+    )
+    assert (
+        stellarator.invessel_build.radial_build.radial_build[
+            component_name_exp
+        ]["vol_id"]
+        == component_volume_id_exp
+    )
+    assert stellarator.magnet_set.volume_ids == magnet_volume_ids_exp
+
     stellarator.export_dagmc(filename=filename_exp)
     stellarator.export_cub5(filename=filename_exp)