Pulse height tally for photons #2452

cfichtlscherer · 2023-04-03T18:05:39Z

This pull request represents a second attempt to provide OpenMC with the capability of a pulse-height tally (see
#1881). This non-Boltzmann tally allows the simulation of all kinds of gamma detector measurements.

Thanks to @gridley and @paulromano for their help and input during the implementation. The feature is extensively validated and produces the same results as MCNP.

This tally only makes sense if photon transport is activated and works only in combination with an EnergyFilter and a CellFilter - checked at the beginning of the simulation in tally.cpp.
Furthermore, the variable vector<int> active_pulse_height_tallies; is introduced, which is used to check if the simulation tallies the pulse-height.

The basic idea for storing the information of the delivered energy over the entire history of a particle is to give the particle class a new attribute, pht_storage (a vector with the length of the number of cells in the system). In this vector, the pulse-height contributions are updated whenever the photon performs a collision or secondary particles are started. At the end of the particle's history, the function score_pulse_height_tally is called. While scoring the pulse-height, we must always iterate over all selected cells. Even if the pulse-height value is zero, this result must be scored. For this reason also the FilterCell class was changed.

I am very happy and thankful for all feedback and hope we can have this feature in the official OpenMC code soon!

When this pull-request is merged, I would happily create a sample page (jupyter notebook) that shows how to quickly and directly simulate detector responses in your system.

Update July 7, 2023
The results of OpenMC pulse-height tally simulations fit very well with the corresponding MCNP6.2 simulations (F8 tally).

For the following plot, photons with the energy of 1MeV were started from a disk source into a sodium-iodine scintillation detector. The identical model was implemented in OpenMC and MCNP6.2. The MCNP simulation used 1e8 particles and ENDFB_VI_8 crosssections, the OpenMC simulation 1e9 particles, and the Lanl_based/mcnp_endfb71 crosssections.

I am happy to share input files and more details with anyone who is interested.

gridley · 2023-04-03T18:25:14Z

I'm impressed by how few lines of code this is! Sweet!

src/particle.cpp

paulromano

Thanks for the updated PR here @cfichtlscherer. Like @gridley I'm impressed by how small the implementation is! There are some design issues we'll need to work through but I'm hopeful we can get to the point of having something that can be merged.

src/particle_data.cpp

src/tallies/tally_scoring.cpp

paulromano · 2023-04-18T19:44:43Z

Also, FYI it looks like you have tests failing in CI which you'll need to look into.

cfichtlscherer · 2023-05-01T19:15:59Z

Also, FYI it looks like you have tests failing in CI which you'll need to look into.

Yes, I think this should be solved by 09de63d

cfichtlscherer · 2023-06-08T08:39:18Z

Thanks a lot @gridley @paulromano Yes we dont need the special estimator anymore which is much better for the general code structure I think. Currently I am still fighting with another bug but I am confident to show you the next version soon :-D

paulromano

@cfichtlscherer Thanks for the updates. Nice to see that the filters are essentially unchanged now! A few more requests but this is getting really close.

One thing I would like to see before we merge: I know that you've previously done comparisons with MCNP. Would it be easy enough to repeat this comparison using the current version in this PR just to make sure things still look good?

src/particle.cpp

paulromano · 2023-06-19T16:54:04Z

src/particle.cpp

+    if (cell_born() == C_NONE)
+      cell_born() = coord(n_coord() - 1).cell;


Ah, sorry for my misunderstanding. Now I see that this code is actually exactly what's being called in event_calculate_xs, so no behavior is changing. That being said, it still feels awkward to me that this block of code is called in pht_secondary_particles. Can you move it up into event_revive_from_secondary?

src/tallies/tally.cpp

paulromano · 2023-06-19T16:55:56Z

src/tallies/tally.cpp

@@ -381,6 +393,10 @@ void Tally::add_filter(Filter* filter)
    energyout_filter_ = filters_.size();
  } else if (dynamic_cast<const DelayedGroupFilter*>(filter)) {
    delayedgroup_filter_ = filters_.size();
+  } else if (dynamic_cast<const CellFilter*>(filter)) {
+    cell_filter_ = filters_.size();
+  } else if (dynamic_cast<const EnergyFilter*>(filter)) {


Still need updates here on the type checks

tests/regression_tests/pulse_height/test.py

paulromano

Thanks for your patience @cfichtlscherer. Hopefully the last round of comments here!

src/tallies/tally.cpp

paulromano · 2023-07-07T13:15:18Z

src/tallies/tally.cpp

@@ -381,6 +393,10 @@ void Tally::add_filter(Filter* filter)
    energyout_filter_ = filters_.size();
  } else if (dynamic_cast<const DelayedGroupFilter*>(filter)) {
    delayedgroup_filter_ = filters_.size();
+  } else if (dynamic_cast<const CellFilter*>(filter)) {
+    cell_filter_ = filters_.size();
+  } else if (dynamic_cast<const EnergyFilter*>(filter)) {


Yes, I would recommend we change all of them here to not use dynamic_cast and instead rely on filt->type() checks.

cfichtlscherer · 2023-07-07T18:58:26Z

All done, think we can merge.

paulromano

Thanks for all the work getting this to the finish line!

Christopher Fichtlscherer added 2 commits April 3, 2023 11:16

pulse-height tally implementation v2.0

974a720

bugfixes

dd18f12

cfichtlscherer requested a review from paulromano as a code owner April 3, 2023 18:05

description of the score pulse-height

199e23b

Christopher Fichtlscherer added 2 commits April 3, 2023 14:25

type in score dic

38b7f70

remove blank line

89820f2

shimwell mentioned this pull request Apr 12, 2023

Pulse height tally #1881

Closed

paulromano requested changes Apr 13, 2023

View reviewed changes

src/particle.cpp Outdated Show resolved Hide resolved

Christopher Fichtlscherer added 5 commits April 13, 2023 22:08

Clang Format

c676ad7

Clang format

eeb8a49

clang

82f7b90

clang format - final

11a7ed2

clang

7c2142e

paulromano requested changes Apr 18, 2023

View reviewed changes

src/particle_data.cpp Outdated Show resolved Hide resolved

src/tallies/tally_scoring.cpp Outdated Show resolved Hide resolved

src/tallies/tally_scoring.cpp Outdated Show resolved Hide resolved

src/tallies/tally_scoring.cpp Outdated Show resolved Hide resolved

Christopher Fichtlscherer added 6 commits April 27, 2023 14:39

adding global variable pulse_height_cells

3ce4b9c

changing the value in the pht_storage at index

7cdbada

calling the right index of the pht_storage

7c9c991

small changes

01e3506

Merge branch 'develop' into pulse-height-2

a746d37

update

158c3af

cfichtlscherer requested a review from pshriwise as a code owner April 30, 2023 23:29

Christopher Fichtlscherer added 2 commits May 1, 2023 10:48

addressing issue openmc-dev#2508

d32174e

missing break statement in tally.cpp

09de63d

Christopher Fichtlscherer added 4 commits May 2, 2023 09:15

revert last changes

e401990

Merge branch 'develop' into pulse-height-2

78bd621

small change comment

deb4196

p.alive in process_advance_particle_events

1712d5c

cfichtlscherer added 6 commits June 8, 2023 14:19

remove estimator

8fab57a

remove estimator py / small changes

fb547a9

when score check if cell in pht_cells

16cbc2a

clang

7c0804a

iterate over all cells in all cell filters

48d2277

update test

cfd219a

paulromano requested changes Jun 19, 2023

View reviewed changes

cfichtlscherer added 12 commits June 20, 2023 11:39

comment photon transport / clang

7aa7c4f

clang update

969ef10

update test

63db229

cell filter check

a630e89

asking birth cell earlier

0bf609d

clang

e8a412d

clang

71b0e17

clang

3eb998f

clang

ff0f6db

clang old version

e4b574d

clange

6f7d8ad

test clang

01a03a6

paulromano requested changes Jul 7, 2023

View reviewed changes

cfichtlscherer added 2 commits July 7, 2023 17:52

changes filt->type() in tally

c796fa1

clang

b10a32e

cfichtlscherer added 3 commits July 8, 2023 15:27

dynamic_cast to FilterType

b43ebae

type

24d8cb5

type 2

ebee016

paulromano changed the title ~~Pulse height tally - second attempt~~ Pulse height tally for photons Jul 10, 2023

paulromano approved these changes Jul 10, 2023

View reviewed changes

paulromano merged commit e3cd406 into openmc-dev:develop Jul 10, 2023
16 checks passed

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Pulse height tally for photons #2452

Pulse height tally for photons #2452

cfichtlscherer commented Apr 3, 2023 •

edited

Loading

gridley commented Apr 3, 2023

paulromano left a comment

paulromano commented Apr 18, 2023

cfichtlscherer commented May 1, 2023

cfichtlscherer commented Jun 8, 2023

paulromano left a comment

paulromano Jun 19, 2023

paulromano Jun 19, 2023

paulromano left a comment

paulromano Jul 7, 2023

cfichtlscherer commented Jul 7, 2023

paulromano left a comment

		if (cell_born() == C_NONE)
		cell_born() = coord(n_coord() - 1).cell;

Pulse height tally for photons #2452

Pulse height tally for photons #2452

Conversation

cfichtlscherer commented Apr 3, 2023 • edited Loading

gridley commented Apr 3, 2023

paulromano left a comment

Choose a reason for hiding this comment

paulromano commented Apr 18, 2023

cfichtlscherer commented May 1, 2023

cfichtlscherer commented Jun 8, 2023

paulromano left a comment

Choose a reason for hiding this comment

paulromano Jun 19, 2023

Choose a reason for hiding this comment

paulromano Jun 19, 2023

Choose a reason for hiding this comment

paulromano left a comment

Choose a reason for hiding this comment

paulromano Jul 7, 2023

Choose a reason for hiding this comment

cfichtlscherer commented Jul 7, 2023

paulromano left a comment

Choose a reason for hiding this comment

cfichtlscherer commented Apr 3, 2023 •

edited

Loading