💨 ♻️ Create Vacuum Vessel class

process/blanket_library.py

@@ -57,25 +57,17 @@ def component_volumes(self):
         blanket_library.dz_blkt_half = self.component_half_height(icomponent=0)
         # Shield
         blanket_library.dz_shld_half = self.component_half_height(icomponent=1)
-        # Vacuum Vessel
-        blanket_library.dz_vv_half = self.component_half_height(icomponent=2)
 
         # D-shaped blanket and shield
         if physics_variables.itart == 1 or fwbs_variables.i_fw_blkt_vv_shape == 1:
-            for icomponent in range(3):
+            for icomponent in range(2):
                 self.dshaped_component(icomponent)
 
         # Elliptical blanket and shield
         else:
-            for icomponent in range(3):
+            for icomponent in range(2):
                 self.elliptical_component(icomponent)
 
-            # This will fail the hts_REBCO and 2D_scan regression tests,
-            # the number of VMCON iterations (nviter) is different.
-            # Seems to be because in the blanket calculations (icomponent=0):
-            # r2 = 1.3836567143743970 rather than old value of r2 = 1.3836567143743972,
-            # r3 = 3.7009701431231936 rather than r3 = 3.7009701431231923.
-
         # Apply coverage factors to volumes and surface areas
         self.apply_coverage_factors()
 
@@ -101,15 +93,8 @@ def component_half_height(self, icomponent: int):
                 + build_variables.dz_xpoint_divertor
                 + divertor_variables.dz_divertor
             )
-        # Vacuum vessel
-        elif icomponent == 2:
-            hbot = (
-                build_variables.z_tf_inside_half
-                - build_variables.dz_shld_vv_gap
-                - build_variables.dz_vv_lower
-            )
         else:
-            raise ProcessValueError(f"{icomponent=} is invalid, it must be either 0,1,2")
+            raise ProcessValueError(f"{icomponent=} is invalid, it must be either 0,1")
 
         # Calculate component internal upper half-height (m)
         # If a double null machine then symmetric
@@ -126,11 +111,6 @@ def component_half_height(self, icomponent: int):
             # Shield
             if icomponent == 1:
                 htop = htop + build_variables.dz_blkt_upper
-            # Vacuum Vessel
-            if icomponent == 2:
-                htop = (
-                    htop + build_variables.dz_blkt_upper + build_variables.dz_shld_upper
-                )
 
         # Average of top and bottom (m)
         return 0.5 * (htop + hbot)
@@ -164,9 +144,6 @@ def dshaped_component(self, icomponent: int):
         # Sheild
         if icomponent == 1:
             r2 = build_variables.dr_blkt_inboard + r2 + build_variables.dr_blkt_outboard
-        # Vaccum Vessel
-        if icomponent == 2:
-            r2 = build_variables.rsldo - r1
 
         # Calculate surface area, assuming 100% coverage
         if icomponent == 0:
@@ -209,19 +186,6 @@ def dshaped_component(self, icomponent: int):
                 build_variables.dr_shld_outboard,
                 build_variables.dz_shld_upper,
             )
-        elif icomponent == 2:
-            (
-                blanket_library.vol_vv_inboard,
-                blanket_library.vol_vv_outboard,
-                fwbs_variables.vol_vv,
-            ) = dshellvol(
-                r1,
-                r2,
-                blanket_library.dz_vv_half,
-                build_variables.dr_vv_inboard,
-                build_variables.dr_vv_outboard,
-                (build_variables.dz_vv_upper + build_variables.dz_vv_lower) / 2,
-            )
 
     def elliptical_component(self, icomponent: int):
         """Calculate component surface area and volume using elliptical scheme
@@ -299,20 +263,6 @@ def elliptical_component(self, icomponent: int):
                 build_variables.dr_shld_outboard,
                 build_variables.dz_shld_upper,
             )
-        if icomponent == 2:
-            (
-                blanket_library.vol_vv_inboard,
-                blanket_library.vol_vv_outboard,
-                fwbs_variables.vol_vv,
-            ) = eshellvol(
-                r1,
-                r2,
-                r3,
-                blanket_library.dz_vv_half,
-                build_variables.dr_vv_inboard,
-                build_variables.dr_vv_outboard,
-                (build_variables.dz_vv_upper + build_variables.dz_vv_lower) / 2,
-            )
 
     def apply_coverage_factors(self):
         """Apply coverage factors to volumes
@@ -383,11 +333,6 @@ def apply_coverage_factors(self):
             blanket_library.vol_shld_inboard + blanket_library.vol_shld_outboard
         )
 
