Heavily fluctuating Radii in dev

[LiekeVanSon]

Describe the bug The radii of in particular massive single stars fluctuate multiple orders of magnitude in later evolution. This causes the maximum stellar radii of stars above 20 solar mass to be a mess and essentially keeps all stars compact at low metallicity. For the highest masses, the star even drops in radius below the size of it's eventual stripped He star

See example HR diagram of a 100 Msun star at different metallicities

Example high mass systems 90 Msun 100 Msun

Label the issue urgency_high severity_moderate

To Reproduce Steps to reproduce the behavior:

1. Run grid of "single" stars (i.e. non-interacting binaries). See grids in this zip file SingleStarGrids.zip
2. Plot either the HR diagram or the maximum radius of each star at every metallicity
3. Notebook I used to create attached gifs and plots is here: SingleStarInspectionPerMetallicity.ipynb.zip

Expected behavior This Issue was not present in v02.35. See expected behavior of max radii in plot below MaxRadius_perZ_V02_35.pdf

And an example of the behavior of a high mass star in the HR diagram at version 02.35

Screenshots Maximal radii

MaxRadius_perZ_v02_42.pdf

Old winds The issue seems not caused by the new wind prescription, the following plots were made with grids using --mass-loss-prescription 'BELCZYNSKI2010' Max radii old winds MaxRadius_perZ_v02_42_oldWinds.pdf

90 Msun old winds 100 Msun old winds

Versioning (please complete the following information):

• OS: Rocky Linux release 8.9
• COMPAS v02.41.06

[jeffriley]

Ok, I ran a few quick tests for 3 different COMPAS versions, thus:

v02.35.00:

./compas -n 1 --random-seed 53 --initial-mass-1 90 --initial-mass-2 0.1 --metallicity 0.002 --semi-major-axis 1E20 --logfile-type csv --detailed-output --mass-loss-prescription VINK

v02.41.06 & v02.42.01:

./compas -n 1 --random-seed 53 --initial-mass-1 90 --initial-mass-2 0.1 --metallicity 0.002 --semi-major-axis 1E20 --logfile-type csv --detailed-output --mass-loss-prescription BELCZYNSKI2010

and plotted the radius of the primary vs time from the detailed output files (I did not filter by record type, but that shouldn't be an issue):

The maximum radius seems to have dropped in later releases, but there's definitely something going on in the CHeB phase that doesn't seem to have been there in v02.35.00.

@ilyamandel if you think the problem described in issue 978 pre-dates v02.35.00, and this is something different, I can try to track down which release introduced this problem, then hopefully determine what the problem is.

[jeffriley]

And the detailed output files if anyone is interested:

v02.35.00.csv v02.41.06.csv v02.42.01.csv

[jeffriley]

This behaviour was introduced in v02.37.03. @ilyamandel that was PR https://github.com/TeamCOMPAS/COMPAS/pull/956 which was the fix for issue https://github.com/TeamCOMPAS/COMPAS/issues/855

I have done no more than identify the release the behaviour was introduced - I have not looked at the code to identify the code that changed the behaviour (I can begin that tomorrow if you like). I suspect some of the behaviour we are seeing (the reduction in the maximum radius) is probably correct - so the PR actually fixed a problem - but I'm not sure about the fluctuating radius ~CHeB - maybe that was an unintended (and unnoticed) ricochet...

v02.37.02.csv v02.37.03.csv

[LiekeVanSon]

Hi Jeff, thanks for tracking this down! I'm sorry I won't have time to dig into the code this week at all due to a workshop we're organizing here at CCA.

I do want to mention that I am highly suspicious of the max Radii as well. In particular if you look at the MaxRadius_perZ_v02_42.pdf plot, you see that (ignoring the obvious problem of the fluctuations) the code now predicts that stars above 30Msun never exceed 200 Rsun for metallicities below ~0.1 Zsolar. For sure these stars might experience LBV winds, but this seems extreme to me (though I'd love to hear others' opinions).

[jeffriley]

Hmmm. You might be right @LiekeVanSon - @ilyamandel, @SimonStevenson ?

I haven't checked yet, but I think the only real code change in PR https://github.com/TeamCOMPAS/COMPAS/pull/956 was to use Mas0 rather than Mass in some age/timescale calculations (which would affect radius etc.). I will look later today.

[ilyamandel]

Thank you for the beautiful illustrations, @LiekeVanSon , and for digging into the code, @jeffriley !

This is indeed the same issue as described in #978 , it was indeed exacerbated by PR #956 , but it predates it.

The evolution of stars in SSE (the Hurley+ code) is continuous in radii. However, the Hurley+ 2000 paper does not match the Hurley+ 2000 code. COMPAS attempted to accurately reproduce the paper, and PR #956 brought us into closer agreement with the paper. However, that only diverged us further from the code and exacerbated the problems.

Fixing the problems is not as trivial as undoing #956 , however. A complete fix requires actually digging through the SSE code and figuring out what needs to be changed in our implementation. I've done the digging through the code, but the required changes are sufficiently non-trivial that I haven't found the time to implement them in many months (my fault!).

@jeffriley kindly agreed to work on the code, and I'll send him my notes with my best recollection of what changes are required and what to refer to in the SSE code. We also agreed with @jeffriley to make comparisons of the HR diagram and time-radius plots for COMPAS vs Hurley code SSE stars (solar metallicity, no mass loss, masses of 1, 5, 20, 50, 100 Msun) to make sure we are really doing the right thing and have a reference going forward.

[jeffriley]

@LiekeVanSon see updates to issue #978 , and PR #1102

[ilyamandel]

Closing after fix in #1102