njoy / NJOY2016

Nuclear data processing with legacy NJOY
https://www.njoy21.io/NJOY2016
Other
98 stars 86 forks source link

Fix bug in PURR related to heating values #18

Closed paulromano closed 7 years ago

paulromano commented 7 years ago

This pull request fixes a subtle bug in purr when generating heating values in the probability table when ihave=2. If a user asks for partial heating values for MT=302 and MT=402 but not MT=318 (which seems reasonable if the isotope in question is not fissionable), this results in ihave=2. purr will then try to use the partial heating values to account for fluctuations in the heating values when generating ptables. However, in rdheat, the fact that MT=318 is missing means that the heating values for it are never set. Furthermore, the heat array that is used in purr is never initialized, so the uninitialized values are then used at the end of purr when calculating heating values. Namely, heat(3,ie,it) will be uninitialized and may possibly be NaN which eventually causes a floating point exception somewhere down the line. For me, this always seemed to happen here.

So far in my testing, this seems to have resolved the random NaNs causing problems for me when building with gfortran (as described in #16). However, I don't know if it will do the same for @brown170.

Semi-related question -- I see that the CMake file for NJOY now (supposedly) requires gfortran 5.1+. Is there a particular reason for this? I say supposedly because it doesn't appear that CMAKE_Fortran_COMPILER_VERSION actually gets set (at least on my machine). I have similar logic in the CMakeLists.txt file for OpenMC that relies on execute_process instead.

jlconlin commented 7 years ago

@paulromano Thanks for submitting this PR. We have been unable to duplicate the NaN issue that @brown170 has been seeing. Do you have an input deck that we can have that reliably reproduces the NaN issue?

paulromano commented 7 years ago

Well, unfortunately there's no such thing as "reliable" uninitialized values -- totally depends on the compiler, optimization, machine, state of memory, etc. I can share my input and my setup though and you're welcome to try to reproduce. On my laptop, which is a Lenovo X1 Carbon running Ubuntu 16.10 (with gfortran 6.2.0), I ran into this bug when processing Fe-58 from ENDF/B-VIII.beta4. My input was as follows:

reconr /
20 21
'ENDF/B-8 PENDF for  26-Fe- 58 '/
2637 2/
0.001 0.0 0.003 /
'ENDF/B-8:  26-Fe- 58 '/
'Processed by NJOY'/
0/
broadr /
20 21 22
2637 1 0 0 0. /
0.001 1.0e6 0.003 /
293.6
0/
heatr /
20 22 23 /
2637 2 /
302 402 /
purr /
20 23 24
2637 1 1 20 64 /
293.6
1.e10
0/
acer /
20 24 0 25 26
1 0 1 .01 /
'ENDF/B-8:  26-Fe- 58  at 293.6'/
2637 293.6
1 1/
/
stop

I only got NaNs and floating point exceptions when I built with CMAKE_BUILD_TYPE=Release.

jlconlin commented 7 years ago

@paulromano I had problems with that input deck too in ACER. I compiled with gfortran 6.3.0_1 via Homebrew. I used your fix and the problem went away.