TeamCOMPAS / COMPAS

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

Misleading detailed output: wrong (post-MT) stellar type stored during mass transfer episodes #774

Closed IsobelMarguarethe closed 1 year ago

IsobelMarguarethe commented 2 years ago

Describe the bug The mass transfer flag happens later than the mass transfer, leading to misleading detailed output.

Label the issue urgency_moderate - This is a moderately urgent issue

To Reproduce Steps to reproduce the behavior:

  1. Run COMPAS --random-seed 23046 --number-of-systems 1 --detailed-output
  2. Compare the pre- and post-mass transfer stellar type of the primary stored in COMPAS_Output.h5 to BSE_Detailed_Output_0.h5:
    
    import h5py

f1 = h5py.File("COMPAS_Output/Detailed_Output/BSE_Detailed_Output_0.h5", "r") f2 = h5py.File("COMPAS_Output/COMPAS_Output.h5", "r")

times_1 = f1['Time'][()] types_primary_1 = f1['Stellar_Type(1)'][()] primary_rl_ratio_1 = f1['Radius(1)|RL'][()] mt_tracking_1 = f1['MT_History'][()]

primary_rl_ratio_2 = f2['BSE_RLOF']['Radius(1)|RL<step'][()] types_primary_2_pre_MT = f2['BSE_RLOF']['Stellar_Type(1)<MT'][()] types_primary_2_post_MT = f2['BSE_RLOF']['Stellar_Type(1)>MT'][()] times_2 = f2['BSE_RLOF']['Time<MT'][()]

stored_mass_transfer_times_map = (mt_tracking_1 == 1)

print("Times of RLOF: {}".format(times_2)) print("Times at which MT_tracking = 1: {}".format(times_1[stored_mass_transfer_times_map])) print("Type of primary pre MT: {}".format(types_primary_2_pre_MT)) print("Type of primary post MT: {}".format(types_primary_2_post_MT)) print("Types of primary when MT_tracking = 1: {}".format(types_primary_1[stored_mass_transfer_times_map]))

Times of RLOF: [4.35815622 6.07168313]

Times at which MT_tracking = 1: [4.35815622 4.36934683]

Type of primary pre MT: [ 2 14]

Type of primary post MT: [ 7 14]

Types of primary when MT_tracking = 1: [7 7]



**Versioning:**
 - MacOS Mojave 10.14.6
 - COMPAS=v02.27.05
reinhold-willcox commented 2 years ago

Hey, sorry for the delay, I haven't kept up with git issues!

So I've run the system, and I can reproduce what you're seeing. There are a couple bugs here, but I think they're not exactly as you've described.

Here is the text output from the detailed plotter:

 Time(Myr),  Event,                                        M1(M_o),  type1,  M2(M_o),  type2,  a(R_o),  e
 0.000000    Zero-age main-sequence, metallicity Z=0.0142  51.940    1       39.768    1       196.638  0.000
 4.357383    Star 1:MS->HG                                 43.920    2       36.714    1       223.696  0.000
 4.358156    Stable mass transfer from 1 to 2              16.015    7       60.916    1       546.751  0.000
 4.358156    Star 1:HG->HeMS                               16.015    7       60.916    1       546.751  0.000
 4.909519    Star 1:HeMS->HeHG                             10.291    8       58.963    1       607.513  0.000
 4.937435    Star 1 undergoes supernova and forms a BH     6.314     14      58.807    1       612.265  0.000
 6.071683    Stellar Merger:BH+MS                          6.314     14      55.101    1       383.468  0.798

Your first output, # Times of RLOF: [4.35815622 6.07168313] is correct. There are 2 MT episodes in this evolution occuring at these times. The first is stable MT from star 1 to 2, where star 1 is a HG initially but becomes a HeMS. The second is (arguably unstable) from star 2->1 and results immediately in a merger. This is correctly flagged as a RLOF event (triggering the BSE_RLOF output), but does not trigger the MT_History tracker. This is likely a bug, though it was possibly done intentionally if whoever wrote it didn't care about tracking MT events that resulted in mergers.

Your second output, # Times at which MT_tracking = 1: [4.35815622 4.36934683] shows a separate bug (I think) because the tracker should be turned off by the second timestep. I thought this might be a case of outputting detailed printing multiple times in one timestep (which can lead to misleading lines like this), but then we shouldn't see two different times in this array, we should just see the same time twice. But this really means we're either triggering the tracker a second time for some reason, or not turning it off soon enough.

However, even if the second MT episode leading to merger did correctly trigger the MT_History tracker, you still would not see it in the code you ran, since you have masked specifically for MT_History == 1, which corresponds to "MASS TRANSFER STABLE STAR1 -> STAR2". To catch all tracked MT events, you should define the mask as MT_History != 0 or something like that. This is from the source code and is admittedly a bit hard to find, so probably we should make this parameter clearer in the documentation.

Your third and fourth output, # Type of primary pre MT: [ 2 14] and # Type of primary post MT: [ 7 14], are correct. The issues you've mentioned about misleading stellar types are because you've only looked at these parameters for star 1. In the first MT, star 1 is HG (2) to HeMS (7). In the second, star 1 is just a BH (14). My guess is you've been treating star 1 as the donor in all cases, but you need to also check RLOF(1)>MT to know if it is actually a donor (and same for star 2). Note that sometimes both stars are donors.

The fifth output, # Types of primary when MT_tracking = 1: [7 7] is correct after accounting for the previously noted bugs. It's just the type of star 1 at the end of the timestep(s) where MT_tracking is on. So the second 7 should not be there, but you would not expect to see anything from the second MT event unless you changed the (mt_tracking == 1) condition.

I hope that answered it, but let me know if anything needs to be clarified. I'll bring up the bugs at the next COMPAS dev meeting so we can iron them out.

IsobelMarguarethe commented 2 years ago

Hi @reinhold-willcox, thanks for this detailed response! Indeed the second RLOF time is not the one I care about here, I should have highlighted that (although I'm pleased that it revealed an additional bug). I did intend to mask specifically for MT_History == 1. The main bug I am concerned about here is that MT_History == 1 at [4.35815622 4.36934683] Myr while Stellar_Type(1) == [7 7] at those times. As I understand it, there should be no mass transfer from the primary once its stellar type changes from 2 to 7.

ilyamandel commented 1 year ago

@IsobelMarguarethe , @reinhold-willcox -- has this particular issue been fixed? Note that there are many DetailedOutput issues (e.g., binary orbit changes happen a step later than SNe which cause them; the RLOF data record a more massive ONeWD progenitor in the case of accretion induced collapse than DetailedOutput does; etc.), but we should probably deal with them one at a time.

reinhold-willcox commented 1 year ago

@ilyamandel Yes, the remaining TODO on this (the reason I didn't close earlier) was the documentation of the Mt_History parameter. But I changed this parameter a few months ago and updated the documentation, so I'll go ahead and close this issue.