Packaged Fuel Array (Universe and Lattice)#

Problem Description#

A three-dimensional fixed-source problem featuring two identical fuel assemblies placed side-by-side inside a water-filled box. Each assembly uses a composite “shooting-star” fuel geometry built from the union of two orthogonal cylinders enclosed by a cladding sphere.

This example demonstrates MC/DC’s universe, translation, and rotation capabilities for constructive solid geometry (CSG) packaging.

Geometry and Materials#

The global domain is a rectangular box: \(x \in [-10,10]\), \(y \in [-5,5]\), \(z \in [-5,5]\) cm, with vacuum boundary conditions.

Each assembly is defined as a universe containing three cells:

  1. Fuel — union of a z-aligned and an x-aligned cylinder (radius 1 cm, half-length 5 cm).

  2. Cladding — spherical shell (radius 3 cm) surrounding the fuel.

  3. Water — region outside the cladding sphere.

The left assembly is translated to \((-5,0,0)\) cm; the right assembly is translated to \((+5,0,0)\) cm and rotated \(10°\) about the \(y\)-axis.

Cross-section data (mono-energetic, cm-1)#

Region

\(\Sigma_c\)

\(\Sigma_s\)

\(\Sigma_f\)

\(\nu\)

Fuel

0.45

—

0.55

2.5

Cladding

0.05

0.95

—

—

Water

0.02

0.08

—

—

Physical Assumptions#

  • Mono-energetic (one-speed) neutron transport.

  • Isotropic scattering.

  • Steady-state fixed-source calculation.

  • No delayed neutrons.

Numerical Setup#

Spatial mesh (tally)

\(201 \times 101\) in the \((x,z)\)-plane

Tally score

Fission rate

Source particles

\(10^{3}\) (demonstration)

Batches

2

Quantities of Interest#

  • Two-dimensional fission rate distribution in the \((x,z)\)-plane.

  • Relative standard deviation map for convergence assessment.

Reference Solution#

No analytical reference. The geometry can be verified using MC/DC’s built-in simulation.visualize_model() method to render the CSG model.

Step-by-Step Walkthrough#

1. Materials (lines 1–27)

 1import numpy as np
 2import mcdc
 3
 4simulation = mcdc.Simulation("Packaged fuel array")
 5
 6# ======================================================================================
 7# Materials
 8# ======================================================================================
 9
10fuel = mcdc.Material.multigroup(
11    capture=np.array([0.45]),
12    fission=np.array([0.55]),
13    nu_p=np.array([2.5]),
14)
15
16cover = mcdc.Material.multigroup(
17    capture=np.array([0.05]),
18    scatter=np.array([[0.95]]),
19)
20
21water = mcdc.Material.multigroup(
22    capture=np.array([0.02]),
23    scatter=np.array([[0.08]]),
24)
25
26# ======================================================================================
27# The assembly

Three mono-energetic materials: fissile fuel, a scattering cladding, and water moderator.

2. Assembly Geometry — Shooting-Star CSG (lines 29–54)

29
30# Surfaces
31cylinder_z = mcdc.Surface.CylinderZ(center=[0.0, 0.0], radius=1.0)
32cylinder_x = mcdc.Surface.CylinderX(center=[0.0, 0.0], radius=1.0)
33
34top_z = mcdc.Surface.PlaneZ(z=2.5)
35bot_z = mcdc.Surface.PlaneZ(z=-2.5)
36top_x = mcdc.Surface.PlaneX(x=2.5)
37bot_x = mcdc.Surface.PlaneX(x=-2.5)
38
39sphere = mcdc.Surface.Sphere(center=[0.0, 0.0, 0.0], radius=3.0)
40
41# Cells
42pellet_z = -cylinder_z & +bot_z & -top_z
43pellet_x = -cylinder_x & +bot_x & -top_x
44shooting_star = pellet_z | pellet_x
45fuel_shooting_star = mcdc.Cell(region=shooting_star, fill=fuel)
46cover_sphere = mcdc.Cell(region=-sphere & ~shooting_star, fill=cover)
47water_tank = mcdc.Cell(region=+sphere, fill=water)
48
49# ======================================================================================
50# Copy the assembly via universe cells
51# ======================================================================================
52
53# Set the universe
54assembly = mcdc.Universe(cells=[fuel_shooting_star, cover_sphere, water_tank])

