Moving Neutron Source

Problem Description

A time-dependent problem in which an isotropic neutron source moves through a three-dimensional box of air-like material along a prescribed piecewise-linear trajectory. This example demonstrates MC/DC’s moving-source capability.

Geometry and Materials

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

A single homogeneous material (air analogue) fills the entire domain:

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

Region

\(\Sigma_c\)

\(\Sigma_s\)

Air

0.002

0.008

Neutron speed: \(v = 2 \times 10^{5}\) cm/s.

Physical Assumptions

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

  • Isotropic scattering.

  • Time-dependent transport with a moving point-like source.

  • No fission.

Numerical Setup

The source moves in three consecutive phases:

Phase

Velocity (cm/s)

Duration (s)

1

\((1, 0, 0)\)

7

2

\((-0.5,\; 2,\; 0)\)

2

3

\((0,\; -3,\; 0)\)

1

Spatial mesh (tally)

\(201 \times 201\) in the \((x,y)\)-plane

Time mesh (tally)

46 equally spaced bins over \(t \in [0,10]\) s

Tally score

Scalar flux

Source particles

\(10^{5}\)

Batches

2

Quantities of Interest

  • Time-resolved 2-D flux distribution \(\phi(x,y,t)\).

  • Animated GIF of the neutron cloud following the moving source.

Reference Solution

No analytical reference. The solution is validated qualitatively by confirming that the flux maximum tracks the prescribed source trajectory.

Step-by-Step Walkthrough

1. Materials (lines 1–15)

 1import numpy as np
 2
 3import mcdc
 4
 5# ======================================================================================
 6# Set model
 7# ======================================================================================
 8
 9# Set materials
10air = mcdc.MaterialMG(
11    capture=np.array([0.002]),
12    scatter=np.array([[0.008]]),
13    speed=np.array([200000.0]),
14)
15

A single air-like material with speed defined for time-dependent transport.

2. Geometry (lines 17–27)

17min_x = mcdc.Surface.PlaneX(x=-5.0, boundary_condition="vacuum")
18max_x = mcdc.Surface.PlaneX(x=5.0, boundary_condition="vacuum")
19min_y = mcdc.Surface.PlaneY(y=-5.0, boundary_condition="vacuum")
20max_y = mcdc.Surface.PlaneY(y=5.0, boundary_condition="vacuum")
21min_z = mcdc.Surface.PlaneZ(z=-10.0, boundary_condition="vacuum")
22max_z = mcdc.Surface.PlaneZ(z=10.0, boundary_condition="vacuum")
23
24# Make cells
25mcdc.Cell(region=+min_x & -max_x & +min_y & -max_y & +min_z & -max_z, fill=air)
26
27# ======================================================================================

A simple box with vacuum boundaries. One cell fills the entire domain.

3. Moving Source (lines 33–49)

33    y=[-0.5, 0.5],
34    z=[-0.5, 0.5],
35    direction=[1.0, 1.0, 0.0],
36    polar_cosine=[-1.0, -0.9],
37    energy_group=0,
38    time=[0.0, 10.0],
39)
40src.move(
41    velocities=[
42        [1.0, 0.0, 0.0],
43        [-0.5, 2.0, 0.0],
44        [0.0, -3.0, 0.0],
45    ],
46    durations=[7.0, 2.0, 1.0],
47)
48
49# ======================================================================================

The source is created with spatial and angular extent, then src.move() assigns a piecewise-linear trajectory: three velocity segments with their durations. The source physically translates through the domain over time.

4. Tallies, Settings, and Run (lines 55–65)

55    x=np.linspace(-5.0, 5.0, 201),
56    y=np.linspace(-5.0, 5.0, 201),
57)
58mcdc.Tally(mesh=mesh, scores=["flux"], time=np.linspace(0, 10, 46))
59
60# Settings
61mcdc.settings.N_particle = 100000
62mcdc.settings.N_batch = 2
63
64# Run
65mcdc.run()

A structured \(201 \times 201\) mesh tally with 46 time bins captures the evolving 2-D flux. The companion process-output.py script generates an animated GIF.

What to try:

  • Change the velocity vectors to create a circular or zigzag path.

  • Add more time-resolution bins for smoother animation.

  • Compare with moving_pellet where the geometry moves instead of the source.

Full Input

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

The complete input used for this example is embedded below:

 1import numpy as np
 2
 3import mcdc
 4
 5# ======================================================================================
 6# Set model
 7# ======================================================================================
 8
 9# Set materials
10air = mcdc.MaterialMG(
11    capture=np.array([0.002]),
12    scatter=np.array([[0.008]]),
13    speed=np.array([200000.0]),
14)
15
16# Set container cell surfaces
17min_x = mcdc.Surface.PlaneX(x=-5.0, boundary_condition="vacuum")
18max_x = mcdc.Surface.PlaneX(x=5.0, boundary_condition="vacuum")
19min_y = mcdc.Surface.PlaneY(y=-5.0, boundary_condition="vacuum")
20max_y = mcdc.Surface.PlaneY(y=5.0, boundary_condition="vacuum")
21min_z = mcdc.Surface.PlaneZ(z=-10.0, boundary_condition="vacuum")
22max_z = mcdc.Surface.PlaneZ(z=10.0, boundary_condition="vacuum")
23
24# Make cells
25mcdc.Cell(region=+min_x & -max_x & +min_y & -max_y & +min_z & -max_z, fill=air)
26
27# ======================================================================================
28# Set source
29# ======================================================================================
30
31src = mcdc.Source(
32    x=[-4.0, -3.0],
33    y=[-0.5, 0.5],
34    z=[-0.5, 0.5],
35    direction=[1.0, 1.0, 0.0],
36    polar_cosine=[-1.0, -0.9],
37    energy_group=0,
38    time=[0.0, 10.0],
39)
40src.move(
41    velocities=[
42        [1.0, 0.0, 0.0],
43        [-0.5, 2.0, 0.0],
44        [0.0, -3.0, 0.0],
45    ],
46    durations=[7.0, 2.0, 1.0],
47)
48
49# ======================================================================================
50# Set tallies, settings, and run MC/DC
51# ======================================================================================
52
53# Tallies
54mesh = mcdc.MeshStructured(
55    x=np.linspace(-5.0, 5.0, 201),
56    y=np.linspace(-5.0, 5.0, 201),
57)
58mcdc.Tally(mesh=mesh, scores=["flux"], time=np.linspace(0, 10, 46))
59
60# Settings
61mcdc.settings.N_particle = 100000
62mcdc.settings.N_batch = 2
63
64# Run
65mcdc.run()

How to Run

From the repository root run:

python examples/moving_source/input.py

Expected Output

An HDF5 mesh tally with time-resolved 2-D flux slices and a GIF animation produced by the example’s post-processing script.