Point Detector Tally

[itay-space]

Description

This repository contains the development and benchmarking of a relativistic point detector within the OpenMC framework. The detector efficiently estimates neutron flux at a point using full relativistic kinematic calculations, making it ideal for high-energy particle simulations. This first implementation of a point detector tally in OpenMC is validated through single and multiple collision experiments, showing excellent agreement with conventional methods. This enhancement significantly improves OpenMC's capabilities for complex shielding analyses and other scenarios where traditional Monte Carlo methods may be inefficient. For further details, refer to the paper "Development and benchmarking of a point detector in OpenMC" here.

Use Example

To configure a point tally in OpenMC, follow the standard procedure for creating a tally but specify the estimator as "point" and provide the coordinates for positioning and the exclusion sphere radius $R_0$ using a list of x, y, z, $R_0$ values. The following code snippet demonstrates how to set up a point detector located at coordinates (0, 0, 1) with an exclusion sphere radius of 0.01 cm, along with an energy filter, using the OpenMC Python API:

point_detector = openmc.Tally(name="pDet")
point_detector.scores = ["flux"]  
point_detector.filters = [energy_filter]  
point_detector.estimator = 'point'  
point_detector.positions = [0, 0, 1 , 0.01]  

Citing

If you use the point detector in your research, please consider giving proper attribution by citing both the following publications:

1. Itay Horin, Or Alon, Arik Kreisel, Tsviki Y. Hirsh, Tamar Dayan, and Hallel Fuks, "Development and Benchmarking of a Point Detector in OpenMC," Annals of Nuclear Energy, 220, 111497 (2025).
   This paper presents the development, implementation, and validation of the point detector tally in OpenMC.

2. Paul K. Romano, Nicholas E. Horelik, Bryan R. Herman, Adam G. Nelson, Benoit Forget, and Kord Smith, "OpenMC: A State-of-the-Art Monte Carlo Code for Research and Development," Ann. Nucl. Energy, 82, 90--97 (2015).
   This citation is for the OpenMC code itself.

Checklist

• I have performed a self-review of my own code
• I have run clang-format (version 15) on any C++ source files (if applicable)
• I have followed the style guidelines for Python source files (if applicable)
• I have made corresponding changes to the documentation (if applicable)
• I have added tests that prove my fix is effective or that my feature works (if applicable)

[gridley]

Wow, there's a lot of new stuff going on here. My main question is: why do you think that point estimates should be interpreted as an estimating method for tallies in general? It seems like in fact, this is a special estimator that can only be applied for the point estimate filter. This special filter can only allow one type of estimator, the point.

Also, there is a substantial amount of undocumented code here, and no accompanying commentary on what's being done here in the PR documentation or docs changes. This is a nontrivial feature, so would need to be described a bit more. For example, I'm unsure what this is talking about when you say this is relativistic.

There is a fair amount of commented out debug printing. That would have to be removed before merging.

There must also be some limitations to this: for example, is this correct for treating resonance upscatter with DBRC or RVS? Getting the correctly normalized PDF value there is kind of nontrivial.

[gridley]

This doesn't really make sense as a way to architect this feature. Semantically, the estimate of a quantity at a point is an entirely different problem from the different approaches to estimating otherwise. To be specific, in this case you are calculating contributions to a tally outside the usual process of filtering the state of a particle.

IMO, this should be classified as an entirely different type of tally similar to how pulse height tallies were done. POINT should not be an option in TallyEstimator, since this would allow someone to make a mesh tally that has the POINT estimator type. That doesn't really make sense, if you ask me. This estimator is fundamentally tied to one specific type of filter, a (set of) dirac delta(s) in space.

[itay-space]

We debated whether to implement this method as a filter or an estimator. Ultimately, we decided on using an estimator because we propose here a different method of estimating flux in space. However, we recognise that this decision might not be optimal given our limited familiarity with the overall structure of OpenMC. We would greatly appreciate your input on this matter.

[itay-space]

