Generate a binary population by hand

The process to generate a synthetic binary population, is similar to the process to evolve a single/multiple binaries by hand: first generate an initial population, then evolve it with the Evolve class.

An initialized binary population consists of a collection of binary systems with assigned primary and secondary masses, orbital periods, eccentricities, metallicities, and star formation histories. These parameters are randomly sampled from observationally motivated distribution functions.

In COSMIC, the initial sample is done through an initial binary sampler which works with the InitialBinaryTable class. There are two samplers available:

1. independent : initialize binaries with independent parameter distributions for the primary mass, mass ratio, eccentricity, separation, and binary fraction

2. multidim : initialize binaries with multidimensional parameter distributions according to Moe & Di Stefano 2017

We consider both cases below.

independent

First import the InitialBinaryTable class and the independent sampler

In [1]: from cosmic.sample.initialbinarytable import InitialBinaryTable

In [2]: from cosmic.sample.sampler import independent

The independent sampler contains multiple models for each binary parameter. You can access the available models using the independent sampler help call:

In [3]: help(independent.get_independent_sampler)
Help on function get_independent_sampler in module cosmic.sample.sampler.independent:

get_independent_sampler(final_kstar1, final_kstar2, primary_model, ecc_model, porb_model, SF_start, SF_duration, binfrac_model, met, size, **kwargs)
    Generates an initial binary sample according to user specified models
    
    Parameters
    ----------
    final_kstar1 : `int or list`
        Int or list of final kstar1
    
    final_kstar2 : `int or list`
        Int or list of final kstar2
    
    primary_model : `str`
        Model to sample primary mass; choices include: kroupa93, kroupa01, salpeter55, custom
        if 'custom' is selected, must also pass arguemts:
        alphas : `array`
            list of power law indicies
        mcuts : `array`
            breaks in the power laws.
        e.g. alphas=[-1.3,-2.3,-2.3],mcuts=[0.08,0.5,1.0,150.] reproduces standard Kroupa2001 IMF
    
    ecc_model : `str`
        Model to sample eccentricity; choices include: thermal, uniform, sana12
    
    porb_model : `str`
        Model to sample orbital period; choices include: log_uniform, sana12
    
    qmin : `float`
        kwarg which sets the minimum mass ratio for sampling the secondary
        where the mass ratio distribution is flat in q
        if q > 0, qmin sets the minimum mass ratio
        q = -1, this limits the minimum mass ratio to be set such that
        the pre-MS lifetime of the secondary is not longer than the full
        lifetime of the primary if it were to evolve as a single star
    
    m_max : `float`
        kwarg which sets the maximum primary and secondary mass for sampling
        NOTE: this value changes the range of the IMF and should *not* be used
            as a means of selecting certain kstar types!
    
    m1_min : `float`
        kwarg which sets the minimum primary mass for sampling
        NOTE: this value changes the range of the IMF and should *not* be used
            as a means of selecting certain kstar types!
    
    m2_min : `float`
        kwarg which sets the minimum secondary mass for sampling
        the secondary as uniform in mass_2 between m2_min and mass_1
    
    msort : `float`
        Stars with M>msort can have different pairing and sampling of companions
    
    qmin_msort : `float`
        Same as qmin for M>msort
    
    m2_min_msort : `float`
        Same as m2_min for M>msort
    
    SF_start : `float`
        Time in the past when star formation initiates in Myr
    
    SF_duration : `float`
        Duration of constant star formation beginning from SF_Start in Myr
    
    binfrac_model : `str or float`
        Model for binary fraction; choices include: vanHaaften or a fraction where 1.0 is 100% binaries
    
    binfrac_model_msort : `str or float`
        Same as binfrac_model for M>msort
    
    met : `float`
        Sets the metallicity of the binary population where solar metallicity is zsun
    
    size : `int`
        Size of the population to sample
    
    zsun : `float`
        optional kwarg for setting effective radii, default is 0.02
    
    Returns
    -------
    InitialBinaryTable : `pandas.DataFrame`
        DataFrame in the format of the InitialBinaryTable
    
    mass_singles : `float`
        Total mass in single stars needed to generate population
    
    mass_binaries : `float`
        Total mass in binaries needed to generate population
    
    n_singles : `int`
        Number of single stars needed to generate a population
    
    n_binaries : `int`
        Number of binaries needed to generate a population

