Hello,
I’m trying to run some Cartesian grid style problems in OpenMC. The problems themselves run quite quickly once the simulation actually begins, however the source initialization is prohibitively slow for the application I had envisioned.
I was wondering if anyone had some insight into what might be causing this problem and what potential workarounds might exist? For instance, the following grid style problem takes about 16 minutes to run with 32 OpenMP threads, despite only being a 5x5 grid and having a relatively small number of total particles run:
import numpy as np
import openmc
#Group structure:
groups = openmc.mgxs.EnergyGroups([ 0.00000000E+00, 1.00000000E+00, 2.00000000E+00, 3.00000000E+00, 4.00000000E+00, 5.00000000E+00
, 6.00000000E+00, 7.00000000E+00])
#Data for Material 0:
mat0_xsdat = openmc.XSdata('mat_0', groups)
mat0_xsdat.order = 0
mat0_xsdat.set_total([ 1.59206000E-01, 4.12970000E-01, 5.90310000E-01, 5.84350000E-01, 7.18000000E-01
, 1.25445000E+00, 2.65038000E+00], temperature=294.)
mat0_xsdat.set_absorption([ 6.01026916E-04, 1.61039000E-05, 3.37260000E-04, 1.94060000E-03, 5.74156300E-03
, 1.50013000E-02, 3.72400000E-02], temperature=294.)
mat0_xsdat.set_fission([ 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00
, 0.00000000E+00, 0.00000000E+00], temperature=294.)
mat0_xsdat.set_nu_fission([ 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00
, 0.00000000E+00, 0.00000000E+00], temperature=294.)
mat0_xsdat.set_chi([ 1.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00
, 0.00000000E+00, 0.00000000E+00], temperature=294.)
scatter_matrix = np.array(\
[[[ 4.44777000E-02, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00],
[ 1.13400000E-01, 2.82334000E-01, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00],
[ 7.23470000E-04, 1.29940000E-01, 3.45256000E-01, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00],
[ 3.74990000E-06, 6.23400000E-04, 2.24570000E-01, 9.10284000E-02, 7.14370000E-05, 0.00000000E+00, 0.00000000E+00],
[ 5.31840000E-08, 4.80020000E-05, 1.69990000E-02, 4.15510000E-01, 1.39138000E-01, 2.21570000E-03, 0.00000000E+00],
[ 0.00000000E+00, 7.44860000E-06, 2.64430000E-03, 6.37320000E-02, 5.11820000E-01, 6.99913000E-01, 1.32440000E-01],
[ 0.00000000E+00, 1.04550000E-06, 5.03440000E-04, 1.21390000E-02, 6.12290000E-02, 5.37320000E-01, 2.48070000E+00]]])
scatter_matrix = np.transpose(scatter_matrix)
mat0_xsdat.set_scatter_matrix(scatter_matrix, temperature=294.)
#Data for Material 1:
mat1_xsdat = openmc.XSdata('mat_1', groups)
mat1_xsdat.order = 0
mat1_xsdat.set_total([ 1.77949000E-01, 3.29805000E-01, 4.80388000E-01, 5.54367000E-01, 3.11801000E-01
, 3.95168000E-01, 5.64406000E-01], temperature=294.)
mat1_xsdat.set_absorption([ 8.02455708E-03, 3.71759686E-03, 2.67688000E-02, 9.62356000E-02, 3.00197400E-02
, 1.11260200E-01, 2.82780200E-01], temperature=294.)
mat1_xsdat.set_fission([ 7.21206000E-03, 8.19301000E-04, 6.45320000E-03, 1.85648000E-02, 1.78084000E-02
, 8.30348000E-02, 2.16004000E-01], temperature=294.)
mat1_xsdat.set_nu_fission([ 2.00599843E-02, 2.02730297E-03, 1.57059918E-02, 4.51830102E-02, 4.33420839E-02
, 2.02090096E-01, 5.25710535E-01], temperature=294.)
mat1_xsdat.set_chi([ 5.87910000E-01, 4.11760000E-01, 3.39060000E-04, 1.17610000E-07, 0.00000000E+00
, 0.00000000E+00, 0.00000000E+00], temperature=294.)
scatter_matrix = np.array(\
[[[ 1.27537000E-01, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00],
[ 4.23780000E-02, 3.24456000E-01, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00],
[ 9.43740000E-06, 1.63140000E-03, 4.50940000E-01, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00],
[ 5.51630000E-09, 3.14270000E-09, 2.67920000E-03, 4.52565000E-01, 1.25250000E-04, 0.00000000E+00, 0.00000000E+00],
[ 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 5.56640000E-03, 2.71401000E-01, 1.29680000E-03, 0.00000000E+00],
[ 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 1.02550000E-02, 2.65802000E-01, 8.54580000E-03],
[ 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 0.00000000E+00, 1.00210000E-08, 1.68090000E-02, 2.73080000E-01]]])
scatter_matrix = np.transpose(scatter_matrix)
mat1_xsdat.set_scatter_matrix(scatter_matrix, temperature=294.)
#Create the cross sections hdf5 file:
mg_cross_sections_file = openmc.MGXSLibrary(groups)
mg_cross_sections_file.add_xsdata(mat0_xsdat)
mg_cross_sections_file.add_xsdata(mat1_xsdat)
mg_cross_sections_file.export_to_hdf5('cross_sections.h5')
#Assign each cross section to a separate material:
materials = {}
for xs in ['mat_0','mat_1']:
materials[xs] = openmc.Material(name=xs)
materials[xs].set_density('macro', 1.)
materials[xs].add_macroscopic(xs)
#Create the materials file for this specification:
materials_file = openmc.Materials(materials.values())
materials_file.cross_sections = 'cross_sections.h5'
materials_file.export_to_xml()
#The following creates the settings and sets the energy mode to multi-group:
settings_file = openmc.Settings()
settings_file.energy_mode = 'multi-group'
#The user should complete the OpenMC input below by specifying geometry, settings, etc.
#grid to be used
xgrid = {}
xgrid[1] = openmc.XPlane(0, boundary_type='reflective')
xgrid[2] = openmc.XPlane(1.2600000000000000E-01)
xgrid[3] = openmc.XPlane(2.5200000000000000E-01)
xgrid[4] = openmc.XPlane(3.7800000000000000E-01)
xgrid[5] = openmc.XPlane(5.0400000000000000E-01)
xgrid[6] = openmc.XPlane(6.3000000000000000E-01, boundary_type='reflective')
ygrid = {}
ygrid[1] = openmc.YPlane(0, boundary_type='reflective')
ygrid[2] = openmc.YPlane(1.2600000000000000E-01)
ygrid[3] = openmc.YPlane(2.5200000000000000E-01)
ygrid[4] = openmc.YPlane(3.7800000000000000E-01)
ygrid[5] = openmc.YPlane(5.0400000000000000E-01)
ygrid[6] = openmc.YPlane(6.3000000000000000E-01, boundary_type='reflective')
bound1 = openmc.ZPlane(0, boundary_type='reflective')
bound2 = openmc.ZPlane(6.3000000000000000E-01, boundary_type='reflective')
#cells grid
cellgrid = {}
#row 1
cellgrid[1,1] = +xgrid[1] & -xgrid[2] & +ygrid[1] & -ygrid[2] & +bound1 & -bound2
cellgrid[1,2] = +xgrid[2] & -xgrid[3] & +ygrid[1] & -ygrid[2] & +bound1 & -bound2
cellgrid[1,3] = +xgrid[3] & -xgrid[4] & +ygrid[1] & -ygrid[2] & +bound1 & -bound2
cellgrid[1,4] = +xgrid[4] & -xgrid[5] & +ygrid[1] & -ygrid[2] & +bound1 & -bound2
cellgrid[1,5] = +xgrid[5] & -xgrid[6] & +ygrid[1] & -ygrid[2] & +bound1 & -bound2
#row 2
cellgrid[2,1] = +xgrid[1] & -xgrid[2] & +ygrid[2] & -ygrid[3] & +bound1 & -bound2
cellgrid[2,2] = +xgrid[2] & -xgrid[3] & +ygrid[2] & -ygrid[3] & +bound1 & -bound2
cellgrid[2,3] = +xgrid[3] & -xgrid[4] & +ygrid[2] & -ygrid[3] & +bound1 & -bound2
cellgrid[2,4] = +xgrid[4] & -xgrid[5] & +ygrid[2] & -ygrid[3] & +bound1 & -bound2
cellgrid[2,5] = +xgrid[5] & -xgrid[6] & +ygrid[2] & -ygrid[3] & +bound1 & -bound2
#row 3
cellgrid[3,1] = +xgrid[1] & -xgrid[2] & +ygrid[3] & -ygrid[4] & +bound1 & -bound2
cellgrid[3,2] = +xgrid[2] & -xgrid[3] & +ygrid[3] & -ygrid[4] & +bound1 & -bound2
cellgrid[3,3] = +xgrid[3] & -xgrid[4] & +ygrid[3] & -ygrid[4] & +bound1 & -bound2
cellgrid[3,4] = +xgrid[4] & -xgrid[5] & +ygrid[3] & -ygrid[4] & +bound1 & -bound2
cellgrid[3,5] = +xgrid[5] & -xgrid[6] & +ygrid[3] & -ygrid[4] & +bound1 & -bound2
#row 4
cellgrid[4,1] = +xgrid[1] & -xgrid[2] & +ygrid[4] & -ygrid[5] & +bound1 & -bound2
cellgrid[4,2] = +xgrid[2] & -xgrid[3] & +ygrid[4] & -ygrid[5] & +bound1 & -bound2
cellgrid[4,3] = +xgrid[3] & -xgrid[4] & +ygrid[4] & -ygrid[5] & +bound1 & -bound2
cellgrid[4,4] = +xgrid[4] & -xgrid[5] & +ygrid[4] & -ygrid[5] & +bound1 & -bound2
cellgrid[4,5] = +xgrid[5] & -xgrid[6] & +ygrid[4] & -ygrid[5] & +bound1 & -bound2
#row 5
cellgrid[5,1] = +xgrid[1] & -xgrid[2] & +ygrid[5] & -ygrid[6] & +bound1 & -bound2
cellgrid[5,2] = +xgrid[2] & -xgrid[3] & +ygrid[5] & -ygrid[6] & +bound1 & -bound2
cellgrid[5,3] = +xgrid[3] & -xgrid[4] & +ygrid[5] & -ygrid[6] & +bound1 & -bound2
cellgrid[5,4] = +xgrid[4] & -xgrid[5] & +ygrid[5] & -ygrid[6] & +bound1 & -bound2
cellgrid[5,5] = +xgrid[5] & -xgrid[6] & +ygrid[5] & -ygrid[6] & +bound1 & -bound2
#defining the cells
cells = {}
#row 1
cells[1,1] = openmc.Cell(region=cellgrid[1,1],fill=materials['mat_1'])
cells[1,2] = openmc.Cell(region=cellgrid[1,2],fill=materials['mat_1'])
cells[1,3] = openmc.Cell(region=cellgrid[1,3],fill=materials['mat_1'])
cells[1,4] = openmc.Cell(region=cellgrid[1,4],fill=materials['mat_1'])
cells[1,5] = openmc.Cell(region=cellgrid[1,5],fill=materials['mat_0'])
#row 2
cells[2,1] = openmc.Cell(region=cellgrid[2,1],fill=materials['mat_1'])
cells[2,2] = openmc.Cell(region=cellgrid[2,2],fill=materials['mat_1'])
cells[2,3] = openmc.Cell(region=cellgrid[2,3],fill=materials['mat_1'])
cells[2,4] = openmc.Cell(region=cellgrid[2,4],fill=materials['mat_1'])
cells[2,5] = openmc.Cell(region=cellgrid[2,5],fill=materials['mat_0'])
#row 3
cells[3,1] = openmc.Cell(region=cellgrid[3,1],fill=materials['mat_1'])
cells[3,2] = openmc.Cell(region=cellgrid[3,2],fill=materials['mat_1'])
cells[3,3] = openmc.Cell(region=cellgrid[3,3],fill=materials['mat_1'])
cells[3,4] = openmc.Cell(region=cellgrid[3,4],fill=materials['mat_1'])
cells[3,5] = openmc.Cell(region=cellgrid[3,5],fill=materials['mat_0'])
#row 4
cells[4,1] = openmc.Cell(region=cellgrid[4,1],fill=materials['mat_1'])
cells[4,2] = openmc.Cell(region=cellgrid[4,2],fill=materials['mat_1'])
cells[4,3] = openmc.Cell(region=cellgrid[4,3],fill=materials['mat_0'])
cells[4,4] = openmc.Cell(region=cellgrid[4,4],fill=materials['mat_0'])
cells[4,5] = openmc.Cell(region=cellgrid[4,5],fill=materials['mat_0'])
#row 5
cells[5,1] = openmc.Cell(region=cellgrid[5,1],fill=materials['mat_0'])
cells[5,2] = openmc.Cell(region=cellgrid[5,2],fill=materials['mat_0'])
cells[5,3] = openmc.Cell(region=cellgrid[5,3],fill=materials['mat_0'])
cells[5,4] = openmc.Cell(region=cellgrid[5,4],fill=materials['mat_0'])
cells[5,5] = openmc.Cell(region=cellgrid[5,5],fill=materials['mat_0'])
#adding all cells to the universe
root_universe = openmc.Universe(cells=(cells[1,1]\
,cells[1,2]\
,cells[1,3]\
,cells[1,4]\
,cells[1,5]\
,cells[2,1]\
,cells[2,2]\
,cells[2,3]\
,cells[2,4]\
,cells[2,5]\
,cells[3,1]\
,cells[3,2]\
,cells[3,3]\
,cells[3,4]\
,cells[3,5]\
,cells[4,1]\
,cells[4,2]\
,cells[4,3]\
,cells[4,4]\
,cells[4,5]\
,cells[5,1]\
,cells[5,2]\
,cells[5,3]\
,cells[5,4]\
,cells[5,5]))
#additional settings info
geometry = openmc.Geometry(root_universe)
geometry.export_to_xml()
# Create a point source
boxspace = openmc.stats.Box(lower_left=[0,0,0], upper_right=[6.3000000000000000E-01,6.3000000000000000E-01,6.3000000000000000E-01])
source = openmc.Source(space=boxspace)
#settings
settings_file.source = source
settings_file.batches = 150
settings_file.inactive = 15
settings_file.generations_per_batch = 15
settings_file.particles = 1000
settings_file.export_to_xml()
#run openmc
openmc.run()
As mentioned, the real issue isn’t necessarily the simulation runtime, but rather the prohibitively long source initialization at the beginning. In fact, this calculation took around 960 seconds, but less than 30 seconds was actually spent running the simulation. The other 15 minutes were mostly spent hung up on the Initializing source particles...portion…
I’ve found I can avoid this issue a bit by greatly reducing the settings_file.particles option to say 10, and then simply ramping up the number of batches to try and still get good statistics, but the problem with this is that with only 10 particles per generation, it will inevitably result in a generation with zero fission and the run will prematurely terminate. I’ve run other problems with far more particles and not had this issue, but they all had many less distinct regions. This leads me to suspect the issue may be with the number of regions that this problem has.
So, is there something that OpenMC is doing with the number of regions at the beginning that is causing this slowdown? If so, is there a way around it? Besides reducing the number of regions that is… In fact I was hoping to do problems on a 50x50 grid, but if this slowdown is unavoidable and potentially made worse by adding more regions, that seems like an unlikely prospect…