The fuel region is the union of a z-cylinder and an x-cylinder (the “shooting star”). The cladding fills the sphere minus the fuel. Water fills outside the sphere. These three cells form a reusable universe.

3. Packaging with Universe, Translation, and Rotation (lines 56–80)

56# Set container cell surfaces
57min_x = mcdc.Surface.PlaneX(x=-10.0, boundary_condition="vacuum")
58mid_x = mcdc.Surface.PlaneX(x=0.0)
59max_x = mcdc.Surface.PlaneX(x=10.0, boundary_condition="vacuum")
60min_y = mcdc.Surface.PlaneY(y=-5.0, boundary_condition="vacuum")
61max_y = mcdc.Surface.PlaneY(y=5.0, boundary_condition="vacuum")
62min_z = mcdc.Surface.PlaneZ(z=-5.0, boundary_condition="vacuum")
63max_z = mcdc.Surface.PlaneZ(z=5.0, boundary_condition="vacuum")
64
65# Make copies via universe cells
66container_left = +min_y & -max_y & +min_z & -max_z & +min_x & -mid_x
67container_right = +min_y & -max_y & +min_z & -max_z & +mid_x & -max_x
68assembly_left = mcdc.Cell(region=container_left, fill=assembly, translation=[-5, 0, 0])
69assembly_right = mcdc.Cell(
70    region=container_right, fill=assembly, translation=[+5, 0, 0], rotation=[0, 10, 0]
71)
72
73# Set model
74simulation.set_model([assembly_left, assembly_right])
75
76# ======================================================================================
77# Set source
78# ======================================================================================
79
80source = mcdc.Source(x=[-0.1, 0.1], isotropic=True, energy=0)

The assembly universe is placed twice using mcdc.Cell(..., fill=assembly):

  • Left — translated to \((-5, 0, 0)\).

  • Right — translated to \((+5, 0, 0)\) and rotated 10° about \(y\).

simulation.set_model() tells MC/DC these are the top-level cells.

4. Source, Tallies, Settings, and Run (lines 82–105)

 82
 83# ======================================================================================
 84# Set tallies, settings, and run MC/DC
 85# ======================================================================================
 86
 87# Tallies
 88mesh = mcdc.MeshStructured(
 89    x=np.linspace(-10, 10, 201),
 90    z=np.linspace(-5, 5, 101),
 91)
 92tally = mcdc.Tally(mesh=mesh, scores=["fission"])
 93simulation.set_tallies([tally])
 94
 95# Settings
 96simulation.settings.N_particle = 1000
 97simulation.settings.N_batch = 2
 98simulation.settings.active_bank_buffer = 1000
 99
100# Run (or visualize)
101visualize = False
102if not visualize:
103    simulation.run()
104else:
105    colors = {

A point-like source near the centre, a structured mesh tally for the \((x,z)\)-plane fission rate, and 1 000 particles in 2 batches. The active_bank_buffer accommodates fission-born particles.

5. Optional Visualization (lines 107–end)

107        cover: "gray",
108        water: "blue",
109    }
110    simulation.visualize_model(
111        vis_plane="xz",
112        y=0.0,
113        x=[-11.0, 11.0],
114        z=[-6, 6],
115        pixels=(400, 400),
116        colors=colors,
117        time=[0.0],
118        save_as=None,
119    )

Set visualize = True to render the CSG geometry with simulation.visualize_model() instead of running the transport.

What to try:

  • Change the rotation angle and observe the effect on the fission map.

  • Add a third assembly copy with a different translation.

  • Use mcdc.Lattice instead of manual universe placement.

Full Input#

Click here to view the input file: examples/fuel_array_packaged/input.py.

