I have a similar problem to another user but the fixes outlined in that post do not solve my issue. I am starting to use meshes to plot data in 2d and 3d forms. For CSG geometry I have had no problems at all but since I have begun working with CAD it has made it much more confusing. I am getting very strange results when plotting scattering, absorption and flux through my basic aluminum sphere (plots will be attached at the bottom) Essentially when using a regular mesh, my plots for scattering and absorption show values only at the source and nowhere else. When I use an unstructured mesh, my tallies are completely empty and I cannot find any data, although the simulation runs fine with a reasonable leakage fraction. I will post my code below and the plots; any help would be greatly appreciated.
I get the following message when running with unstructured mesh:
WARNING: Output for a MOAB mesh (mesh 24) was requested but will not be
written. Please use the Python API to generated the desired VTK
tetrahedral mesh.
import cadquery as cq
from cadquery import Assembly
from cad_to_dagmc import CadToDagmc
%matplotlib inline
import os
import openmc
import matplotlib.pyplot as plt
import openmc.lib
from stl_to_h5m import stl_to_h5m
assert(openmc.lib._dagmc_enabled())
Aluminum = openmc.Material(name='aluminum')
Aluminum.add_element('Al', 1, percent_type='ao')
#Aluminum.set_density('g/cm3', 2.7)
mats = openmc.Materials()
mats.append(Aluminum)
mats.append(air)
mats.export_to_xml()
from cad_to_dagmc import CadToDagmc
my_model = CadToDagmc()
my_model.add_stp_file("Sphere.step", material_tags=["aluminum"])
my_model.export_dagmc_h5m_file()
dag_univ = openmc.DAGMCUniverse(filename = 'dagmc.h5m')
# creates an edge of universe boundary at a large radius
vac_surf = openmc.Sphere(r=10000, surface_id=9999, boundary_type="vacuum")
# specifies the region as below the universe boundary
region = -vac_surf
# creates a cell from the region and fills the cell with the dagmc geometry
containing_cell = openmc.Cell(cell_id=9999, region=region, fill=dag_univ)
model = openmc.Model()
model.geometry = openmc.Geometry(root=[containing_cell])
model.geometry.export_to_xml()
settings = openmc.Settings()
settings.batches = 10
settings.inactive = 0
settings.particles = 500
settings.run_mode = 'fixed source'
source = openmc.Source()
source.space = openmc.stats.Point((0,9,0))
source.angle = openmc.stats.Isotropic()
source.energy = openmc.stats.Discrete([14e3], [1])
source.particle = 'neutron'
settings.source = source
settings.export_to_xml()
mesh = openmc.RegularMesh().from_domain(
model.geometry, # the corners of the mesh are being set automatically to surround the geometry
dimension=[100, 1, 100] # only 1 cell in the Y dimension (change to cube for vtk file)
)
#mesh = openmc.UnstructuredMesh("sphere.h5m", 'moab')
tallies = openmc.Tallies()
# Create mesh filter for tally
mesh_filter = openmc.MeshFilter(mesh)
mesh_tally_1 = openmc.Tally(name='absorption')
mesh_tally_1.filters = [mesh_filter]
mesh_tally_1.scores = ['absorption']
tallies.append(mesh_tally_1)
mesh_tally_2 = openmc.Tally(name='scattering')
mesh_tally_2.filters = [mesh_filter]
mesh_tally_2.scores = ['scatter']
tallies.append(mesh_tally_2)
tallies.export_to_xml()
mesh_tally_3 = openmc.Tally(name='flux')
mesh_tally_3.filters = [mesh_filter]
mesh_tally_3.scores = ['flux']
tallies.append(mesh_tally_3)
tallies.export_to_xml()
model = openmc.model.Model(model.geometry, mats, settings, tallies)
# deletes old files
!rm summary.h5
!rm statepoint.*.h5
!rm plot_14.png
# runs the simulation
sp_filename = model.run()
# open the results file
results = openmc.StatePoint(sp_filename)
my_tally = results.get_tally(scores=['absorption'])
my_slice = my_tally.get_slice(scores=['absorption'])
my_slice.mean.shape = (mesh.dimension[0], mesh.dimension[2]) # setting the resolution to the mesh dimensions
fig = plt.subplot()
plt.show(fig.imshow(my_slice.mean, extent=[-250,250,-250,250]))