The final_kstar1 and final_kstar2 parameters are lists that contain the kstar types that you would like the final population to contain.

The final kstar is the final state of the binary system we are interested in and is based on the BSE kstar naming convention, see Evolutionary State of the Star for more information.

Thus, if you want to generate a population containing double white dwarfs with CO and ONe WD primaries and He-WD secondaries, the final kstar inputs would be:

In [4]: final_kstar1 = [11, 12]

In [5]: final_kstar2 = [10]

Similar to the help for the sampler, the different models that can be used for each parameter to be sampled can be accessed by the help function for the argument. The syntax for each parameter sample is always: sample_`parameter`. See the example for the star formation history (SFH) below:

In [6]: help(independent.Sample.sample_SFH)
Help on function sample_SFH in module cosmic.sample.sampler.independent:

sample_SFH(self, SF_start=13700.0, SF_duration=0.0, met=0.02, size=None)
    Sample an evolution time for each binary based on a user-specified
    time at the start of star formation and the duration of star formation.
    The default is a burst of star formation 13,700 Myr in the past.
    
    Parameters
    ----------
    SF_start : float
        Time in the past when star formation initiates in Myr
    SF_duration : float
        Duration of constant star formation beginning from SF_Start in Myr
    met : float
        metallicity of the population [Z_sun = 0.02]
        Default: 0.02
    size : int, optional
        number of evolution times to sample
        NOTE: this is set in cosmic-pop call as Nstep
    
    Returns
    -------
    tphys : array
        array of evolution times of size=size
    metallicity : array
        array of metallicities

Using the final kstar inputs above, the initial binary population is sampled as:

In [7]: InitialBinaries, mass_singles, mass_binaries, n_singles, n_binaries = InitialBinaryTable.sampler('independent', final_kstar1, final_kstar2, binfrac_model=0.5, primary_model='kroupa01', ecc_model='sana12', porb_model='sana12', qmin=-1, m2_min=0.08, SF_start=13700.0, SF_duration=0.0, met=0.02, size=10000)

In [8]: print(InitialBinaries)
       kstar_1  kstar_2     mass_1    mass_2          porb       ecc  metallicity   tphysf    mass0_1   mass0_2  rad_1  rad_2  lum_1  ...  B_2  bacc_1  bacc_2  tacc_1  tacc_2  epoch_1  epoch_2  tms_1  tms_2  bhspin_1  bhspin_2  tphys  binfrac
0          1.0      0.0   1.545164  0.527492    508.888756  0.358912         0.02  13700.0   1.545164  0.527492    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0      0.5
1          1.0      1.0   3.157127  2.300120  17362.885312  0.410910         0.02  13700.0   3.157127  2.300120    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0      0.5
2          1.0      1.0   1.810039  0.828554      2.076867  0.364290         0.02  13700.0   1.810039  0.828554    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0      0.5
3          1.0      1.0   1.407712  1.334419    872.238666  0.364673         0.02  13700.0   1.407712  1.334419    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0      0.5
4          1.0      1.0   1.397477  0.839581     32.245373  0.087933         0.02  13700.0   1.397477  0.839581    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0      0.5
...        ...      ...        ...       ...           ...       ...          ...      ...        ...       ...    ...    ...    ...  ...  ...     ...     ...     ...     ...      ...      ...    ...    ...       ...       ...    ...      ...
10469      1.0      1.0   3.110464  1.084257  72553.613837  0.708145         0.02  13700.0   3.110464  1.084257    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0      0.5
10470      1.0      0.0   0.847130  0.562958      2.968470  0.047779         0.02  13700.0   0.847130  0.562958    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0      0.5
10471      1.0      1.0   1.620213  0.786274      6.293226  0.116765         0.02  13700.0   1.620213  0.786274    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0      0.5
10472      1.0      1.0  11.050187  7.539196   7669.899401  0.872235         0.02  13700.0  11.050187  7.539196    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0      0.5
10473      1.0      0.0   1.046963  0.506145     15.377534  0.388461         0.02  13700.0   1.046963  0.506145    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0      0.5

[10474 rows x 38 columns]

NOTE: the length of the initial binary data set, InitialBinaries, does not always match the size parameter provided to InitialBinaryTable.sampler. This is because the sampler accounts for a binary fraction specified by the user with the binfrac_model parameter, which is either a fraction between 0 and 1 or mass dependent following the prescription in van Haaften+2013.

