Regression test stdev discrepancy when DD is active

CEMeNT-PSAAP / MCDC

MC/DC: Monte Carlo Dynamic Code

https://mcdc.readthedocs.io/en/latest/

BSD 3-Clause "New" or "Revised" License

25 stars 24 forks source link

Regression test stdev discrepancy when DD is active #231

Open alexandermote opened 3 months ago

alexandermote commented 3 months ago

When developing the dd_cooper 3D domain decomposition test for #230, I noticed that the regression test succeeded for mean flux, but failed for standard deviation. The values it fails by are generally very small, under 10^-16, but it's off by enough to trip the failure condition in run.py. It's entirely likely that this is some problem with the input I developed, but I was asked to create an issue to track it just to be safe.

My best guess was that it had something to do with implicit capture, which is turned on for this problem. I turned it off and it passes at 50 particles, but fails at 1000 in a similar fashion to when implicit capture is turned on.

alexandermote commented 3 months ago

Adding some details that we tracked down in the comments for #230:

The issue is not unique to 3D domain decomposition problems. Running dd_slab_reed without domain decomposition gives a different standard deviation value than we get when running with domain decomposition.
The issue is not caused by the tally recombination introduced in #230. Running dd_slab_reed without this recombination algorithm (which it does not need, being a 1D problem) gives us the same issue.
The issue is worsened when more particles are added. With a low enough number of particles, the stdev values will align. This paired with the number of particles in a tally cell being used to calculate stdev values in tally_closeout() makes me think the number of particles in a cell may be changing slightly when domain decomposition is active.

alexandermote commented 2 weeks ago

Two updates to this issue that I've noticed while writing my M&C paper:

In cases where mean flux is accurate and stdev isn't, I noticed that the stdev values were off by a factor of exactly the number of processors in the subdomain. Multiplying this stdev value by the number of processors seemed to fix this issue.
In cases with 3D domain decomposition, the stdev values seemed to be off in areas where a domain boundary existed with no other geometry boundary. Changing the geometry mesh to align with the DD mesh almost completely removed these issues, but didn't solve the problem outright.