Disrupted/detached step unable to get final values

[maxbriel]

The following binary fails in step_detached due to failing in the get_star_final_values(secondary, secondary.htrack, m01) function on line 1929, see details below.

Error Traceback

Example binary

STAR1 = SingleStar(**{'mass': 9,
                      'state': 'H-rich_Core_H_burning'})
STAR2 = SingleStar(**{'mass':0.8,
                      'state': 'H-rich_Core_H_burning'})

BINARY = BinaryStar(STAR1, STAR2,  
                    **{'time': 0.0, 'state': 'detached', 'event': 'ZAMS', 'orbital_period':4513.150157, 'eccentricity': 0.0},
                    properties = sim_pop)

[maxbriel]

The issue originates from line 605 in interpolation.py.

         while (kvalue_low is None or np.isnan(kvalue_low)):
                mass_low = np.max(self.grid_mass[mass_low > self.grid_mass])
                try:
                    kvalue_low = self.grid_final_values[mass_low][key]
                except KeyError:
                    self.load_grid(mass_low)
                    kvalue_low = self.grid_final_values[mass_low][key]

kvalue_low = self.grid_final_values[mass_low][key] stays nan.

The key isS1_avg_c_in_c_core_at_He_depletion and it cannot find any masses.

@mkruckow could it be that these models do not reach He depletion?

[mkruckow]

For low mass stars it makes sense that there is no well defined value at He depletion, because they never get to He depletion. I think, we should stop the loop as soon as it reached the lowest mass in the grid and at this point keep the nan value.

[celiotine]

Looks like this bug is being caused by the same code described in Issue #286

[maxbriel]

I've tried exiting the loop when the mask becomes empty.

if np.sum(mass_low > self.grid_mass) == 0:
    break

This solves the crash, but results in weird binaries. The above example of $M2=0.77 M\odot$ should have a main-sequence lifetime longer than the age of the Universe.

The evolution of the binary with the loop fix does not indicate this. Specifically, the binary does not reach t_max_simulation_time.

The star is correctly matched to a nearby star in the single star grid, but the ODE solver in step_detached ends with s.t_events[2] as its output, indicating $t_{max}$ of the track (~1e9 yr)., but this is shorter than the the age of the Universe.

The 2 things to consider:

1. Why are these single star models returning a shorter than expected age?
2. How do we handle non-existent properties in the final_values of the system?

Matching details:

Matching attributes and their normalizations : ['mass', 'center_h1', 'log_R', 'he_core_mass'] [20.0, 1.0, 2.0, 10.0]
minimization in matching considered acceptable, with 0.00000004 < 0.01 tolerance
matching  H-rich_Core_H_burning  star with track of intial mass m0, at time t0: 0.799  [Msun], 38.170 [Myrs] 
 with m(t0), log10(R(t0), center_he(t0), surface_he4(t0), surface_h1(t0), he_core_mass(t0), center_c12(t0) = 
 0.799 -0.146 0.2719 0.2703 0.7155 0.000 0.0022
 The same values of the secondary at the end of the previous step was = 
 0.800 -0.147 0.2718 0.2703 0.7155 0.000 0.0019

[astroJeff]

Update (Aug 22): We think the issue may be related to the He core being "np.nan". @maxbriel can you look into this a bit more and find exactly what the problem is?

[maxbriel]

@mkruckow and I had a long discussion and debug session on this. We traced it down to the EEPS.

1. The matched PSyGrid track has a max_age of about ~10Gyr, but is said to reach central hydrogen exhaustion for low-mass stars, ie reaching termination flag: Single, low-mass star depleted hydrogen, terminating However, looking at the final entry in the PSyGrid track the star has not exhausted its hydrogen!

2. The raw MESA track continues evolving to ~20 Gyr and reaches central hydrogen exhaustion.

3. During the EEP creation the remaining evolution of the model is cut-off and not included. The PSyGrid object misses the rest of the evolution and steps after ~10 Gyr. It's unclear to me why this is cut off.

Next step: After the low-mass single star reruns, we check the EEPs for the low-mass stars to see if it's able to include the last evolutionary part.

[astroJeff]

Update (Sept 12): Needs to be rechecked now that @mkruckow has re-run the single star grids for infinite time. @maxbriel to check.

[maxbriel]

I've implemented the fix in 2 locations (for high and low-mass searches in the grid). This works with the above example and the models hit maxtime and are directed to END. I've made a PR.