If you want to have separate binary fractions and mass pairings for low and high mass stars, you can by supplying the msort kwarg to the sampler. This sets the mass above which an alternative mass pairing (specified by kwargs qmin_msort and m2_min_msort) and binary fraction model (specified by kwarg binfrac_model_msort) are used. This is handy if you want, for example, a higher binary fraction and more equal mass pairings for high mass stars.

Since we are interested in binaries, we only retain the binary systems that are likely to produce the user specified final kstar types. However, we also keep track of the total mass of the single and binary stars as well as the number of binary and single stars so that we can scale our results to larger populations. If you don’t want to filter the binaries, you can supply final kstars as

In [9]: final_kstars = np.linspace(0, 14, 15)

In [10]: InitialBinaries, mass_singles, mass_binaries, n_singles, n_binaries = InitialBinaryTable.sampler('independent', final_kstars, final_kstars, binfrac_model=0.5, primary_model='kroupa01', ecc_model='sana12', porb_model='sana12', qmin=-1, m2_min=0.08, SF_start=13700.0, SF_duration=0.0, met=0.02, size=10000)

Below we show the effect of different assumptions for the independent initial sampler. The standard assumptions are shown in orange and the results of Sana et al. 2012 are shown in orange.

(Source code, png, hires.png, pdf)

../_images/index-11.png

multidim

COSMIC implements multidimensionally distributed initial binaries according to Moe & Di Stefano 2017. The python code used in COSMIC to create this sample was written by Mads Sorenson, and is based on the IDL codes written to accompany Moe & Di Stefano 2017.

The multidimensional initial binary data is sampled in COSMIC as follows:

In [11]: from cosmic.sample.initialbinarytable import InitialBinaryTable

In [12]: from cosmic.sample.sampler import multidim

To see the arguments necessary to call the multidimensional sampler use the help function:

In [13]: help(multidim.get_multidim_sampler)
Help on function get_multidim_sampler in module cosmic.sample.sampler.multidim:

get_multidim_sampler(final_kstar1, final_kstar2, rand_seed, nproc, SF_start, SF_duration, met, size, **kwargs)
    adapted version of Maxwell Moe's IDL code that generates a population of single and binary stars
    
    Below is the adapted version of Maxwell Moe's IDL code
    that generates a population of single and binary stars
    based on the paper Mind your P's and Q's
    By Maxwell Moe and Rosanne Di Stefano
    
    The python code has been adopted by Mads Sørensen
    
    Version history:
    V. 0.1; 2017/02/03
    By Mads Sørensen
    - This is a pure adaption from IDL to Python.
    - The function idl_tabulate is similar to
    the IDL function int_tabulated except, this function seems to be slightly
    more exact in its solution.
    Therefore, relative to the IDL code, there are small numerical differences.
    
    Comments below beginning with ; is the original nodes by Maxwell Moe.
    Please read these careful for understanding the script.
    ; NOTE - This version produces only the statistical distributions of
    ;        single stars, binaries, and inner binaries in hierarchical triples.
    ;        Outer tertiaries in hierarchical triples are NOT generated.
    ;        Moreover, given a set of companions, all with period P to
    ;        primary mass M1, this version currently uses an approximation to
    ;        determine the fraction of those companions that are inner binaries
    ;        vs. outer triples. Nevertheless, this approximation reproduces
    ;        the overall multiplicity statistics.
    ; Step 1 - Tabulate probably density functions of periods,
    ;          mass ratios, and eccentricities based on
    ;          analytic fits to corrected binary star populations.
    ; Step 2 - Implement Monte Carlo method to generate stellar
    ;          population from those density functions.
    
    Parameters
    ----------
    final_kstar1 : `list` or `int`
        Int or list of final kstar1
    
    final_kstar2 : `list` or `int`
        Int or list of final kstar2
    
    rand_seed : `int`
        Int to seed random number generator
    
    nproc : `int`
        Number of processors to use to generate population
    
    SF_start : `float`
            Time in the past when star formation initiates in Myr
    
    SF_duration : `float`
            Duration of constant star formation beginning from SF_Start in Myr
    
    met : `float`
        Sets the metallicity of the binary population where solar metallicity is 0.02
    
    size : `int`
        Size of the population to sample
    
    **porb_lo : `float`
        Lower limit in days for the orbital period distribution
    
    **porb_hi: `float`
        Upper limit in days for the orbital period distribution
    
    Returns
    -------
    InitialBinaryTable : `pandas.DataFrame`
        DataFrame in the format of the InitialBinaryTable
    
    mass_singles : `float`
        Total mass in single stars needed to generate population
    
    mass_binaries : `float`
        Total mass in binaries needed to generate population
    
    n_singles : `int`
        Number of single stars needed to generate a population
    
    n_binaries : `int`
        Number of binaries needed to generate a population

The random seed is used to reproduce your initial sample, since there are several stochastic processes involved in the muldimensional sample. As in the independent sampler, the final_kstar1 and final_kstar2 inputs are lists containing the kstar types that the evolved population should contain.

The multidimensional sample is generated as follows:

In [14]: InitialBinaries, mass_singles, mass_binaries, n_singles, n_binaries = InitialBinaryTable.sampler('multidim', final_kstar1=[11], final_kstar2=[11], rand_seed=2, nproc=1, SF_start=13700.0, SF_duration=0.0, met=0.02, size=10)

In [15]: print(InitialBinaries)
   kstar_1  kstar_2     mass_1    mass_2          porb       ecc  metallicity   tphysf    mass0_1   mass0_2  rad_1  rad_2  lum_1  ...  B_2  bacc_1  bacc_2  tacc_1  tacc_2  epoch_1  epoch_2  tms_1  tms_2  bhspin_1  bhspin_2  tphys   binfrac
0      1.0      1.0   2.441146  0.973536  6.063103e+01  0.816199         0.02  13700.0   2.441146  0.973536    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0  0.528636
1      1.0      1.0   2.198825  1.245038  5.712584e+04  0.486927         0.02  13700.0   2.198825  1.245038    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0  0.506297
2      1.0      1.0   1.644136  1.411615  1.777354e+02  0.749940         0.02  13700.0   1.644136  1.411615    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0  0.453639
3      1.0      1.0   4.203177  1.976197  4.249118e+02  0.593713         0.02  13700.0   4.203177  1.976197    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0  0.652902
4      1.0      1.0   2.513813  1.417578  6.293644e+06  0.543736         0.02  13700.0   2.513813  1.417578    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0  0.528636
5      1.0      1.0   0.878763  0.833142  2.497349e+00  0.011429         0.02  13700.0   0.878763  0.833142    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0  0.404277
6      1.0      1.0   3.730858  0.863017  4.735979e+01  0.235551         0.02  13700.0   3.730858  0.863017    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0  0.614692
7      1.0      1.0  15.081145  9.660231  7.947508e+02  0.748707         0.02  13700.0  15.081145  9.660231    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0  0.891344
8      1.0      1.0   0.942717  0.863542  3.450768e+02  0.716432         0.02  13700.0   0.942717  0.863542    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0  0.405617
9      1.0      1.0   0.975778  0.962636  1.539399e+03  0.405061         0.02  13700.0   0.975778  0.962636    0.0    0.0    0.0  ...  0.0     0.0     0.0     0.0     0.0      0.0      0.0    0.0    0.0       0.0       0.0    0.0  0.406506

[10 rows x 38 columns]

Note

NOTE that in the multidimensional case, the binary fraction is a parameter in the sample. This results in the size of the initial binary data matching the size provided to the sampler. As in the independent sampling case, we keep track of the total sampled mass of singles and binaries as well as the total number of single and binary stars to scale thesimulated population to astrophysical populations.

(Source code)

Evolving an initial binary population

As in Using COSMIC to evolve binaries, once an initial binary population is sampled, it is evolved using the Evolve class. Note that the same process used in Using COSMIC to evolve binaries applies here as well: the BSEDict must be supplied, but only need be resupplied if the flags in the dictionary change.

The syntax for the Evolve class is as follows:

In [16]: from cosmic.evolve import Evolve

In [17]: BSEDict = {'xi': 1.0, 'bhflag': 1, 'neta': 0.5, 'windflag': 3, 'wdflag': 1, 'alpha1': 1.0, 'pts1': 0.001, 'pts3': 0.02, 'pts2': 0.01, 'epsnov': 0.001, 'hewind': 0.5, 'ck': 1000, 'bwind': 0.0, 'lambdaf': 0.0, 'mxns': 3.0, 'beta': -1.0, 'tflag': 1, 'acc2': 1.5, 'grflag' : 1, 'remnantflag': 4, 'ceflag': 0, 'eddfac': 1.0, 'ifflag': 0, 'bconst': 3000, 'sigma': 265.0, 'gamma': -2.0, 'pisn': 45.0, 'natal_kick_array' : [[-100.0,-100.0,-100.0,-100.0,0.0], [-100.0,-100.0,-100.0,-100.0,0.0]], 'bhsigmafrac' : 1.0, 'polar_kick_angle' : 90, 'qcrit_array' : [0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0], 'cekickflag' : 2, 'cehestarflag' : 0, 'cemergeflag' : 0, 'ecsn' : 2.25, 'ecsn_mlow' : 1.6, 'aic' : 1, 'ussn' : 0, 'sigmadiv' :-20.0, 'qcflag' : 1, 'eddlimflag' : 0, 'fprimc_array' : [2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0,2.0/21.0], 'bhspinflag' : 0, 'bhspinmag' : 0.0, 'rejuv_fac' : 1.0, 'rejuvflag' : 0, 'htpmb' : 1, 'ST_cr' : 1, 'ST_tide' : 1, 'bdecayfac' : 1, 'rembar_massloss' : 0.5, 'kickflag' : 0, 'zsun' : 0.019, 'bhms_coll_flag' : 0, 'don_lim' : -1, 'acc_lim' : -1}

In [18]: bpp, bcm, initC, kick_info  = Evolve.evolve(initialbinarytable=InitialBinaries, BSEDict=BSEDict)

In [19]: print(bcm.iloc[:10])
     tphys  kstar_1   mass0_1    mass_1       lum_1     rad_1        teff_1   massc_1    radc_1        menv_1        renv_1  ...        RRLO_2          porb            sep       ecc  B_1  B_2  SN_1  SN_2  bin_state  merger_type  bin_num
0      0.0      1.0  2.441146  2.441146   35.232889  1.784981  10572.795481  0.000000  0.000000  1.000000e-10  1.000000e-10  ...  1.583882e-01  6.063103e+01      97.768428  0.816199  0.0  0.0   0.0   0.0          0         -001        0
0  13700.0     11.0  0.709369  0.709369    0.000004  0.011381   2375.330102  0.709369  0.011381  1.000000e-10  1.000000e-10  ... -1.000000e+00  0.000000e+00       0.000000 -1.000000  0.0  0.0   0.0   0.0          1         0301        0
1      0.0      1.0  2.198825  2.198825   23.108512  1.693523   9768.255374  0.000000  0.000000  1.000000e-10  1.000000e-10  ...  7.513324e-04  5.712584e+04    9422.994116  0.486927  0.0  0.0   0.0   0.0          0         -001        1
1  13700.0     11.0  0.635782  0.636342    0.000004  0.012297   2308.729961  0.636342  0.012297  1.000000e-10  1.000000e-10  ...  2.669142e-06  4.688219e+05   26820.943009  0.480578  0.0  0.0   0.0   0.0          0         -001        1
2      0.0      1.0  1.644136  1.644136    6.800611  1.501786   7640.162387  0.000000  0.000000  1.000000e-10  1.000000e-10  ...  7.842567e-02  1.777354e+02     192.977710  0.749940  0.0  0.0   0.0   0.0          0         -001        2
2  13700.0     11.0  0.640720  0.640720    0.000008  0.012241   2774.448733  0.640720  0.012241  1.000000e-10  1.000000e-10  ... -1.000000e+00  0.000000e+00       0.000000 -1.000000  0.0  0.0   0.0   0.0          1         0110        2
3      0.0      1.0  4.203177  4.203177  280.998633  2.406241  15302.981667  0.000000  0.000000  1.000000e-10  1.000000e-10  ...  2.880054e-02  4.249118e+02     436.314107  0.593713  0.0  0.0   0.0   0.0          0         -001        3
3  13700.0     15.0  0.836652  0.000000    4.544572  0.008210  93428.548967  0.986530  0.008210  1.000000e-10  1.000000e-10  ...  1.000000e-04  1.821547e-02       0.129137 -1.000000  0.0  0.0   0.0   0.0          0         -001        3
4      0.0      1.0  2.513813  2.513813   39.600085  1.812518  10803.200211  0.000000  0.000000  1.000000e-10  1.000000e-10  ...  4.069085e-05  6.293644e+06  226326.490102  0.543736  0.0  0.0   0.0   0.0          0         -001        4
4  13700.0     11.0  0.686027  0.686028    0.000003  0.011668   2306.395953  0.686028  0.011668  1.000000e-10  1.000000e-10  ...  1.120193e-07  6.321749e+07  718868.035718  0.543713  0.0  0.0   0.0   0.0          0         -001        4

[10 rows x 39 columns]

In [20]: print(bpp)
           tphys    mass_1    mass_2  kstar_1  kstar_2         sep         porb       ecc    RRLO_1    RRLO_2  evol_type          aj_1  ...  omega_spin_2  B_1  B_2  bacc_1  bacc_2  tacc_1  tacc_2    epoch_1    epoch_2  bhspin_1  bhspin_2  bin_num
0       0.000000  2.441146  0.973536      1.0      1.0   97.768428    60.631028  0.816199  0.215872  0.158388        1.0      0.000000  ...    554.100227  0.0  0.0     0.0     0.0     0.0     0.0   0.000000   0.000000       0.0       0.0        0
0     666.708170  2.441146  0.973536      2.0      1.0   97.768028    60.630656  0.816198  0.497026  0.160097        2.0    666.708170  ...    542.236341  0.0  0.0     0.0     0.0     0.0     0.0   0.000000   0.000000       0.0       0.0        0
0     670.696861  2.441021  0.973537      2.0      1.0   17.973508     4.779189  0.000000  1.000511  0.160073        3.0    670.792519  ...    480.179965  0.0  0.0     0.0     0.0     0.0     0.0  -0.095658   0.002708       0.0       0.0        0
0     671.325357  2.437581  0.976959      3.0      1.0   17.887815     4.745063  0.000000  1.121076  0.161241        2.0    673.226666  ...    936.783251  0.0  0.0     0.0     0.0     0.0     0.0  -1.901310   9.088475       0.0       0.0        0
0     671.325357  2.437581  0.976959      3.0      1.0   17.887815     4.745063  0.000000  1.121076  0.161241        7.0    673.226666  ...    936.783251  0.0  0.0     0.0     0.0     0.0     0.0  -1.901310   9.088475       0.0       0.0        0
..           ...       ...       ...      ...      ...         ...          ...       ...       ...       ...        ...           ...  ...           ...  ...  ...     ...     ...     ...     ...        ...        ...       ...       ...      ...
9   12274.884934  0.975778  0.962636      2.0      1.0  699.307674  1539.399391  0.405061  0.009971  0.008788        2.0  12274.884934  ...    195.795345  0.0  0.0     0.0     0.0     0.0     0.0   0.000000   0.000000       0.0       0.0        9
9   12921.550046  0.975252  0.962636      3.0      1.0  699.497439  1540.235065  0.405061  0.013824  0.009830        2.0  12948.031011  ...    156.362675  0.0  0.0     0.0     0.0     0.0     0.0 -26.480965   0.002465       0.0       0.0        9
9   12940.183607  0.975232  0.962636      3.0      2.0  699.504606  1540.266629  0.405061  0.014034  0.009871        2.0  12966.664572  ...    174.363440  0.0  0.0     0.0     0.0     0.0     0.0 -26.480965   0.002592       0.0       0.0        9
9   13621.915879  0.969860  0.962107      3.0      3.0  701.633743  1549.665737  0.405060  0.060484  0.013700        2.0  13648.396845  ...    192.407385  0.0  0.0     0.0     0.0     0.0     0.0 -26.480965 -28.604601       0.0       0.0        9
9   13700.000000  0.928905  0.962729      3.0      3.0  707.627279  1586.209059  0.396045  0.312854  0.014122       10.0  13726.480965  ...    173.807788  0.0  0.0     0.0     0.0     0.0     0.0 -26.480965 -28.604601       0.0       0.0        9

[110 rows x 44 columns]

ClusterMonteCarlo

New in COSMIC 3.4, you can now use COSMIC to sample initial conditions that can be used in the simulation of a Globular Cluster (GC), using the ClusterMonteCarlo (CMC) software package. To create these initial conditions, and save them in a format readable by CMC, you can do the following.

In [21]: from cosmic.sample.initialcmctable import InitialCMCTable

In [22]: from cosmic.sample.sampler import cmc

To see the arguments necessary to call the CMC sampler use the help function:

In [23]: help(cmc.get_cmc_sampler)
Help on function get_cmc_sampler in module cosmic.sample.sampler.cmc:

get_cmc_sampler(cluster_profile, primary_model, ecc_model, porb_model, binfrac_model, met, size, **kwargs)
    Generates an initial binary sample according to user specified models
    
    Parameters
    ----------
    cluster_profile : `str`
        Model to use for the cluster profile (i.e. sampling of the placement of objects in the cluster and their velocity within the cluster)
        Options include:
        'plummer' : Standard Plummer sphere.
            Additional parameters:
            'r_max' : `float`
                the maximum radius (in virial radii) to sample the clsuter
        'elson' : EFF 1987 profile.  Generalization of Plummer that better fits young massive clusters
            Additional parameters:
            'gamma' : `float`
                steepness paramter for Elson profile; note that gamma=4 is same is Plummer
            'r_max' : `float`
                the maximum radius (in virial radii) to sample the clsuter
        'king' : King profile 
            'w_0' : `float`
                King concentration parameter
            'r_max' : `float`
                the maximum radius (in virial radii) to sample the clsuter
    
    primary_model : `str`
        Model to sample primary mass; choices include: kroupa93, kroupa01, salpeter55, custom
        if 'custom' is selected, must also pass arguemts:
        alphas : `array`
            list of power law indicies
        mcuts : `array`
            breaks in the power laws.
        e.g. alphas=[-1.3,-2.3,-2.3],mcuts=[0.08,0.5,1.0,150.] reproduces standard Kroupa2001 IMF
    
    ecc_model : `str`
        Model to sample eccentricity; choices include: thermal, uniform, sana12
    
    porb_model : `str`
        Model to sample orbital period; choices include: log_uniform, sana12
    
    msort : `float`
        Stars with M>msort can have different pairing and sampling of companions
    
    pair : `float`
        Sets the pairing of stars M>msort only with stars with M>msort
    
    binfrac_model : `str or float`
        Model for binary fraction; choices include: vanHaaften or a fraction where 1.0 is 100% binaries
    
    binfrac_model_msort : `str or float`
        Same as binfrac_model for M>msort
    
    qmin : `float`
        kwarg which sets the minimum mass ratio for sampling the secondary
        where the mass ratio distribution is flat in q
        if q > 0, qmin sets the minimum mass ratio
        q = -1, this limits the minimum mass ratio to be set such that
        the pre-MS lifetime of the secondary is not longer than the full
        lifetime of the primary if it were to evolve as a single star
    
    m2_min : `float`
        kwarg which sets the minimum secondary mass for sampling
        the secondary as uniform in mass_2 between m2_min and mass_1
    
    qmin_msort : `float`
        Same as qmin for M>msort; only applies if qmin is supplied
    
    met : `float`
        Sets the metallicity of the binary population where solar metallicity is zsun 
    
    size : `int`
        Size of the population to sample
    
    zsun : `float`
        optional kwarg for setting effective radii, default is 0.02
    
    Optional Parameters
    -------------------
    virial_radius : `float`
        the initial virial radius of the cluster, in parsecs
        Default -- 1 pc
    
    tidal_radius : `float`
        the initial tidal radius of the cluster, in units of the virial_radius
        Default -- 1e6 rvir
    
    seed : `float`
        seed to the random number generator, for reproducability
    
    Returns
    -------
    Singles: `pandas.DataFrame`
        DataFrame of Single objects in the format of the InitialCMCTable
    
    Binaries: `pandas.DataFrame`
        DataFrame of Single objects in the format of the InitialCMCTable
In [24]: from cosmic.sample import InitialCMCTable

In [25]: Singles, Binaries = InitialCMCTable.sampler('cmc', binfrac_model=0.2, primary_model='kroupa01', ecc_model='sana12', porb_model='sana12', qmin=-1.0, cluster_profile='plummer', met=0.014, size=40000, params='../examples/Params.ini', gamma=4, r_max=100)

In [26]: InitialCMCTable.write(Singles, Binaries, filename="input.hdf5")

In [27]: InitialCMCTable.write(Singles, Binaries, filename="input.fits")