search

TeamCOMPAS / COMPAS

COMPAS rapid binary population synthesis code
http://compas.science
MIT License
64 stars 66 forks source link

NaN Values for 'Mass_Core@CO(SN)' with Latest COMPAS (version 03.07.00) #1245

Open lee0517-hub opened 6 days ago

lee0517-hub commented 6 days ago

Describe the bug When using the latest version of COMPAS (version 03.07.00), I encountered an issue where the data for 'Mass_Core@CO(SN)' is NaN for some binaries.

Label the issue urgency_high

To Reproduce Steps to reproduce the behavior: ./COMPAS --critical-mass-ratio-prescription NONE --kick-magnitude-distribution ZERO --random-seed 1724542090 --common-envelope-alpha 1.6265542767303343 --wolf-rayet-multiplier 0.8084912681292984 --number-of-systems 1 --initial-mass-function KROUPA --mass-ratio-distribution DUQUENNOYMAYOR1991 --semi-major-axis-distribution SANA2012 --minimum-secondary-mass 5 --metallicity-distribution LOGUNIFORM --tides-prescription PERFECT --remnant-mass-prescription FRYER2022 --common-envelope-allow-radiative-envelope-survive False --common-envelope-allow-immediate-RLOF-post-CE-survive False

Expected behavior A clear and concise description of what you expected to happen.

Screenshots 2024-10-18

Versioning (please complete the following information):

Additional context Add any other context about the problem here.

jeffriley commented 1 day ago

The immediate cause of this problem is that in HeGB::CalculateCoreMassOnPhase_Static(), p_time is greater than tx, but it is also greater than tinf2, so we end up taking the sqrt of a negative number, hence the NaN.

The calculation is given by Hurley et al., 2000 eq 39, modified per the discussion following eq 84. I believe our implementation of the equation is correct, so we need to determine why p_Time is > tinf2. Hurley doesn't mention caveats, so I assume he expected tinf2 would always be greater than p_Time (t in eq 39). Preventing the NaN naively is easy, but I'm not sure that would be physical... Needs more investigation.

lee0517-hub commented 1 day ago

It seems that this bug has only recently appeared. I am not certain from which version it started, but I can confirm that the issue does not exist in version v03.00.05. I hope this information helps and saves you some time in the debugging process.

jeffriley commented 23 hours ago

Problem was introduced in v03.05.02.

@ilyamandel

Prior to v03.05.02, running the command @lee0517-hub gives above produces a single entry in the BSE_Supernovae file (no NaNs).

Running the same command with v03.05.02 results in a second SN event, with both Mass_Core@CO(SN) and Mass_CO_Core@CO(SN) set to Nan.

I ran the same command with v02.50.01, and there was only a single entry in the BSE_Supernovae file, so I don't believe this problem was introduced by reinstating IsCCSN() etc. to pre-v03.00.00 states.

I removed the call to EvolveOnPhase(0.0) from the Initialise() function of all evolved stellar types in v03.05.02 and re-ran the command shown above, and there was only a single entry in the BSE_Supernovae file, so it would seem the addition of the call to EvolveOnPhase(0.0) in Initialise() for all evolved stellar types has introduced the problem (or maybe just exposed it, but I think we need to determine if the call toEvolveOnPhase(0.0) in Initialise() is producing valid results - possibly making attribute values a step ahead of where they should be?).

EDIT: changelog.h indicates that you added a call to EvolveBinary(0.0) - it was actually EvolveOnPhase(0.0) - we should probably fix changelog.h.

EDIT: @ilyamandel should we be updating timescales and gbparams after calling EvolveOnPhase(0.0) in Initialise()?