Computing Long-Duration Gamma-Ray Bursts from Double Compact Object Populations 🌌

If you haven’t done so yet, export the path POSYDON environment variables. For example:

%env PATH_TO_POSYDON=/Users/simone/Google Drive/github/POSYDON-public/
%env PATH_TO_POSYDON_DATA=/Volumes/T7/

Creating the GRB TransientPopulation

We will use the same BBH population as in the Analyzing Merging BBH Populations: Rates & Observations tutorial to compute the long gamma-ray burstss (LGRB) rate associated with the formation of merging BBHs.

Make sure the BBH_contact.h5 file is present to continue with this tutorial.

from posydon.popsyn.synthetic_population import Population

file = 'BBH_contact.h5'

BBH_population = Population(file)
BBH_population.mass_per_metallicity

BBH_population.calculate_formation_channels(mt_history=True)

We had already selected only binaries that result in a merger in the BBH_contact.h5 file, so we don’t need to perform any selection on the binaries.

Thus, our next step is to create a selection function for LGRBs based on the properties in the history and onlines of the binaries.

We have already included a basic selection function in posydon.popsyn.transient_select_func, which select if a GRB occurs based on the presence of m_disk_radiated for either star. It will also output some pre and post supernova properties.

As you can see below, the GRB_selection requires an additional input parameter of S1_S2. The Population.create_transient_population function cannot pass these additional arguments to sub-function. Thus, we will be required to wrap the GRB_selection functions.

This is useful, if you would like to create similar population, but use different model parameters.

from posydon.popsyn.transient_select_funcs import GRB_selection

GRB_selection(BBH_population.history[10], BBH_population.oneline[10], BBH_population.formation_channels.loc[[10]], S1_S2='S2')

import pandas as pd

def GRB_wrapper(transient_chunk, oneline_chunk, formation_channels_chunk):
    """calculate the GRBs from star 1 and star 2 for a given chunk of the population"""
    df_1 = GRB_selection(transient_chunk, oneline_chunk, formation_channels_chunk, S1_S2='S1')
    df_2 = GRB_selection(transient_chunk, oneline_chunk, formation_channels_chunk, S1_S2='S2')

    # combine the two dataframes
    if df_1 is not None and df_2 is not None:
        GRB_df = pd.concat([df_1, df_2])
    elif df_1 is not None:
        GRB_df = df_1
    elif df_2 is not None:
        GRB_df = df_2
    else:
        return None

    return GRB_df

LGRB = BBH_population.create_transient_population(GRB_wrapper, 'LGRB')

You’ve now created a LGRB population, where either the first and/or second star has some amount of radiative disk.

This is, of course, just an example and the actual amount of disk required to power a LGRB will be higher than some of the values included in our “LGRB” rate.

Moreover, the number of binaries in this population is insufficient to actually see the all unique events.

LGRB.get_efficiency_over_metallicity()

LGRB.plot_efficiency_over_metallicity(channels=True)

After this selection, we follow similar steps to the BBH analysis for creating and plotting a cosmic star formation history weighted rate.

If you’ve followed the previous tutorial on the BBH analaysis, you will still have access to those populations too.

LGRB_rates = LGRB.calculate_cosmic_weights('IllustrisTNG')

LGRB_rates.calculate_intrinsic_rate_density(mt_channels=True)

from posydon.popsyn.synthetic_population import Rates
BBH_rates  = Rates('BBH_contact.h5', 'BBH', 'IllustrisTNG')

BBH_rates.intrinsic_rate_density

Below we plot the BBH rate and the un-normalised LGRB rate together, showing how you can access multiple event rates, and transient populations from the same file.

import matplotlib.pyplot as plt
plt.plot(BBH_rates.centers_redshift_bins, BBH_rates.intrinsic_rate_density['total'], label='BBH')
plt.plot(LGRB_rates.centers_redshift_bins, LGRB_rates.intrinsic_rate_density['total'], label='LGRB')

plt.yscale("log")
plt.xlim(0,10)
plt.legend()
plt.show()

LGRB_rates.population

fig, ax = plt.subplots(1,1)
LGRB_rates.plot_hist_properties('S1_spin_postSN', intrinsic=True, label='S1', color='blue', ax=ax, show=False)
LGRB_rates.plot_hist_properties('S2_spin_postSN', intrinsic=True, label='S2', color='orange', ax=ax, show=False)
ax.legend()
ax.set_xlabel('BH spin post LGRB')
ax.set_ylabel('\# Events in bin')
plt.show()

It’s possible to calculate more properties for the LGRBs using the posydon.popsyn.GRB module.

It contains the get_GRB_properties that can provide a more detailed GRB properties based on empirical or simulation relations for the LGRB energy and beaming factor.

Enjoy exploring this function and creating a more detailed GRB population! Keep in mind that it’s not possible to add additional columns to the transient population without overwritting the population.

The next tutorial will focus on rerunning specific binaries in a population or setting up an unique binary yourself.