TeamCOMPAS / COMPAS

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

Are we correctly aging stars on mass loss? #1134

Closed ilyamandel closed 1 month ago

ilyamandel commented 1 month ago

Hurley+ (2000), section 7.1, tells us that the fractional age on the MS should stay constant on mass loss. However, when I change the mass of a MS star, this does not seem to happen:

BaseStar *clone = Clone(OBJECT_PERSISTENCE::EPHEMERAL, false); double deltaMass = -m_Mass/1.0E5; clone->UpdateAttributesAndAgeOneTimestep(deltaMass, deltaMass, 0.0, true, false); std::cout<<"zetaEq:Tau"<<Tau()<<"new donor Tau"<Tau()<<std::endl; The old and new tau are not the same.

I tried updates such as clone->UpdateInitialMass(); // update initial mass (MS, HG & HeMS) clone->UpdateAgeAfterMassLoss(); // update age (MS, HG & HeMS) but these do not change tau.

jeffriley commented 1 month ago

This might just be an issue of looking at attributes during a timestep when things are still in a state of change. The age of the star (m_Age) is updated after mass loss, but tau (m_Tau) is not updated until the start of the next timestep.

I don't think tau is changed after massloss in Hurley's sse either - I think tau and age (of the star) are treated separately. Even in the paper, I think Hurley differentiates between t and tau. Given eq 11 I think tau probably should be updated when the age of the star is modified. As noted above, in COMPAS it is updated at the start of the next timestep when we call CalculateTauOnPhase(). In COMPAS SSE mode, and I think Hurley's sse, updating the age after massloss is the last thing done before starting the next timestep - so tau is not used until it is updated at the start of the next timestep. Maybe we need to look at that in BSE mode - I haven't checked (yet), but maybe we are using tau after aging the star and before tau is updated at the beginning of the next timestep.

I'll look in more detail a bit later

jeffriley commented 1 month ago

@ilyamandel if you need tau to be updated when we age the star after mass loss, we can change the code to do that - there should be no impact other than in any code you write to use tau (I don't think existing code uses tau before it is recalculated at the start of the next timestep).

ilyamandel commented 1 month ago
I added this to BaseStar::CalculateZetaEquilibrium() for testing.
Command line: 

./COMPAS -n 1 --mode BSE --initial-mass-1 2.1 --initial-mass-2 2 --semi-major-axis 0.05
/ILYA/

BaseStar *clone = Clone(OBJECT_PERSISTENCE::EPHEMERAL, false);
clone->UpdateAttributesAndAgeOneTimestep(deltaMass, deltaMass, 0.0, true, true);
std::cout<<std::fixed<<std::setprecision(15)<<"old mass"<<Mass()<<"new mass"<<clone->Mass()<<std::endl;
std::cout<<"old Tau"<<Tau()<<"new donor Tau"<<clone->Tau()<<std::endl;
std::cout<<"old Age"<<Age()<<"new Age"<<clone->Age()<<std::endl;
std::cout<<"old Time"<<Time()<<"new Time"<<clone->Time()<<std::endl;
std::cout<<"old type"<<(int)StellarType()<<"new type"<<(int)clone->StellarType()<<std::endl;

clone->CalculateTimescales();
//clone->UpdateInitialMass(); // update initial mass (MS, HG & HeMS)
clone->UpdateAgeAfterMassLoss(); // update age (MS, HG & HeMS)
std::cout<<"old Tau"<<Tau()<<"new donor Tau"<<clone->Tau()<<"new donor age/tMS"<<clone->Age()/clone->Timescale(TIMESCALE::tMS)<<std::endl;
std::cout<<"old Age"<<Age()<<"new Age"<<clone->Age()<<std::endl;
std::cout<<"old Time"<<Time()<<"new Time"<<clone->Time()<<std::endl;

clone->CalculateTimescales();
std::cout<<"old Tau"<<Tau()<<"new donor Tau"<<clone->Tau()<<"new donor age/tMS after recalc"<<clone->Age()/clone->Timescale(TIMESCALE::tMS)<<std::endl;
jeffriley commented 1 month ago

@ilyamandel I believe the answer to the eponymous question (if I can use that here...) is "Yes" (or maybe it would be better as "Probably, yes - in most cases...")).

The issue for the clone here is that in BaseStar::UpdateAttributesAndAgeOneTimestepPreamble() we call CalculateGBParams() and CalculateTimescales() - which we need to do because BaseStar::UpdateAttributesAndAgeOneTimestep() is an interface for the BSE code. Updating the timescales there is what's causing the issue described above.

The question is, should we recalculate timescales if dT = 0.0? I think we probably shouldn't, so we should probably change BaseStar::UpdateAttributesAndAgeOneTimestepPreamble() to only call CalculateTimescales() if p_DeltaTime > 0.0. Similarly, I wonder if we should only call CalculateGBParams() if one of the deltaMs is > 0.0?

I said "Probably, yes - in most cases..." above because there are a couple of places in the code we call BaseStar::UpdateAttributesAndAgeOneTimestep() with p_DeltaTime explicitly set to 0.0, so we should look at those - I think even in those cases it would be correct to not recalculate timescales. The calls to BaseStar::UpdateAttributesAndAgeOneTimestep() with p_DeltaTime = 0.0 are in BaseStar::ResolveMassLoss() and Star::UpdateAttributes().

If we are going to change BaseStar::UpdateAttributesAndAgeOneTimestepPreamble(), we should pass through the p_ForceRecalculate parameter from BaseStar::UpdateAttributesAndAgeOneTimestep() - or maybe just remove BaseStar::UpdateAttributesAndAgeOneTimestepPreamble() entirely and move the functionality into BaseStar::UpdateAttributesAndAgeOneTimestep().

Aside: the signature of BaseStar::UpdateAttributesAndAgeOneTimestep() in the description is wrong - it's missing p_ResolveEnvelopeLoss.

ilyamandel commented 1 month ago

@jeffriley :

I tested what happens on fast case A mass transfer and confirm that we are doing the right thing there -- the donor is aged in a way that conserves Tau.

The issue with my previous test was indeed that BaseStar::UpdateAttributesAndAgeOneTimestep() with p_DeltaTime set to 0.0 was already aging the star in a way that did not preserve Tau.

I am happy to close this issue since the main concern is really a non-issue, but will leave it to you to do so in case you'd like to keep it open as a reminder to address the changes you proposed in your previous message -- your call.

jeffriley commented 1 month ago

Changed BaseStar::UpdateAttributesAndAgeOneTimestepPreamble() so timescales are not recalculated when we know dT = 0 in v02.48.01. Not sure yet whether we can do a similar thing for GBParams and mass/mass0 - will look at that later. Closing this now.