Skip to content

Commit

Permalink
Merge pull request CEMeNT-PSAAP#120 from ilhamv/main
Browse files Browse the repository at this point in the history
update c5g7 examples
  • Loading branch information
ilhamv authored Sep 12, 2023
2 parents 60e0c22 + f4c79df commit e6bed85
Show file tree
Hide file tree
Showing 3 changed files with 16 additions and 34 deletions.
6 changes: 3 additions & 3 deletions examples/c5g7/3d/TDX/input.py
Original file line number Diff line number Diff line change
Expand Up @@ -283,7 +283,7 @@ def set_mat(mat):
# At highest energy

energy = np.zeros(7)
energy[-1] = 1.0
energy[0] = 1.0

source = mcdc.source(
point=[pitch * 17 / 2, -pitch * 17 / 2, 0.0], time=[0.0, 15.0], energy=energy
Expand All @@ -295,10 +295,10 @@ def set_mat(mat):

# Tally
t_grid = np.linspace(0.0, 20.0, 201)
mcdc.tally(scores=["fission"], t=t_grid, g=g_grid)
mcdc.tally(scores=["fission"], t=t_grid)

# Setting
mcdc.setting(N_particle=1e5, active_bank_buff=10000)
mcdc.setting(N_particle=1e4, active_bank_buff=1000)

# Run
mcdc.run()
44 changes: 13 additions & 31 deletions examples/c5g7/3d/TDX/process.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,53 +4,35 @@
from matplotlib import cm
from matplotlib import colors


# =============================================================================
# Plot results
# =============================================================================

# Results
# Get results
with h5py.File("output.h5", "r") as f:
fis_avg = f["tally/fission/mean"][:]
fis_sd = f["tally/fission/sdev"][:]
t = f["tally/grid/t"][:]
x = f["tally/grid/x"][:]
z = f["tally/grid/z"][:]
dx = x[1] - x[0]
dz = z[1] - z[0]
dV = dx * dx * dz

t_mid = 0.5 * (t[:-1] + t[1:])
dt = t[1:] - t[:-1]

# Normalize
norm = np.sum(fis_avg[0]) / dt[0]
fis_avg /= norm
fis_sd /= norm
fis_avg /= dt
fis_sd /= dt

# Reference
# Get reference
with h5py.File("reference.h5", "r") as f:
fis_avg_ref = f["tally/fission/mean"][:]
fis_sd_ref = f["tally/fission/sdev"][:]
fis_ref = f["tally/fission/mean"][:]
fis_ref_sd = f["tally/fission/sdev"][:]

# Normalize
norm = np.sum(fis_avg_ref[:, 0]) / dt[0]
fis_avg_ref /= norm
fis_sd_ref /= norm
fis__ref = np.zeros_like(fis_avg_ref[0])
fis__sd_ref = np.zeros_like(fis_avg_ref[0])
for i in range(7):
fis__ref += fis_avg_ref[i]
fis__sd_ref += np.square(fis_sd_ref[i])
fis__sd_ref = np.sqrt(fis__sd_ref)
fis__ref /= dt
fis__sd_ref /= dt

norm = fis_ref[0] / dt[0]
fis_avg /= norm * dt
fis_sd /= norm * dt
fis_ref /= norm * dt
fis_ref_sd /= norm * dt

plt.plot(t_mid, fis__ref, "-m", fillstyle="none", label="Ref. (MC)")
# Plot
plt.plot(t_mid, fis_ref, "-m", fillstyle="none", label="Ref. (MC)")
plt.fill_between(
t_mid, fis__ref - fis__sd_ref, fis__ref + fis__sd_ref, alpha=0.2, color="m"
t_mid, fis_ref - fis_ref_sd, fis_ref + fis_ref_sd, alpha=0.2, color="m"
)
plt.plot(t_mid, fis_avg, "ok", fillstyle="none", label="MC")
plt.fill_between(t_mid, fis_avg - fis_sd, fis_avg + fis_sd, alpha=0.2, color="k")
Expand Down
Binary file modified examples/c5g7/3d/TDX/reference.h5
Binary file not shown.

0 comments on commit e6bed85

Please sign in to comment.