Merge pull request #14 from davidepettinari/support

Add HDPE source supports

Author: RemDelaporteMathurin

analysis/neutron/helpers.py

@@ -25,35 +25,4 @@ def get_exp_cllif_density(temp, LiCl_frac=0.695):
 
     rho = a + b * C + c * temp + d * C**2 + e * C**3 + f * temp * C**2 + g * C * temp**2
 
-    return rho
-
-
-def calculate_cylinder_volume(radius, height):
-    volume = math.pi * radius**2 * height
-    return volume
-
-
-def translate_surface(
-    surface: (
-        openmc.XPlane | openmc.YPlane | openmc.ZPlane | openmc.Plane | openmc.Sphere
-    ),
-    dx: float,
-    dy: float,
-    dz: float,
-) -> openmc.Surface:
-    """Translate an OpenMC surface by dx, dy, dz."""
-    if isinstance(surface, openmc.XPlane):
-        surface.x0 += dx
-    elif isinstance(surface, openmc.YPlane):
-        surface.y0 += dy
-    elif isinstance(surface, openmc.ZPlane):
-        surface.z0 += dz
-    elif isinstance(surface, openmc.Plane):
-        surface.d += surface.a * dx + surface.b * dy + surface.c * dz
-    elif isinstance(surface, openmc.Sphere):
-        surface.x0 += dx
-        surface.y0 += dy
-        surface.z0 += dz
-    else:
-        raise TypeError(f"Unsupported surface type: {type(surface)}")
-    return surface
+    return rho

analysis/neutron/openmc_model.py

@@ -1,6 +1,6 @@
 import openmc
-import vault
 from libra_toolbox.neutronics.neutron_source import A325_generator_diamond
+from libra_toolbox.neutronics import vault
 import helpers
 
 
@@ -167,15 +167,15 @@ def baby_geometry(x_c: float, y_c: float, z_c: float):
     sphere = openmc.Sphere(x0=x_c, y0=y_c, z0=z_c, r=50.00)  # before r=50.00
 
     ######## Lead bricks positioned under the source #################
-    positions = [
+    lead_positions = [
         (x_c - 13.50, y_c, z_c - z_tab),
-        (x_c - 4.50, y_c, z_c - z_tab),
+        (x_c - 2.96, y_c, z_c - z_tab),
         (x_c + 36.50, y_c, z_c - z_tab),
-        (x_c + 27.50, y_c, z_c - z_tab),
+        (x_c + 25.96, y_c, z_c - z_tab),
     ]
 
     lead_blocks = []
-    for position in positions:
+    for position in lead_positions:
         lead_block_region = openmc.model.RectangularParallelepiped(
             position[0] - lead_width / 2,
             position[0] + lead_width / 2,
@@ -186,6 +186,25 @@ def baby_geometry(x_c: float, y_c: float, z_c: float):
         )
         lead_blocks.append(lead_block_region)
 
+    src_supp_length = 40.00
+    src_supp_height = 20.00
+    src_supp_width = 2.54
+    src_supp_position = [
+        (x_c - 13.50 + lead_width / 2, y_c, z_c - z_tab),
+        (x_c + 25.96 + lead_width / 2, y_c, z_c - z_tab),
+    ]
+    src_supports = []
+    for position in src_supp_position:
+        source_support = openmc.model.RectangularParallelepiped(
+            position[0],
+            position[0] + src_supp_width,
+            position[1] - src_supp_length / 2,
+            position[1] + src_supp_length / 2,
+            position[2],
+            position[2] + src_supp_height,
+        )
+        src_supports.append(source_support)
+        
     # regions
     source_wall_region = -ext_cyl_source & +source_region
     source_region = -source_region
@@ -209,6 +228,8 @@ def baby_geometry(x_c: float, y_c: float, z_c: float):
     lead_block_2_region = -lead_blocks[1]
     lead_block_3_region = -lead_blocks[2]
     lead_block_4_region = -lead_blocks[3]
+    src_supp_1_region = -src_supports[0] & ~source_wall_region & ~source_region
+    src_supp_2_region = -src_supports[1] & ~source_wall_region & ~source_region 
     he_region = (
         +z_plane_5
         & -z_plane_12
@@ -228,6 +249,8 @@ def baby_geometry(x_c: float, y_c: float, z_c: float):
         & ~lead_block_2_region
         & ~lead_block_3_region
         & ~lead_block_4_region
+        & ~src_supp_1_region
+        & ~src_supp_2_region
     )
     sphere_region = (
         -sphere
@@ -248,6 +271,8 @@ def baby_geometry(x_c: float, y_c: float, z_c: float):
         & ~lead_block_2_region
         & ~lead_block_3_region
         & ~lead_block_4_region
+        & ~src_supp_1_region
+        & ~src_supp_2_region
     )
 
     # cells
@@ -287,6 +312,10 @@ def baby_geometry(x_c: float, y_c: float, z_c: float):
     lead_block_3_cell.fill = lead
     lead_block_4_cell = openmc.Cell(region=lead_block_4_region)
     lead_block_4_cell.fill = lead
+    src_supp_1_cell = openmc.Cell(region=src_supp_1_region)
+    src_supp_1_cell.fill = HDPE
+    src_supp_2_cell = openmc.Cell(region=src_supp_2_region)
+    src_supp_2_cell.fill = HDPE
 
     cells = [
         source_wall_cell_1,
@@ -307,6 +336,8 @@ def baby_geometry(x_c: float, y_c: float, z_c: float):
         lead_block_2_cell,
         lead_block_3_cell,
         lead_block_4_cell,
+        src_supp_1_cell,
+        src_supp_2_cell,
     ]
 
     return sphere, cllif_cell, cells
@@ -330,6 +361,7 @@ def baby_model():
         air,
         epoxy,
         he,
+        HDPE,
     ]
 
     # BABY coordinates
@@ -486,6 +518,12 @@ def baby_model():
 lead.add_nuclide("Pb207", 0.221, "ao")
 lead.add_nuclide("Pb208", 0.524, "ao")
 
+# High Density Polyethylene
+# Reference:  PNNL Report 15870 (Rev. 1)
+HDPE = openmc.Material(name="HDPE")
+HDPE.set_density("g/cm3", 0.95)
+HDPE.add_element("H", 0.143724, "wo")
+HDPE.add_element("C", 0.856276, "wo")
 
 if __name__ == "__main__":
     model = baby_model()
@@ -519,4 +557,4 @@ def baby_model():
     with open(processed_data_file, "w") as f:
         json.dump(existing_data, f, indent=4)
 
-    print(f"Processed data stored in {processed_data_file}")
+    print(f"Processed data stored in {processed_data_file}")