The complete input used for this example is embedded below:

  1import numpy as np
  2import mcdc
  3
  4simulation = mcdc.Simulation("Packaged fuel array")
  5
  6# ======================================================================================
  7# Materials
  8# ======================================================================================
  9
 10fuel = mcdc.Material.multigroup(
 11    capture=np.array([0.45]),
 12    fission=np.array([0.55]),
 13    nu_p=np.array([2.5]),
 14)
 15
 16cover = mcdc.Material.multigroup(
 17    capture=np.array([0.05]),
 18    scatter=np.array([[0.95]]),
 19)
 20
 21water = mcdc.Material.multigroup(
 22    capture=np.array([0.02]),
 23    scatter=np.array([[0.08]]),
 24)
 25
 26# ======================================================================================
 27# The assembly
 28# ======================================================================================
 29
 30# Surfaces
 31cylinder_z = mcdc.Surface.CylinderZ(center=[0.0, 0.0], radius=1.0)
 32cylinder_x = mcdc.Surface.CylinderX(center=[0.0, 0.0], radius=1.0)
 33
 34top_z = mcdc.Surface.PlaneZ(z=2.5)
 35bot_z = mcdc.Surface.PlaneZ(z=-2.5)
 36top_x = mcdc.Surface.PlaneX(x=2.5)
 37bot_x = mcdc.Surface.PlaneX(x=-2.5)
 38
 39sphere = mcdc.Surface.Sphere(center=[0.0, 0.0, 0.0], radius=3.0)
 40
 41# Cells
 42pellet_z = -cylinder_z & +bot_z & -top_z
 43pellet_x = -cylinder_x & +bot_x & -top_x
 44shooting_star = pellet_z | pellet_x
 45fuel_shooting_star = mcdc.Cell(region=shooting_star, fill=fuel)
 46cover_sphere = mcdc.Cell(region=-sphere & ~shooting_star, fill=cover)
 47water_tank = mcdc.Cell(region=+sphere, fill=water)
 48
 49# ======================================================================================
 50# Copy the assembly via universe cells
 51# ======================================================================================
 52
 53# Set the universe
 54assembly = mcdc.Universe(cells=[fuel_shooting_star, cover_sphere, water_tank])
 55
 56# Set container cell surfaces
 57min_x = mcdc.Surface.PlaneX(x=-10.0, boundary_condition="vacuum")
 58mid_x = mcdc.Surface.PlaneX(x=0.0)
 59max_x = mcdc.Surface.PlaneX(x=10.0, boundary_condition="vacuum")
 60min_y = mcdc.Surface.PlaneY(y=-5.0, boundary_condition="vacuum")
 61max_y = mcdc.Surface.PlaneY(y=5.0, boundary_condition="vacuum")
 62min_z = mcdc.Surface.PlaneZ(z=-5.0, boundary_condition="vacuum")
 63max_z = mcdc.Surface.PlaneZ(z=5.0, boundary_condition="vacuum")
 64
 65# Make copies via universe cells
 66container_left = +min_y & -max_y & +min_z & -max_z & +min_x & -mid_x
 67container_right = +min_y & -max_y & +min_z & -max_z & +mid_x & -max_x
 68assembly_left = mcdc.Cell(region=container_left, fill=assembly, translation=[-5, 0, 0])
 69assembly_right = mcdc.Cell(
 70    region=container_right, fill=assembly, translation=[+5, 0, 0], rotation=[0, 10, 0]
 71)
 72
 73# Set model
 74simulation.set_model([assembly_left, assembly_right])
 75
 76# ======================================================================================
 77# Set source
 78# ======================================================================================
 79
 80source = mcdc.Source(x=[-0.1, 0.1], isotropic=True, energy=0)
 81simulation.set_sources([source])
 82
 83# ======================================================================================
 84# Set tallies, settings, and run MC/DC
 85# ======================================================================================
 86
 87# Tallies
 88mesh = mcdc.MeshStructured(
 89    x=np.linspace(-10, 10, 201),
 90    z=np.linspace(-5, 5, 101),
 91)
 92tally = mcdc.Tally(mesh=mesh, scores=["fission"])
 93simulation.set_tallies([tally])
 94
 95# Settings
 96simulation.settings.N_particle = 1000
 97simulation.settings.N_batch = 2
 98simulation.settings.active_bank_buffer = 1000
 99
100# Run (or visualize)
101visualize = False
102if not visualize:
103    simulation.run()
104else:
105    colors = {
106        fuel: "red",
107        cover: "gray",
108        water: "blue",
109    }
110    simulation.visualize_model(
111        vis_plane="xz",
112        y=0.0,
113        x=[-11.0, 11.0],
114        z=[-6, 6],
115        pixels=(400, 400),
116        colors=colors,
117        time=[0.0],
118        save_as=None,
119    )

How to Run#

From inside examples/fuel_array_packaged run:

python input.py

Expected Output#

An HDF5 mesh tally and optional visualization images produced by the simulation.visualize_model() method when run with visualization enabled.