Plotting axial fission distribution

User Support
1

Hello. Thanks for any help in advance. I am trying to plot axial fission distribution for a MSR. I have successfully plotted a 2D radial fission distribution.

Here is my code for axial view:

tallies=openmc.Tallies()
mesh_axial = openmc.RegularMesh()
mesh_axial.dimension = [10, 1, 10]
mesh_axial.lower_left = [-90, -90, -100]
mesh_axial.upper_right = [90, +90, 100]
mesh_f_axial = openmc.MeshFilter(mesh_axial)

tally_axial = openmc.Tally(name=‘axial’)
tally_axial.filters = [mesh_f_axial]
tally_axial.scores = [‘fission’]

tallies.append(tally_axial)

settings=DefineSetting()
model = openmc.Model(geometry=geometry,materials=materials, settings=settings)
model.tallies = tallies

statepoint_filename = model.run()
sp = openmc.StatePoint(statepoint_filename)

the above code runs a successful openmc run and generates that settings, statepoint, tallies, etc.

import numpy as np
import matplotlib.pyplot as plt

Load the last statepoint file and keff value

with openmc.StatePoint(statepoint_filename) as sp:
# Get the OpenMC pin power tally data
mesh_tally = sp.get_tally(name=‘axial’)
fission_rates = mesh_tally.get_values(scores=[‘fission’])

Reshape array to 2D for plotting

fission_rates.shape = mesh_axial.dimension

Normalize to the average pin power

fission_rates /= np.mean(fission_rates[fission_rates > 0.])

Force zeros to be NaNs so their values are not included when matplotlib calculates

the color scale

fission_rates[fission_rates == 0.] = np.nan

Plot the pin powers and the fluxes

plt.figure()
plt.imshow(fission_rates, interpolation=‘none’, cmap=‘jet’, origin=‘lower’)
plt.colorbar()
plt.title(‘Fission distribution axial’)
plt.savefig(“FissionAxial”)
plt.show()

the above code generates the following error:

TypeError: Invalid shape (10, 1, 10) for image data

I’ve tried some other variations, but I’m clearing flailing - I don’t understand what the mesh definitions are providing to the plotting routine.

Thanks again for help! Brian

1 Like
2

Hi Brian

Thanks for posting

It might be worth taking a look at this pull request to openmc where I am attempting to make plotting RegularMesh tallies easier.

If you are super keen you could install this branch and try it out on your model

Also if you have any feedback for the pull request that would also be welcome

1 Like
3

Hi Brian,

Neat application! I think the problem you’re running into is the fission_rates variable doesn’t have a shape that’s valid for imshow. It expects a 2D array of data values or a 3D array where the third axis contains RGB values for color. Documentation on that is here.

A minimal change to get a plot from your code would be:

# removes all axes of size 1 from the array
fission_rates = np.squeeze(fission_rates)

plt.figure()
plt.imshow(fission_rates, interpolation=‘none’, cmap=‘jet’, origin=‘lower’)
plt.colorbar()
plt.title(‘Fission distribution axial’)
plt.savefig(“FissionAxial”)
plt.show()