openmc-dev / openmc

OpenMC Monte Carlo Code
https://docs.openmc.org
Other
771 stars 496 forks source link

Tally filter for time #1932

Closed cfichtlscherer closed 2 years ago

cfichtlscherer commented 2 years ago

For different applications (e.g., the simulation of neutron coincident measurement), it could be helpful to have a global time variable in OpenMC. For example, this could tell you after which period of time neutrons got absorbed (from an observer ). To the best of my knowledge, OpenMC is at the moment not capable of this functionality.

I added a TimeFilter which will return the time the particle spent in the simulation since its birth. The idea is to add a time attribute to the particle class. This attribute stores the time since when the particle is alive (from the observer perspective). It is updated in the Particle::event_advance() class: whenever the neutron advances, we calculate the speed from its energy (relativistic) and then check how long the particle took for the distance. This period is then added to the time attribute of the particle.

Is this the right way to do this, or am I missing something here?

This will be much more difficult for photons, since the speed of light changes with the refractive index of that material. I am not really sure how to proceed here. Maybe it could be a good idea to add an attribute for the refractive index in the material class? Probably having a database would be good here as well.

Implementation vise, the TimeFilter is very similar to the EnergyFilter.

I have the feeling it is too early for a pull request. But I wanted to check if perhaps somebody is already working on this? Has an idea on the photons or some general comments; if I am making a mistake here? Or would it be better to implement it differently?

Of course I am always ready to share the code :-) if someone is more interested in the details.

Thanks a lot already!

gridley commented 2 years ago

Have you seen this?

https://github.com/paulromano/openmc/tree/time-filter

You make a really interesting point about photons, but I think for most of the applications of OpenMC, we should probably default to using the speed of light in vacuum unless otherwise specified. What is an application where time-dependent photon transport matters?

Aside from that, I think your intuitition is correct about how to implement this. That is pretty much what Paul did, although IIRC he did not use the relativistic velocity for neutrons. Actually, we should definitely use the relativistic model since time-of-flight experiments with 14.1 MeV neutrons (beta = 0.17) are pretty commonplace in the security world, from what I understand. So relativity should def be taken into account...

Do you know of any open validation tests we could do for transient stuff like this? I know it would be applicable to time-of-flight experiments, but am unaware of any validation-grade experimental data on this. It would be tricky since we'd have to add detector response stuff, which means pulse height tallies, which are not yet in the main code base. IMO we shouldn't merge pretty advanced physics features like this without validating against a few things..

I think once this is merged into the main codebase, we would want to build Sjenitzer's dynamic MC method on top of it. You'd have some timestepper that takes a reference to a time filter, and you do this forced precursor decay process inside each time bin...

paulromano commented 2 years ago

@cfichtlscherer As Gavin pointed out, I have a branch in progress with exactly what you've described. It works pretty well so far and allows you to apply the time filter to volumetric tallies. I actually do use the relativistic formulas so any time dilation effect on high-energy neutrons is accounted for. The point about refractive index is indeed interesting -- to my knowledge other codes do not account for this effect. I really don't know much about optics, but perhaps the refractive index is very close to 1 for high energy photons and hence can be ignored?

In any event, I'm hoping to submit a PR for my time-filter branch soon. One issue I need to figure out is how to apply it to surface tallies because it currently doesn't work correctly in that case.

ojschumann commented 2 years ago

@cfichtlscherer In addition to the Particle class, you need to update the MPI-code that sends the fission bank. That is quite straightforward and documented somewhere in the code. I also implemented a time-variable for OpenMC in some internal code, so I might share that as well. I use this feature e.g. for a time-stepping simulation to calculate transients, Rossi-alpha (I'm aware of #1596 and hope this one will be merged in the future) or probability of initiation.

Some remarks:

I would be very interested in this topic and could help with coding or benchmarking.

gridley commented 2 years ago

Just to clarify: I had indeed remembered Paul's code incorrectly. It did use relativistic velocities the whole time.

cfichtlscherer commented 2 years ago

Thanks a lot for your answers! And no, I haven't seen Pauls's branch. Next time I will check better what's already there. I searched for TimeFilter and couldn't find it. That's why I started working on it.

I am, at the moment, only interested in neutrons. But maybe this could become a unique feature of OpenMC - probably we should check first if this is needed and cannot be neglected in a specific application. Implementation should not be too complicated if the refraction index is a material attribute.

@gridley For the validation, I propose running an example in MCNP, Geant4, and OpenMC. Unfortunately, I don't have an MCNP license at the moment, so I could only do the OpenMC and Geant4 models.

We could otherwise maybe simplify the reactions (the entire physics of the simulation) and see if we can think of an analytical validation model (a la [The Verification of Monte Carlo Codes in Middle Earth, 1994]). But I think this kind of validation is only meaningful to a limited extent.

I agree with @ojschumann experimental validation will be very complicated as all possible sources of error must be addressed.

@ojschumann Yes! I also did the MPI-code part. You're right.

Having a benchmark for calculating multiplicity counter results would be very interesting. Probably I should finish pull request #1881 first.

@ojschumann could you maybe send me the link to the ESARDA workshop paper? I would be interested.

I think https://tuprints.ulb.tu-darmstadt.de/5621/1/thesis-print-defended.pdf provides a good overview of the needed data (@ohnemax is my supervisor, so I might be biased here ;-) ).

I will spend more time on this feature in 2022.

I wish you a Merry Christmas and a Happy New Year.

ojschumann commented 2 years ago

@cfichtlscherer You can find my paper here.

paulromano commented 2 years ago

@cfichtlscherer I didn't publicize my branch very well so it's definitely not your fault for missing it! I've been thinking that it may be beneficial to have periodic videocalls for people in the development community to raise better awareness of what everyone is working on and so that we're not duplicating each other's work. Perhaps in the new year we can start organizing something like that. Would you all be interested in that?

In any event, this issue motivated me to stop dragging my feet on this feature -- just opened a pull request (#1934)

cfichtlscherer commented 2 years ago

Yes, I would be very interested in periodic video calls 😀

py1sl commented 2 years ago

In terms of validation with experimental data, something like the pulsed neutron die away/ differential die away experiments or the Livermore pulsed spheres might be an option. I believe there are both fissile and non fissile experiments. I am more interested in the pulsed neutron die away, but I will look into others as well

ilhamv commented 2 years ago

Hi everyone,

I was planning to add such time-dependent features too! Thank you @ojschumann for nudging my PR so that I'm aware of this ongoing development.

Besides TimeFilter, we may want to consider adding TimeEdgeFilter as well, where we use a "time-edge estimator" that scores the product of particle's speed and weight whenever a time edge is crossed. I think this time-edge filter may be useful whenever time-edge quantities are desired instead of the time-average ones. Note that such TimeEdgeFilter would not work in combination with a surface tally (unless there is a clever way to force surface and time-edge crossings at the same time?).

As for verification of the time-dependent features, we can use Ganapol's semi-analytical benchmark AZURV1 (https://inis.iaea.org/search/search.aspx?orig_q=RN:41070601). With this, we can verify three scalar flux quantities: (1) time-average spatial-average, (2) time-average spatial-edge, and (3) time-edge spatial-average. Attached is a time-edge spatial-average solution of the slab problem with scattering ratio set to 1.1 (supercritical) and run with 1000 histories; we can also perform a convergence test by increasing the number of histories and observe the N^-0.5 behavior of the errors to be even more confident (also attached).

paulromano commented 2 years ago

The original intent of the this issue was addressed in #1934, so I'm going to go ahead and close this. If there are further feature requests (e.g., the TimeEdgeFilter discussed above), those should be made into separate issues.