Tally shaping error for post processing using statepoint file

Hi OpenMC community-

I am a senior undergrad working on a capstone project and would like to ultimately display our flux and heating tallies for our project. We have been able to successfully run OpenMC using CSG geometry but are confused as how to display tally data. We are running into an error when trying to shape our tally results. It is our understanding that all tallies must be the same size prior to plotting using matplot lib, especially as we are pulling the same result from different cells.

We have 16 cells defined in our geometry, with one being a void. We are trying to gather tally results on a cell-by-cell basis, with the end goal of plotting an image of our geometry for each individual tally (one image using flux on all components/cells, one image using heating on all components/cells, etc.)

We have set up our tallies as follows:

flibe_1_filter = openmc.CellFilter([1])
flibe_1_tally = openmc.Tally()
flibe_1_tally.filters = [flibe_1_filter]
flibe_1_tally.scores = ['flux','heating' ,'absorption', '(n,a)', '(n,2a)', '(n,3a)', '(n,gamma)', 'He3-production']

The issue we are having is after running openMC. When trying to get flux data for plotting purposes, we would like to shape our data for all cells into a 100 x 100 array. This way, we can plot them and display an image of the flux, heating, etc.

steel_1_tally = sp.get_tally(scores=['flux'], id=steel_1_tally.id)
steel_1_flux = steel_1_tally.mean
steel_1_flux.shape = (100, 100)
print(steel_1_flux.shape)
steel_1_flux = steel_1_flux[::-1, :]

On line one above, you can see we pull in the flux score and extract the mean value. We then attempt to shape the tally. However, we are receiving an error when trying to shape the results. When printing the size of the array, we get (1, 1, 2). When we try to force all tallies to be (1, 1, 2) and plot them, we get this error:

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

Must we shape the data to plot the results correctly? Does this indicate a computational issue with the model, or is there another way to get plottable results?

Hi there. Cell tallies are great for getting a result on a cell. You might be interesting in the mesh tallies. In particular a 2d RegularMesh tallies are particularly easy to plot with matplotlib.

Here is an example

Hi @Shimwell -

Thanks for your response! This is a great help. We are giving it a try. We have also looked at some previous posts for guidance, and have formed our mesh as follows:

mesh = openmc.RegularMesh().from_domain(geometry)
#mesh = openmc.RegularMesh()
mesh.lower_left = [-120, -40, 0]
mesh.upper_right = [250,100,250]
dimension=[100, 1, 100]
mesh_filter = openmc.MeshFilter(mesh)
mesh_tally = openmc.Tally(name='tallies_on_mesh')
mesh_tally.filters = [mesh_filter]
mesh_tally.scores = ['flux','heating' ,'absorption', '(n,a)', '(n,2a)', '(n,3a)', '(n,gamma)', 'He3-production', 'scatter', 'elastic']  # you can have multiple scores on a single tally, this is more memory efficient that having new tallies each score
tallies = openmc.Tallies([mesh_tally])

We define lower left and upper right coordinates to account for the fact that we do have geometry which references a negative sense.

openmc.run() does seem to work, but then it does produce a similar issue to what we had before upon post processing when executing a reshaping command:

my_slice.mean.shape = (mesh.dimension[0], mesh.dimension[2]) # setting the resolution to the mesh dimensions

which produces:

ValueError: cannot reshape array of size 1000 into shape (10,10)

the size of our array produced is (1000, 1, 10) similarly to the cell tally structure. Any way to resolve this?

1 Like

This problem has been solved. The issue with the meshing was an incorrect placement of bottom left and top right coordinates for mesh volumes. The mesh must encompass the entire model. Make sure to also use parenthesis in place of square brackets for the mesh.lower, upper, and dimensions. The corrected code is as follows:

mesh = openmc.RegularMesh()
mesh.lower_left = (-125, -205, 0)
mesh.upper_right = (250,400,370)
mesh.dimension = (350, 1, 350)
mesh_filter = openmc.MeshFilter(mesh)
mesh_tally = openmc.Tally(name='tallies_on_mesh')
mesh_tally.filters = [mesh_filter]
mesh_tally.scores = ['flux','heating' ,'absorption', '(n,a)', '(n,t)', '(n,gamma)', 'scatter', 'elastic'] 
tallies = openmc.Tallies([mesh_tally])

Slicing and plotting the mesh is successful.