-nan in GiantBranch::CalculateCoreMassAtBGB()

[jeffriley]

Running COMPAS with:

./compas -n1 --random-seed 0

results in -nan in GiantBranch::CalculateCoreMassAtBGB() in the following statement:

return std::min((0.95 * gbParams(McBAGB)), std::sqrt(std::sqrt(c + (MC_L_C1 * PPOW(p_Mass, MC_L_C2)))));      // sqrt is much faster than PPOW()

because the value of c is negative enough to cause the parameter to sqrt() to be negative

double GiantBranch::CalculateCoreMassAtBGB(const double p_Mass, const DBL_VECTOR &p_GBParams) {
#define gbParams(x) p_GBParams[static_cast<int>(GBP::x)]                // for convenience and readability - undefined at end of function
#define massCutoffs(x) m_MassCutoffs[static_cast<int>(MASS_CUTOFF::x)]  // for convenience and readability - undefined at end of function
   
    double luminosity = GiantBranch::CalculateLuminosityAtPhaseBase_Static(massCutoffs(MHeF), m_AnCoefficients);
    double Mc_MHeF    = BaseStar::CalculateCoreMassGivenLuminosity_Static(luminosity, p_GBParams);
    double c          = (Mc_MHeF * Mc_MHeF * Mc_MHeF * Mc_MHeF) - (MC_L_C1 * PPOW(massCutoffs(MHeF), MC_L_C2));         // pow() is slow - use multiplication

    return std::min((0.95 * gbParams(McBAGB)), std::sqrt(std::sqrt(c + (MC_L_C1 * PPOW(p_Mass, MC_L_C2)))));            // sqrt is much faster than PPOW()

#undef massCutoffs
#undef gbParams
}

To Reproduce
Run COMPAS with:

./compas -n1 --random-seed 0

Expected behavior
Production code does not cause -nan

**Versioning

OS: Ubuntu 22.04
COMPAS v03.07.01

EDIT: at first I thought this only occurs when p_Mass is 0.0, but that is not the case - the problem is evident for other values of p_Mass (e.g. 0.87)

[ilyamandel]

It says above Eq. (44) in Hurley+2000 that we should only apply this equation when the mass is greater than MHeF -- but we don't make that check in the code. We can, of course, put in such a check -- but there's no such thing as M_c,BGB for stars below the helium flash mass in the Hurley prescription. So, for such stars, we shouldn't be calling this function at all...

[ilyamandel]

I added a check to only compute McBGB for stars with mass above MHeF, following text above Eq. 44 in Hurley+, 2000 in #1298 .