Polyethylene (HDPE) cross section for neutron moderation

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?

HDPE_material.add_elements_from_formula(“C2H4”, percent_type=“ao”)
HDPE_material.set_density(“g/cm3”, 0.95)

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

HDPE_material.add_element(‘C’, 2)
HDPE_material.add_element(‘H’, 4)
HDPE_material.set_density(“g/cm3”, 0.95)

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

Here is a collection of material definitions which includes a polyethylene neutronics_material_maker/pnnl_materials.json at 55fa45b6811c19a34874ab9f501c253a56df6051 · fusion-energy/neutronics_material_maker · GitHub

1 Like

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.

I quite like this calculate_cexs function for getting the cross section data from materials.

The plot_xs function in that same file is also handy if you don’t want the data and just want a plot.

1 Like

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:

2 Likes

@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?

@nickschw: because the bound scattering cross section of H1 is ~82 b. That is a know fact, but you can also see it here:

https://www.ncnr.nist.gov/resources/n-lengths/list.html

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.

1 Like

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

Hello,

indeed this needs to be followed-up on.

I have done some testing and can confirm that the OpenMC Official Lucite TSL is faulty.

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)

1 Like

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 will check it and start a pull request.

2 Likes

Hello,

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.

mat_water = openmc.Material(name=“Water”)
mat_water.add_element(“H”, 2.0, percent_type=“ao”)
mat_water.add_element(“O”, 1.0, percent_type=“ao”)
mat_water.set_density(“g/cc”, 1.0)

mat_hdpe = openmc.Material(name=“HDPE”)
mat_hdpe.add_element(“C”, 2)
mat_hdpe.add_element(“H”, 4)
mat_hdpe.set_density(“g/cm3”, 0.95)

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!

Hello pamk,

the water material will need:
c_H_in_H2O
mat_water.add_s_alpha_beta('c_H_in_H2O')

the polyethylene will need:
c_H_in_CH2
mat_hdpe.add_s_alpha_beta('c_H_in_CH2')

If those are not used, simulations of thermal neutron fluxes will be incorrect.
The effects of S(a,b) for these two materials is significant.

For now you should use the LANL ENDF S(a,b) libraries:
https://openmc.org/lanl-data-libraries/

As there is problems with the OpenMC native official libraries:
https://openmc.org/official-data-libraries/

This is being actively addressed.

Cheers,
Perry

Dear Perry,

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!