Default to high resolution random number generators

[mainegra]

There are situations were ranmar gives wrong results. It is a better approach to have someone switch to the low resolution ranmar rng intentionally, than having users obtaining wrong results unaware of this issue. Make ranlux (Mortran apps) or the high-resolution version (C++ apps) the default random number generator in EGSnrc. This is especially important as more simulations move to low energies and small scales.

[ftessier]

(see discussion also in #292)

[mainegra]

Hi-res vs low-res ranmar rng timing for isotropic photon source:

• Kerma in 5 cm radius air sphere @ 1m 40 keV (egs_kerma):
  • ph transport only: 8% slowdown
• DD in 5 cm radius water sphere for 5 keV (egs_app)
  • ph only: 4% slowdown
  • ph+e- : 8% slowdown
• DD in 15 cm radius water sphere for 200 keV (egs_app)
  • ph+e- : 7.24% slowdown
• DD in 15 cm radius water sphere for 2 MeV (egs_app)
  • ph+e- : 7.5% slowdown
• DD in 15 cm radius lead sphere for 2 MeV (egs_app)
  • ph+e- : 4% slowdown!!!
• DD in 15 cm radius lead sphere for 20 MeV (egs_app)
  • ph+e- : 4% slowdown!!!

In low Z materials the slowdown is about 8% from keV energies up to 2 MeV. While in high Z materials the penalty is only 4%.

ranlux vs low-res ranmar rng timing for isotropic photon source (Mortran):

• DD in 5 cm radius water sphere for 5 keV (EDKnrc)
  • ph only: 14.4% slowdown
  • ph+e- : 20.5% slowdown

@dworogers it would be interesting to check what happens for Linac simulations!

[MartinMartinov]

I don't still have the hard numbers, but in my gold nanoparticle simulations when I first started using the high resolution option years ago, I found slowdowns ranging 10-14%

[ftessier]

Hmm, that's not negligible. I wonder if we make the low-res version the default, with a prominent notice to turn the hi-res option when "needed" (with a clearer definition), or default to hi-res with a prominent notice stating that one could gain ~10% with the low-res option, at your own risk! The latter is the safer. Will look at ranlux timing, and of course the real answer is xoshiro256+ 😉

[mainegra]

@MartinMartinov thanks for sharing that info! I'm running a case for a Pb sphere and will update above. @ftessier a proper definition on when to use the hires rng might be harder to find. Perhpaps example cases such as those encountered by Martin, Rowan and Pat. And now the case of low energy photons I have encountered. Seems like we could use a co-op student to work on this in the Summer .... 😉

[mainegra]

@ftessier if you are interested here is a bit of a soap opera I discover last night in the world of RNGs researchers:

https://www.pcg-random.org/posts/on-vignas-pcg-critique.html

Spoiler alert: They might be arguing about minutiae in our MC world!

[ftessier]

Wow, thanks for that perspective; perhaps pcg64 is the better contender. Will have to do some more reading on this entertaining drama! Then again, if the code is in the public domain, we can include "all of them" and let the users choose... In terms of hires/lowres, I'd side with making hires the default with some clear messaging about the potential ~10% saving because indeed, it is hard to define broadly in which cases hires is needed.

[rtownson]

I agree with highres default with instructions to switch to lowres when the time savings are needed. I think the general philosophy in EGSnrc is that the default settings are accurate in the vast majority of conditions.

[mainegra]

@ftessier see above the result ranlux vs ranmar in Mortran. It is actually getting worse!

[ftessier]

Darn! Skip all this and implement xoshiro256+ right away? 🤓

[ftessier]

Posted in the thread again. Here are my replies regarding the pcg vs xoshiro debate:

Interesting and entertaining feud between the xoshiro and pcg generator authors! From there I was led to a talk by pcg author Melissa O'Neill with an great historical perspective and a simple explanation of what pcg does.

[ftessier]

On pcg vs xoshiro: after additional reading, I side with the xoshiro generator, for its dead simple code, its speed, and its ability to jump ahead to generate independent random sequences for parallel simulations. In terms of randomness quality, I defer to the mathematicians; to me it suffices that it passes the TestU01 BigCrush test. Note that neither the pcg paper nor the xoshiro papers have been peer-reviewed (yet), if that matters.

[ftessier]

Only thing left as an implementation detail is how to sample doubles uniformly on [0,1). Do we just take the high bits from the xoshiro256+ 64-bit integer result and plug them into the mantissa of a double? Is that correct? Other suggestions welcome.

[ftessier]

Using the 53 most significant bits as the mantissa is sufficient for our purpose in EGSnrc, although strictly speaking it is not uniform: values in (0, 2e-53) are never picked. I think we can live with that!

double random_real_53(void) {
	return ldexp((double)(random64() & ((1ULL << 53) - 1)), -53);
}

This is from http://mumble.net/~campbell/tmp/random_real.c, which also proposes a neat solution to the the non-uniformity issue.

[mainegra]

@ftessier looking forward to your xoshiro256+ implementation! 🤓

[mxxo]

Perhaps it's worth mentioning ranlux++ as another potential choice? Greatly improved performance over plain ranlux but still "suitable for use in the most demanding Monte Carlo calculations". There exists a portable C++ implementation that was recently added to ROOT: https://arxiv.org/abs/2106.02504

[ftessier]

I have not looked at the ranlux prng recently, because they are (were) painfully slow generators. But an implementation is useful to eventually check other generators. Thanks for the reference, which is discussing some interesting optimizations, so I will look into this again. I quickly note in Figure 2 performance above about 10 ns/number (with conversion), whereas xoshiro is sub 1 ns/number according to https://prng.di.unimi.it.