Wow, there's a lot of new stuff going on here. My main question is: why do you think that point estimates should be interpreted as an estimating method for tallies in general? It seems like in fact, this is a special estimator that can only be applied for the point estimate filter. This special filter can only allow one type of estimator, the point.

Also, there is a substantial amount of undocumented code here, and no accompanying commentary on what's being done here in the PR documentation or docs changes. This is a nontrivial feature, so would need to be described a bit more. For example, I'm unsure what this is talking about when you say this is relativistic.

There is a fair amount of commented out debug printing. That would have to be removed before merging.

There must also be some limitations to this: for example, is this correct for treating resonance upscatter with DBRC or RVS? Getting the correctly normalized PDF value there is kind of nontrivial.

Hello @gridley ,

Thank you for your thoughtful review and feedback on the code. We have added partial documentation to some of the functions and are actively working on documenting the rest. Regarding your comment about relativistic computations, the added documentation should make this aspect clear now. Please let us know if anything is still unclear or needs further elaboration.

We have recently submitted a paper to the Annals of Nuclear Energy detailing the development and benchmarking of the code.

Regarding your question about velocity sampling: we introduced v_t, which represents the most recently sampled target velocity within the particle data class. This allows us to calculate the angular probability density function (PDF) for scattering towards the detector based on a specific sampled velocity.

Thank you again for your time and guidance.

Best regards,
Itay

[egor1abs]

Thanks for this work @itay-space! This would be good to have.
Are there any problems that you need help with to merge into the main branch?

[itay-space]

I’ve resolved all the merge conflicts and confirmed that the code compiles cleanly. I know the documentation is still lacking in parts, if there’s any specific area you’d like to see documented or clarified, let me know and I’ll add it right away.

[itay-space]

The PR is now a single commit, so the previous noise/conflict history is gone, and it’s been brought up to date with upstream via a merge to resolve the divergence. If there’s any missing documentation or gaps, just flag them and I’ll add or clarify them right away.

[GuySten]

I am interested in helping to review this feature.

This feature is really big so I suggest breaking it down into logical chunks.
Im my opinion a natural point of separation is the ability to calculate pdfs correctly ( On the cpp side ).
So I suggest starting with a PR that only implements that functionality.

@gridley, as the reviewer, what are your thoughts on the matter?

[shimwell]

I see this PR has been split into multiple PRs. I just have a couple of small question on the feature itself if you have time @itay-space

I notice in the paper there is a sentence "Future work will extend the detector’s application to gamma rays". Can I ask if that includes just a gamma if it is from the source term or does that include gammas from coupled n p transport (when settings.photon_transport=True)?

I see that the tested tallies were for neutron flux. I just want to ask if this works for neutron dose? Dose has a score of 'flux' but also has a openmc.EnergyFunctionFilter as a filter so perhaps this complicates things?

Many thanks for the work and the extra effort taken to get the feature into the OpenMC main branch

[shimwell]

Also one more question, when trying this feature I guess this branch is useful as it is the only one with all the parts included. Or would you recommend the slightly newer deploy-clean branch

[shimwell]

When running openmc I get this error with both this deploy and the deploy-clean branch

free(): corrupted unsorted chunks
Aborted (core dumped)

[itay-space]

Thanks for the questions @shimwell! The extension wasn’t limited to source gammas, it should also handle secondary gammas from coupled n,γ transport. Regarding the dose tally, there’s no issue since the point tally is implemented as an estimator, all filters (including EnergyFunctionFilter) should work the same.
For the runtime error, please try the deploy-old branch, it’s stable, though not synced with the latest OpenMC version.

[shimwell]

thanks , oh that error message might be me making a mistake somewhere else. Appears to be working, I shall post any findings. Thanks

[shimwell]

I notice some conflicts have crept into this PR, perhaps this is not important as it has been broken down the PR into separate parts.

However I think this one is the only branch (or deply-clean or deploy-old) which can run the point detector for testing at the moment.