I have an experiment with a He3 detector to detect D-D neutrons. The detector is encased in a HDPE block that I’m trying to model. Is it valid to use the following to define the material?
I guess as a broader statement more to do with the structure of materials, one wouldn’t use plain ole C12 to model a material made from graphite. So since the cross section data in openmc doesn’t include “polyethylene,” is it valid to use it’s chemical formula instead?
Yep that would make a material from the specified ratio of elements. The formula method you have there is the most concise way to do it. You could also use the material.add_element method to build up a material but I think in this case the formula is more convenient. You get the same end result as
I looked to see if s alpha beta cross sections exist for polyethylene but couldn’t see any in ENDF/B8.0. So I think what you have is the best you can do
Here are some examples of making a material from elements if you want another option
That’s very helpful thank you. I’ve also installed the neutronics material maker – what a great resource! I do see there’s thermal scattering data c_H_in_CH2 which could be what I’m looking for. Is there a way in openmc to view the thermal scattering cross sections of .h5 files to verify it is indeed correct?
I would reference this paper to confirm its accuracy.
Hi, if you are using the standard ENDF/B-VIII.0 nuclear data files, there is a normalization error in the file for polyethylene. But actually it is really easy to fix it using openmc.data:
ch2_path = '/home/student/tsl_software/endfb-viii.0-hdf5/neutron/c_H_in_CH2.h5'
data_h5 = openmc.data.ThermalScattering.from_hdf5(ch2_path)
for T in data_h5.temperatures:
sigma_bound = data_h5.elastic.xs[T].bound_xs
data_h5.elastic.xs[T].bound_xs = sigma_bound / 2.0
print(f'Original sigma_bound: {sigma_bound} -> updated to: {(sigma_bound/2.0)}')
ch2_path_bak = ch2_path + '.bak'
os.system(f'mv {ch2_path} {ch2_path_bak}')
data_h5.export_to_hdf5(ch2_path)
This comes from one of the examples from the thermal scattering school at ESS:
@Shimwell that function is very helpful – thank you! @Ignacio wow what a great find – I wonder if anyone else has found similar discrepancies between thermal cross sections of NCrystals and ENDF? Are we sure that it’s not the NCrystal results that are off by a factor of two in the other direction?
Aside of that, it is not just a comparison between codes: I am comparing with experimental data. There I am comparing with EXFOR entry # 30748002, measured by Dawidowski and Granada.
This is not an unusual normalization problem with these materials. I need to check if this is a problem in the source data used to generate the library, or in the processing scripts.
As far as I can tell, the data that shows up in the HDF5 file (at least for the bound cross section) is exactly reproduced from the ENDF file itself, so it does not appear to be a processing issue (unless you’re talking about LEAPR→ENDF?)
Hi Paul. Sorry for hijacking the thread, maybe we should move the discussion somewhere else. I dig a bit and the problem is that some of the TSL evaluations with incoherent elastic component are not standard, and sigma_bound is multiplied by natom.
In particular, for ENDF/B-VIII.0 this affects polyethylene and lucite. A simple workaround is to check that the free atom scattering cross section derived from the bound atom scattering cross section in the elastic component matches the free atom scattering cross section from the inelastic component (or, if you have to divide by natom).
import openmc.data
import os
from io import StringIO
def get_file(url):
user_agent = '--user-agent="Mozilla/5.0"'
zipfile = 'file.zip'
file = url.split('/')[-1].replace('zip', 'dat')
os.system(f'wget --quiet -O {zipfile} {user_agent} "{url}"')
os.system(f'unzip -o {zipfile}')
return file
def read_file(file):
data_endf = openmc.data.endf.Evaluation(file)
file_obj = StringIO(data_endf.section[7, 4])
items = openmc.data.endf.get_head_record(file_obj)
items, values = openmc.data.endf.get_list_record(file_obj)
msigma = values[0]
awr = values[2]
energy_max = values[3]
natom = int(values[5])
sigma_free = msigma/natom
file_obj = StringIO(data_endf.section[7, 2])
items = openmc.data.endf.get_head_record(file_obj)
items, values = openmc.data.endf.get_tab1_record(file_obj)
sigma_bound = items[0]
return sigma_free, sigma_bound, awr, natom
endf8_urls = ['https://www-nds.iaea.org/public/download-endf/ENDF-B-VIII.0/tsl/tsl_0005_H(YH2).zip',
'https://www-nds.iaea.org/public/download-endf/ENDF-B-VIII.0/tsl/tsl_0007_H(ZrH).zip',
'https://www-nds.iaea.org/public/download-endf/ENDF-B-VIII.0/tsl/tsl_0010_H(ice-Ih).zip',
'https://www-nds.iaea.org/public/download-endf/ENDF-B-VIII.0/tsl/tsl_0037_H(CH2).zip',
'https://www-nds.iaea.org/public/download-endf/ENDF-B-VIII.0/tsl/tsl_0039_H(Lucite).zip',
'https://www-nds.iaea.org/public/download-endf/ENDF-B-VIII.0/tsl/tsl_0034_s-ch4.zip']
files = [get_file(url) for url in endf8_urls]
tol = 0.5 # tolerance in barns
print('\n Checking files:\n')
for fn in files:
sigma_free, sigma_bound, awr, natom = read_file(fn)
if (abs(sigma_free - sigma_bound/(1+1/awr)**2) < tol):
print(f'{fn} is standard: {sigma_free} b aprox. eq. {(sigma_bound/(1+1/awr)**2)} b with natom = {natom}')
elif (abs(sigma_free - sigma_bound/natom*1.0/(1+1/awr)**2) < tol):
print(f'{fn} is non standard: {sigma_free} b aprox. eq. {(sigma_bound/(1+1/awr)**2)}/{natom} b')
else:
print(f'Something is wrong with {fn}: sigma_free= {sigma_free} b, sigma_bound = {sigma_bound} b, awr = {awr}, natom = {natom} ')
which gives:
tsl_0005_H(YH2).dat is standard: 20.43634 b aprox. eq. 20.43634000403785 b with natom = 1
tsl_0007_H(ZrH).dat is standard: 20.43634 b aprox. eq. 20.47800060779665 b with natom = 1
tsl_0010_H(ice-Ih).dat is standard: 20.505 b aprox. eq. 20.065000186502324 b with natom = 1
tsl_0037_H(CH2).dat is non standard: 20.43634 b aprox. eq. 40.872679858200634/2 b
tsl_0039_H(Lucite).dat is non standard: 20.4363375 b aprox. eq. 163.49071943280254/8 b
tsl_0034_s-ch4.dat is standard: 20.43634 b aprox. eq. 20.41903460041051 b with natom = 4
This was in a solid block of lucite in my own, closed benchmark.
OpenMC w/OpenMC ENDF80 calculates thermal neutron fluxes 8% higher than they should be (compared against MCNP w/LANL ENDF80 and OpenMC /wLANL ENDF80).
I am now a little worried about the polyethylene TSL… (that one is very important, has much wider use and importance)
I will be creating an issue, hopefully I can include some benchmarks (HEU-MET-THERM-004, PU-MET-THERM-004)
Hi @yrrepy. Yes, there must be a bug somewhere in openmc.data.njoy.make_ace_thermal() because the normalization of the thermal elastic cross section for several solid S(a,b) tables is off. The lower energy limit for thermal elastic XS of hydrogen in solid moderators should be ~82 b.
I wanted to kindly follow up on this discussion. I am using OpenMC to simulate neutron shielding with water and High-Density Polyethylene (HDPE) materials, which I have defined as outlined below.
Since I am conducting a shielding simulation, I believe I should be cautious with the thermal scattering data that I’m using. Do you have any recommendations for me, by any chance? Also, should I use a function like c_H_in_H2O in order to force the thermal scattering calculations in the simulations?
Thanks for your help. I will try that and see how much it changes my results.
Generally, for shielding and dose calculations involving a fixed neutron source with energies between 10 MeV and 20 MeV, will OpenMC results be close to the real values? Should variance reduction methods be used in these types of calculations? Thanks!
This is getting off track and should be it’s own post.
S(a,b) will not be important at 10-20 MeV, but when those neutrons are thermalized (<1eV), they will be. So in shielding and dose the S(a,b) is needed.
Generally speaking OpenMC will produce results just as good as the best MC codes out there.
There exists numerous substantial uncertainties in the nuclear data at 10-20 MeV, so all of these codes will not accurately reproduce reality (the real values). It depends on the particular situation and quantity in question by how much these will deviate.
Variance reduction is often needed for shielding, but depends on the situation. Notably for thick lead or concrete.
Thanks for your response.
You are right. I will create another post right now to ask my follow-up question regarding the variance reduction method. Thanks again for your help!