Your First Binary Simulations with POSYDON 🌠

Tutorial goal:

In this tutorial, we will run a small population of 10 binaries locally, and explore how to manipulate the output data from your population.

New concepts:

Population ini file

If you haven’t done so yet, export the path POSYDON environment variables. Set these parameters in your .bash_profile or .zshrc if you use POSYDON regularly.

[ ]:

%env PATH_TO_POSYDON=/YOUR/POSYDON/PATH/
%env PATH_TO_POSYDON_DATA=/YOUR/POSYDON_DATA/PATH/

Creating the Initialisation File

To run population synthesis with POSYDON, a population_params.ini file is required. This file described how the stellar population is created and what prescriptions and parameters are implemented in specific steps.

POSYDON comes with a default population_params_default.ini file found at PATH_TO_POSYDON/posydon/popsyn or have a look here.

The file is split in three main parts. You can find more details about their properties by clicking on their links. 1. `SimulationProperties <>`__: - these describe the properties and parameters of different steps in the evolution of a binary systems. 2. `BinaryPopulation <>`__: - parameters of the initial sampling of the binary population, such as initial mass function, period distribution, and metallicity. - Also contains parameters on how the population is ran, such as how many binaries are kept in memory. 3. `SavingOutput <>`__ - Describes the data from the binary and each individual star to the output files.

We will copy the default population run parameter file to the current folder.

[ ]:

import os
import shutil
from posydon.config import PATH_TO_POSYDON

path_to_params = os.path.join(PATH_TO_POSYDON, "posydon/popsyn/population_params_default.ini")
shutil.copyfile(path_to_params, './population_params.ini')

Reprocessed POSYDON v1 data

If you’re using the reprocessed POSYDON v1 dataset, you will only have solar metallicity available! This means you will only be able to run populations with metallicity = [1]!

You should be able to follow along with this tutorial and can follow along with the “One metallicity notebook”, but the multi-metallicity component is not yet available.

Creating and running a binary population

The copied population_params.ini contains the parameters to run 10 binaries at a metallicity of \(Z=1 Z_\odot\).

If you open the file and scroll down to the BinaryPopulation section, you will see how they’re defined:

metallicity = [1] # [2., 1., 0.45, 0.2, 0.1, 0.01, 0.001, 0.0001]
# In units of solar metallicity
...
number_of_binaries = 10
# int

If you like to run a small population in a notebook, you can use the PopulationRunner to do this. If you want to run a specific binary instead, have a look at the `Binary Tutorial <>`__.

The PopulationRunner class takes the .ini file and sets-up a population run.

This will create `BinaryPopulations <>`__ for each metallicity defined in the .ini file. In this case, we can check and see that a single BinaryPopulation is created and contains 10 binaries.

[ ]:

from posydon.popsyn.synthetic_population import PopulationRunner
poprun = PopulationRunner('./population_params.ini', verbose=True)

[ ]:

print('Number of binary populations:',len(poprun.binary_populations))
print('Metallicity:', poprun.binary_populations[0].metallicity)
print('Number of binaries:', poprun.binary_populations[0].number_of_binaries)

For this tutorial, we set verbose=True, which shows you the progress of the population run. This overwrites the population verbose set inside the population_params.ini file.

Now we are ready to evolve the binary population. This should take about 30 seconds, but depends on your machine.

[ ]:

poprun.evolve()

Inspecting the population: Population class

When you ran the population, you might have seen that a temporary folder with the name 1e+00_Zsun_batches was created while the binaries were evolved. This is a temporary folder in which populations are temporarly saved. After the binary evolution has finished, the binaries in the folder are moved to a single file named 1e+00_Zsun_popululation.h5. This is done automatically, when you run a population using the PopulationRunner class.

The created file contains 3 main components:

history: the evolution of an individual binary in a pandas DataFrame
oneline: a single line to describe the initial and final conditions and some one-of parameters, such as the metallicity.
mass_per_metallicity: some metadata on the population, such as the total simulated mass, the actual underlying mass of the population, and the number of binaries in the file.

The Population provides an interface with these components in the file, such that you’re able to share the populations runs and can work with large population that do not fit in memory.

Older Population Files

If you’re using older population files, you can make them compatible with the Population class by calling Population(pop_file, metallicity,ini_file), where the metallicity is in solar units. You will only need to do this once; afterwards you can initialise the class like normal.

[ ]:

from posydon.popsyn.synthetic_population import Population

[ ]:

pop = Population('1e+00_Zsun_population.h5', verbose=True)

The pop.mas_per_met shows some basic information about the population you’ve just created.

The index (metallicity) is the metallicity of your population in solar units.
simulated mass is the total ZAMS mass that has been evolved in the population.
underlying mass is the actual mass of the population if one integrates the IMF and period distribution fully.
number_of_systems shows the 10 systems in the file.

[ ]:

pop.mass_per_metallicity

There are some additional metadata properties available, such as: - metallicities: the metallicity in absolute metallicity - solar_metallicities: the metallicities in the file in solar metallicity - number_of_systems: the total number of systems in the Population file - indices: the indices of the binaries in the file - columns: the columns available in the history and oneline dataframes

[ ]:

# you can also access the total number of systems in the file with
print(pop.number_of_systems)

Population.history

[ ]:

pop.history.lengths

pop.history loads in the full history of all the binaries into memory. You can access individual or a selections of the population using several methods:

pop.history[5]
pop.history[[0,4]]
pop.history[‘time]
pop.history.select()

The select function is the most powerfull way to access the binaries, because it allows you to perform selections based on the specific columns available in the history dataframe. For example, below we can select on state == 'RLO1', which gives us all the rows with RLO1 occuring.

The available identifiers are limited to string columns (state, event, step_names, S1_state, S2_state), index, and columns names.

[ ]:

# select only binary_index 5
pop.history[5]

[ ]:

# select binary 0 and 4
pop.history[[0,4]]

[ ]:

pop.history['time'].head()

You can also check what columns are available in the history file:

[ ]:

pop.history.columns

[ ]:

# using the select function
pop.history.select(where='index == 9')

[ ]:

# selecting all RLO1 states and only time and state columns
pop.history.select(where='state == RLO2', columns=['time', 'state'])

[ ]:

# selecting rows 10 to 16
pop.history.select(start=10, stop=16)

You might have notices while using the above functions that not all the binaries will have the same length in the history. You can access these with pop.history_lengths. This information is also stored in the population file.

[ ]:

pop.history_lengths

Population.oneline

Population.oneline provides a similar interface to accessing the DataFrame in the population file as Population.history, with similar functionality being available.

The select function only has acces to: - index - column names - string columns: state_i, state_f, event_i, event_f, step_names_i, step_names_f, S1_state_i, S1_state_f, S2_state_i, S2_state_f, S1_SN_type, S2_SN_type, interp_class_HMS_HMS, interp_class_CO_HeMS, interp_class_CO_HMS_RLO, interp_class_CO_HeMS_RLO, mt_history_HMS_HMS, mt_history_CO_HeMS, mt_history_CO_HMS_RLO, mt_history_CO_HeMS_RLO

[ ]:

pop.oneline[5]
pop.oneline[[0,4]]
pop.oneline.select(where='index == 9')
pop.oneline.select(where='index == [0,9]')

[ ]:

pop.oneline.columns

Population.formation_channels

You might be interested in figuring out what sort of formation pathways/channels a binary has followed through its evolution.

This is not a standard output of the population synthesis, but you can include it into the population file by calculating it. If you would like more detail on the initial mass transfer, you can set mt_history=True.

This will write the formation channels to the Population file, which can be accessed by Population.formation_channels.

[ ]:

pop.calculate_formation_channels(mt_history=True)

[ ]:

# the formation channels are loaded in with pop.formation_channels
pop.formation_channels

Next time you open this population file, the formation_channels will be available without having to be recalculated.

[ ]:

pop = Population('1e+00_Zsun_population.h5')
pop.formation_channels

Selecting a sub-population

You might just want a small sub-selection of the full population, especially if you’re working with large population and multi-metallicity runs.

The Population.export_selection() function will export just the indices of the binaries you’re interested in into a new file. The simulated and underlying mass will remain the same, since they are dependent on the population run.

If we select just 2 binaries and export them, we create a new population of just the binaries you’re interested in. In the `BBH analysis <>`__ and `GRB analysis <>`__ tutorials, we show how to perform a selection with multiple criteria and metallicities.

[ ]:

indices = [0,9]
pop.export_selection(indices, 'selected.h5')

[ ]:

selected = Population('selected.h5')
selected.mass_per_metallicity

If you would like to know the simulated mass of just your population, you can calulate this using the online ZAMS values.

[ ]:

import numpy as np

print('selected simulated mass: ', np.sum(selected.oneline[['S1_mass_i', 'S2_mass_i']].to_numpy()))

By default this export-selection will not overwrite nor append if the output file is already present. You have to explicitly state what you would like to append to or overwrite the population file.

With the append=True you are able to combine multiple stellar populations into a single file. This is especially useful when creating multi-metallicity populations.

[ ]:

# this will overwrite the existing file
pop.export_selection(indices, 'selected.h5', overwrite=True)

selected = Population('selected.h5')
selected.mass_per_metallicity

[ ]:

# This will add to the file and add the extra simulated mass
pop.export_selection(indices, 'selected.h5', append=True)
selected = Population('selected.h5')
selected.mass_per_metallicity

Feel free to explore the small binary population you’ve just created!

If you want to learn more about population synthesis and how to build more complex models it is advised to continue with the remaining tutorials and consult the POSYDON documentation.

Local MPI runs

To speed up population synthesis runs, you can run on a computing cluster, as described in HPC Facilities, or you can distribute the population synthesis across multiple cores on your local machine using MPI.

To enable local MPI runs, go into the population_params.ini and change use_MPI to True.

It’s important to note that you cannot run have this option enabled for cluster runs!

We create a binary population simulation script to run the population:

[ ]:

%%writefile script.py
from posydon.popsyn.synthetic_population import PopulationRunner

if __name__ == "__main__":
    synth_pop = PopulationRunner("./population_params.ini")
    synth_pop.evolve()

This script can be initiated using a local where NR_processors is the number of processors you would like to us.

[ ]:

mpiexec -n ${NR_processors} python script.py

This will create a folder for each metallicity in the population and store output of the parallel runs in it.

You will have to concatenate these runs manually into a single population file per metallicity, which can be achieved using the following code:

[ ]:

from posydon.popsyn.synthetic_population import PopulationRunner

synth_pop = PopulationRunner("./population_params.ini")
for pop in synth_pop.binary_populations:
    synth_pop.merge_parallel_runs(pop)