-        # Apply vacuum vessel coverage factor
-        # moved from dshaped_* and elliptical_* to keep coverage factor
-        # changes in the same location.
-        fwbs_variables.vol_vv = fwbs_variables.fvoldw * fwbs_variables.vol_vv
-
     def primary_coolant_properties(self, output: bool):
         """Calculates the fluid properties of the Primary Coolant in the FW and BZ.
         Uses middle value of input and output temperatures of coolant.

process/caller.py

@@ -296,6 +296,8 @@ def _call_models_once(self, xc: np.ndarray) -> None:
             # DCLL model
             self.models.dcll.run(output=False)
 
+        self.models.vacuum_vessel.run()
+
         self.models.divertor.run(output=False)
 
         self.models.cryostat.run()

process/hcpb.py

@@ -474,11 +474,6 @@ def nuclear_heating_magnets(self, output: bool):
         if build_variables.dr_vv_outboard > d_vv_all:
             d_vv_all = build_variables.dr_vv_outboard
 
-        if d_vv_all > 1.0e-6:
-            ccfe_hcpb_module.vv_density = fwbs_variables.m_vv / fwbs_variables.vol_vv
-        else:
-            ccfe_hcpb_module.vv_density = 0.0
-
         # Calculation of average blanket/shield thickness [m]
         if physics_variables.itart == 1:
             # There is no inner blanket for TART design [m]

process/main.py

@@ -105,7 +105,7 @@
 from process.structure import Structure
 from process.superconducting_tf_coil import SuperconductingTFCoil
 from process.tf_coil import TFCoil
-from process.vacuum import Vacuum
+from process.vacuum import Vacuum, VacuumVessel
 from process.water_use import WaterUse
 
 os.environ["PYTHON_PROCESS_ROOT"] = os.path.join(os.path.dirname(__file__))
@@ -664,6 +664,7 @@ def __init__(self):
         self.availability = Availability()
         self.buildings = Buildings()
         self.vacuum = Vacuum()
+        self.vacuum_vessel = VacuumVessel()
         self.water_use = WaterUse()
         self.pulse = Pulse()
         self.ife = IFE(availability=self.availability, costs=self.costs)

process/vacuum.py

@@ -5,8 +5,13 @@
 
 from process import constants
 from process import process_output as po
+from process.blanket_library import dshellvol, eshellvol
+from process.data_structure import blanket_library as blanket_library
 from process.data_structure import build_variables as buv
+from process.data_structure import ccfe_hcpb_module as ccfe_hcpb_module
 from process.data_structure import divertor_variables as dv
+from process.data_structure import fwbs_variables as fwbs_variables
+from process.data_structure import physics_variables as physics_variables
 from process.data_structure import physics_variables as pv
 from process.data_structure import tfcoil_variables as tfv
 from process.data_structure import times_variables as tv
@@ -684,3 +689,166 @@ def vacuum(
             )
 
         return pumpn, nduct, dlscalc, mvdsh, dimax
+
+
+class VacuumVessel:
+    """Class containing vacuum vessel routines"""
+
+    def __init__(self) -> None:
+        pass
+
+    def run(self) -> None:
+        blanket_library.dz_vv_half = self.calculate_vessel_half_height(
+            z_tf_inside_half=buv.z_tf_inside_half,
+            dz_shld_vv_gap=buv.dz_shld_vv_gap,
+            dz_vv_lower=buv.dz_vv_lower,
+            n_divertors=dv.n_divertors,
+            dz_blkt_upper=buv.dz_blkt_upper,
+            dz_shld_upper=buv.dz_shld_upper,
+            z_plasma_xpoint_upper=buv.z_plasma_xpoint_upper,
+            dr_fw_plasma_gap_inboard=buv.dr_fw_plasma_gap_inboard,
+            dr_fw_plasma_gap_outboard=buv.dr_fw_plasma_gap_outboard,
+            dr_fw_inboard=buv.dr_fw_inboard,
+            dr_fw_outboard=buv.dr_fw_outboard,
+        )
+        # D-shaped blanket and shield
+        if physics_variables.itart == 1 or fwbs_variables.i_fw_blkt_vv_shape == 1:
+            (
+                blanket_library.vol_vv_inboard,
+                blanket_library.vol_vv_outboard,
+                fwbs_variables.vol_vv,
+            ) = self.calculate_dshaped_vessel_volumes(
+                rsldi=buv.rsldi,
+                rsldo=buv.rsldo,
+                dz_vv_half=blanket_library.dz_vv_half,
+                dr_vv_inboard=buv.dr_vv_inboard,
+                dr_vv_outboard=buv.dr_vv_outboard,
+                dz_vv_upper=buv.dz_vv_upper,
+                dz_vv_lower=buv.dz_vv_lower,
+            )
+        else:
+            (
+                blanket_library.vol_vv_inboard,
+                blanket_library.vol_vv_outboard,
+                fwbs_variables.vol_vv,
+            ) = self.calculate_elliptical_vessel_volumes(
+                rmajor=pv.rmajor,
+                rminor=pv.rminor,
+                triang=pv.triang,
+                rsldi=buv.rsldi,
+                rsldo=buv.rsldo,
+                dz_vv_half=blanket_library.dz_vv_half,
+                dr_vv_inboard=buv.dr_vv_inboard,
+                dr_vv_outboard=buv.dr_vv_outboard,
+                dz_vv_upper=buv.dz_vv_upper,
+                dz_vv_lower=buv.dz_vv_lower,
+            )
+
+        # Apply vacuum vessel coverage factor
+        # moved from dshaped_* and elliptical_* to keep coverage factor
+        # changes in the same location.
+        fwbs_variables.vol_vv = fwbs_variables.fvoldw * fwbs_variables.vol_vv
+
+        ccfe_hcpb_module.vv_density = fwbs_variables.m_vv / fwbs_variables.vol_vv
+
+    def calculate_vessel_half_height(
+        self,
+        z_tf_inside_half: float,
+        dz_shld_vv_gap: float,
+        dz_vv_lower: float,
+        n_divertors: int,
+        dz_blkt_upper: float,
+        dz_shld_upper: float,
+        z_plasma_xpoint_upper: float,
+        dr_fw_plasma_gap_inboard: float,
+        dr_fw_plasma_gap_outboard: float,
+        dr_fw_inboard: float,
+        dr_fw_outboard: float,
+    ) -> float:
+        """Calculate vacuum vessel internal half-height (m)"""
+
+        z_bottom = z_tf_inside_half - dz_shld_vv_gap - dz_vv_lower
+
+        # Calculate component internal upper half-height (m)
+        # If a double null machine then symmetric
+        if n_divertors == 2:
+            z_top = z_bottom
+        else:
+            z_top = z_plasma_xpoint_upper + 0.5 * (
+                dr_fw_plasma_gap_inboard
+                + dr_fw_plasma_gap_outboard
+                + dr_fw_inboard
+                + dr_fw_outboard
+            )
+
+            z_top = z_top + dz_blkt_upper + dz_shld_upper
+
+        # Average of top and bottom (m)
+        return 0.5 * (z_top + z_bottom)
+
+    def calculate_dshaped_vessel_volumes(
+        self,
+        rsldi: float,
+        rsldo: float,
+        dz_vv_half: float,
+        dr_vv_inboard: float,
+        dr_vv_outboard: float,
+        dz_vv_upper: float,
+        dz_vv_lower: float,
+    ) -> tuple[float, float, float]:
+        """Calculate volumes of D-shaped vacuum vessel segments"""
+
+        r_1 = rsldi
+        r_2 = rsldo - r_1
+
+        (
+            vol_vv_inboard,
+            vol_vv_outboard,
+            vol_vv,
+        ) = dshellvol(
+            rmajor=r_1,
+            rminor=r_2,
+            zminor=dz_vv_half,
+            drin=dr_vv_inboard,
+            drout=dr_vv_outboard,
+            dz=(dz_vv_upper + dz_vv_lower) / 2,
+        )
+
+        return vol_vv_inboard, vol_vv_outboard, vol_vv
+
+    def calculate_elliptical_vessel_volumes(
+        self,
+        rmajor: float,
+        rminor: float,
+        triang: float,
+        rsldi: float,
+        rsldo: float,
+        dz_vv_half: float,
+        dr_vv_inboard: float,
+        dr_vv_outboard: float,
+        dz_vv_upper: float,
+        dz_vv_lower: float,
+    ) -> tuple[float, float, float]:
+        # Major radius to centre of inboard and outboard ellipses (m)
+        # (coincident in radius with top of plasma)
+        r_1 = rmajor - rminor * triang
+
+        # Calculate distance between r1 and outer edge of inboard ...
+        # ... section (m)
+        r_2 = r_1 - rsldi
+        r_3 = rsldo - r_1
+
+        (
+            vol_vv_inboard,
+            vol_vv_outboard,
+            vol_vv,
+        ) = eshellvol(
+            r_1,
+            r_2,
+            r_3,
+            dz_vv_half,
+            dr_vv_inboard,
+            dr_vv_outboard,
+            (dz_vv_upper + dz_vv_lower) / 2,
+        )
+        return vol_vv_inboard, vol_vv_outboard, vol_vv