Hℹ guys,
I use the cylindrical mesh to tally the fission rate and I want to see its distribution. Are there a method
to get it ?
Any help will be thanks!
Best regards!
Zhuangli
Hi @lz152,
I’ve been meaning to finish up a feature that will apply data (from a mesh tally) to any 3D mesh and write a VTK file from the Python API. It’s still in a working state, but if you’re up for working with a source installation, here’s the branch if you want to give it a try!
To apply data from a flux mesh tally you would do something like:
# setup everything other than tallies
cyl_mesh = openmc.CylindricalMesh()
cyl_mesh.r_grid = np.linspace(0, 5, 8)
cyl_mesh.phi_grid = np.linspace(0, 2 * np.pi, 4)
cyl_mesh.z_grid = np.linspace(0, 2, 6)
mesh_filter = openmc.MeshFilter(cyl_mesh)
tally = openmc.Tally(name='my tally')
tally.filters = [mesh_filter]
tally.scores = ['flux']
tallies = openmc.Tallies([tally])
tallies.export_to_xml()
openmc.run()
Run simulation
sp = openmc.StatePoint('statepoint.100.h5')
tally = sp.get_tally(name='my tally')
flux_mean = tally.mean.reshape(*cyl_mesh.dimension)
data = {'flux': flux_mean}
mesh.write_vtk_mesh('cyl_mesh.vtk', data=data)
Note that volume and source normalizations will be required to get the fission rate correct.
Hope this helps!
-Patrick
EDIT 5-02-22: correcting value of tally variable when getting the tally from the statepoint file.
Hi Shriwise,
Thank you. it is a nice code for cylindrical-mash-tally to get a visualization.
I will try to use it!
Anything else:
Dose the cylindrical-mesh-tool available with a source installation?
Best regards!
ZhuangLi
No problem! Let me know how it goes!
Dose the cylindrical-mesh-tool available with a source installation?
If you do a source installation with that branch, cylindrical meshes will be available, yes. I hope I understood your question correctly.
Yeah.
You are right! 
Thank you 
~
1 Like
Hello @pshriwise,
Sorry to revive this topic. The mesh filter documentation says that only regular, rectangular meshes can be used together with the mesh filter. And i indeed found strange results while inadvertently using a cylindrical mesh in a mesh filter. I would like to ask if that is really the case and If it is, how i could define a cylindrical mesh tally?
Best regards,
Artur Souza
@Artur You should be able to use a cylindrical mesh in conjunction with the MeshFilter class. The documentation there should be updated so thanks for pointing that out! Feel free to elaborate on what you believe is strange about the results you’re seeing.
2 Likes
Noting here that the documentation on this class has now been updated in this PR thanks to @paulromano.
1 Like
Hello @paulromano,
I found strange results when trying to reproduce the radial power distribution on a simple pin cell model using a cylindrical mesh. I defined the cylindrical mesh tally as follows:
cylindrical_mesh = openmc.CylindricalMesh()
cylindrical_mesh.r_grid = np.linspace(0.0, fuel_outer_radius.r, 20)
cylindrical_mesh.phi_grid = [0.0, 2*np.pi]
cylindrical_mesh.z_grid = [-height/2, height/2]
cylindrical_mesh_filter = openmc.MeshFilter(cylindrical_mesh)
T4 = openmc.Tally(4, name = 'Cylindrical Mesh Tally')
T4.filters = [cylindrical_mesh_filter]
T4.scores = ['flux', 'fission', 'kappa-fission']
And extracted the results with the following code:
power = cylindrical_mesh_results.get_values(scores = ['kappa-fission'], value = 'mean').flatten()
power_std_dev = cylindrical_mesh_results.get_values(scores = ['kappa-fission'], value = 'std_dev').flatten()
radii_list = np.linspace(0, fuel_radius.r, 20)
middle_points = (np.diff(radii_list)/2) + radii_list[:len(power)]
mesh_volumes = cylindrical_mesh_results.filters[0].mesh.volumes.flatten()
RADIAL = plt.figure(figsize = (8, 6))
GRADIAL = RADIAL.add_subplot(111)
GRADIAL.errorbar(middle_points, power/mesh_volumes, yerr = power_std_dev/mesh_volumes, fmt = 'ro', markersize = 2, ecolor = 'gray', elinewidth = 1, capsize = 2.5)
GRADIAL.set_title('Radial Power Distribution', size = 16, pad = 12)
GRADIAL.set_xlabel('Radial position $(cm)$', math_fontfamily = 'cm', size = 12)
GRADIAL.set_ylabel('Power density $(eV/cm^{3}*particle)$', math_fontfamily = 'cm', size = 12)
That returns the following graph:
This radial distribution, however, disagrees with the well known behavior of a fuel pellet having a higher power density at its periphery than at its center. It also disagrees with the cartesian mesh results from the same simulation:
I know its hard to see the different shades of red in the figure, but i have numerically checked and the power indeed increases as we move from the center to the boundary of the fuel pellet.
Best regards,
Artur Souza
@Artur It looks like there may be a bug with cylindrical mesh tallies when the tracklength estimator is used. If you change the estimator to collision:
T4.estimator = 'collision'
and rerun your problem, you should see the correct behavior. I’ll try to do some digging on this to see if I can figure out the root cause and get it fixed. Thanks for the detailed message reporting this!
1 Like
I have just rerun the simulation using the collision estimator and i got the correct results. Thank you very much for the help.
@Artur Just a heads up – I submitted a bug fix for the underlying issue here, so in our next version of OpenMC you should be able to use the default tracklength estimator with no issues.
1 Like
Noticed this thread and wanted to mention I’ve been tinkering with a cylindrical mesh slicer that can be used with matplotlib to plot R-Z and Phi-Z axis slices of a cylindrical mesh.
1 Like
Thank you for the heads up @paulromano. And thank you for the suggestion @Shimwell. I will definitely be checking your cylindrical mesh plotter.