openmc
openmc copied to clipboard
Adding support for NCrystal materials in OpenMC
After some discussion with @paulromano we cleaned up this version of OpenMC with support for NCrystal. NCrystal is a software library that provides low energy neutron physics to Monte Carlo codes. It was initially designed with Geant4 in mind, but it currently has interfaces with many other codes, and it can be called from Python, C and C++. NCrystal is developed mainly at the European Spallation Source. You can read more about NCrystal in these papers:
X.-X. Cai and T. Kittelmann, NCrystal: A library for thermal neutron transport, Computer Physics Communications 246 (2020) 106851, https://doi.org/10.1016/j.cpc.2019.07.015
T. Kittelmann and X.-X. Cai, Elastic neutron scattering models for NCrystal, Computer Physics Communications 267 (2021) 108082, https://doi.org/10.1016/j.cpc.2021.108082
X.-X. Cai, T. Kittelmann, et. al., "Rejection-based sampling of inelastic neutron scattering", Journal of Computational Physics 380 (2019) 400-407, https://doi.org/10.1016/j.jcp.2018.11.043
or watch the presentation from @tkittel, at ND2022:
https://indico.frib.msu.edu/event/52/contributions/954/
We have used NCrystal together with NJOY to generate a large library of new materials in ACE format (the NJOY-NCrystal library), but calling NCrystal directly allows to introduce in OpenMC additional physics that is not supported by the ACE format.
In particular, NCrystal allows to generate and sample on the fly scattering kernels for solid polycrystalline materials at any temperature. This is done at initialization, based on a human readable description stored in text form in the NCMAT format. For existing ENDF-6 scattering kernels, it allows the direct sampling of S(alpha, beta), which allows for continuous angular distributions (which are still not supported by the ACE format). NCrystal also support oriented materials (single crystals, highly-oriented pyrolytic graphite) and has plug-in support for other physics, currently being used to include small angle neutron scattering and inelastic magnetic scattering. All this would make OpenMC more useful to the neutron scattering community, while adding advanced nuclear data capabilities for the reactor community.
This pull request is planned as an optional feature that could be activated with cmake, using the OPENMC_USE_NCRYSTAL flag after installing and initializing NCrystal. A new method is added to the Material class in Python, to create materials from NCrystal, which takes the composition and density, creates the OpenMC material and associates it with NCrystal.
After reading the xml files, an NCrystal object is created in memory, which provides thermal scattering cross sections and scattering events (for the calculation of the scattering events, the direction of the neutron is a parameter to be able to consider oriented materials). When the thermal scattering cross section is computed, the nuclear elastic scattering component is removed, similarly to what is done in the standard S(alpha, beta) treatment.
There a couple of open topics that we could not decide, and are marked with "fixme" in the code. For instance, NCrystal uses a cache that could be implemented in OpenMC to reduce the computation time. Also, new tags for NCrystal events could be added for tallies. But, overall, we think it is ready for merging upstream.
I can add to what @marquezj wrote, that since last week NCrystal is also available on conda-forge (since @paulromano suggested this during our meeting).
Thanks for the PR @marquezj! Really excited to see this. Just wanted to give you a heads up I probably won't get to reviewing this until early next week.
@paulromano I added comments concerning the two non-obvious fixmes in the PR. Let us know if you have any specific wishes concerning this (or anything else), and we will try to update the PR.
So @paulromano, unless I overlooked something I believe @marquezj took care of all the changes requested so far (including applying clang-format).
FYI we haven't forgotten this, but some other work came in the way. We hope to spend some time on this on Friday.
I just pushed a few updates, mainly to use a new py API method I have just added in NCrystal 3.4.1, which helps getting the basic material composition in Python even for complicated (multiphase/enriched/...) materials. In the rare and somewhat academic case of materials where the same element now appears both as a natural element and as a explicit isotope, we are now using the openmc natural abundancies to break apart the natural element. This is like what we have been doing all along in Geant4 (with the C++ API), where we instead of course were asking Geant4 to provide the natural abundancies.
I also added material_id
and name
optional parameters on the Material.from_ncrystal
method, which then simply pass these two parameters on to the Material constructor.
Outstanding issues are tests, documentation and CI as requested. We will have a look at those tomorrow.
@paulromano I added a tools/ci/gha-install-ncrystal.sh
script now (can most likely be improved a bit in the future). Should we also update the gha-install.sh
and gha-install.py
files? Also, will it be possible to enable the workflows on the current PR, so we can try to get the new tests for NCrystal to work?
Update: I guess there will be a few more files to update, like perhaps gha-script.sh, ci.yml, ...). Do you prefer to take care of that, or should we try?
Hi, sorry for the delay. I added documentation and tests. @tkittel check the documentation: is brief but more can be found in the NCrystal repo and papers (I put links). I added one simple test to check it works. It will be easy to add tests with more complex config strings.
@marquezj it looks fine to me. Once we hopefully have some even better NCrystal documentation, we can go back to the openmc documentation and see if we need to update various links, or add new links. But with the current pile of items to do, it might not happen for quite some time :grin:
So the tests are unsurprisingly failing with:
import NCrystal
ModuleNotFoundError: No module named 'NCrystal'
@paulromano could you kindly suggest where you would like changes, to enable NCrystal being built? (btw. NCrystal is on conda-forge now, but I guess this CI setup is not conda-based?).
@paulromano : Is it normal that some of the checks do not start? All the started ones have now passed. Are there any remaining things we need to look at, or are you happy to merge this now?
Cheers, Thomas
@tkittel No, you don't need to worry about that. GitHub is expecting a set of test configurations that does not include the new ncrystal option, so that's why those extra ones are shown. Once this is merged, we'll just change the list of expected test configurations to include the ones with ncrystal included.
Thanks @paulromano , all changes committed as requested.
Anything else? Let me know if you want me to do a rebase of the branch on main now, or do any squashing or whatever :-)
Hi @paulromano
Despite the status saying "1 change requested", I think we did everything. Let us know if you think more is needed :smile:
Cheers, Thomas
@tkittel Apologies this fell off my radar -- will take another look soon!
@tkittel One final question -- should ncrystal_max_energy
be configurable? Right now it defaults to 5 eV and it doesn't appear that there's any way to change this.
I would postpone that for another time, since there is a common issue here with e.g. the Geant4 hooks. Once we eventually look into it, it might anyway be that we go for another solution altogether for this threshold (e.g. a per-material threshold rather than a global one).
Ok, in that case let's change it to be a global constant as opposed to a global variable in the settings::
namespace. I'll send a quick PR to your branch updating that momentarily.
Very exciting to us as well, thanks a lot @paulromano for the support and shepherding! :smile: