Dynamic time resolution¶
COSMIC has the ability to set time resolution of the bcm array depending on the current state of the evolution. It does this using the timestep_conditions parameter, which we’ll get into below.
How it works¶
Simply put, you tell COSMIC what time resolution to use during different phases of evolution based on a series of conditions. These conditions can cover everything from the stellar type and masses of the stars, to whether or not they are undergoing mass transfer, to the binary state, etc. Below we’ll talk through how to specify those conditions, and then demonstrate their usage in a few worked examples.
Timestep condition syntax¶
The timestep resolution for the bcm array can be set using timestep_conditions parameter. The format of this parameter is a list of up to 15 lists, where each sublist contains zero or more conditions followed by the desired resolution. The conditions are in the form of a string with the format 'column_name operator condition'. The resolution is in the form of a string with the format 'dtp=resolution'. The columns that can be used are listed in cosmic.checkstate.CHECKSTATE_COLUMNS.
The conditions are evaluated at each time step and if the condition is met, the resolution is set to the desired value. Within a sublist, all of the conditions must be met for its resolution to apply. A sublist with no conditions (just a dtp) always applies. If more than one sublist applies at the same time, the last one in the list sets the resolution, so put general conditions first and more specific ones after. When no sublist applies, the resolution goes back to the dtp used for the binary.
The resolution can be any non-negative float value given in units of Myr. If the resolution is set to 0.0, every time step that COSMIC uses internally will be outputted in the bcm array. If the resolution is set to None then time resolution ends and only the final state is outputted.
Conditions can use the operators ==, !=, <, <=, > and >=. A column can be given both a lower and an upper bound, either as a chained comparison (e.g. '1<kstar_1<4') or as two separate conditions (e.g. 'kstar_1>1', 'kstar_1<4'), as well as one !=. Values can be negative (e.g. 'sep==-1') and spaces around the operators are ignored. For example, the following uses a coarse resolution on the main sequence, a fine resolution on the giant branch, every time step during core helium burning and a medium resolution until the star becomes a remnant
timestep_conditions = [['kstar_1 < 2', 'dtp=5.0'], ['1 < kstar_1 < 4', 'dtp=0.01'],
['kstar_1 == 4', 'dtp=0.0'], ['4 < kstar_1 < 10', 'dtp=0.5']]
Setting the same bound twice for a column in one sublist (e.g. 'kstar_1 > 1', 'kstar_1 > 2'), using an unsupported operator or leaving out the dtp raises a ValueError.
Syntax examples¶
Before we dive in to using these with COSMIC, here are some examples of valid syntax for the timestep_conditions parameter. Each example describes the condition and then shows the corresponding syntax. Note that the order of the conditions in a sublist does not matter, but the order of the sublists does.
Basic
# increase the timestep resolution during the Hertzsprung gap (``kstar_1==2``)
timestep_conditions = [['kstar_1 == 2', 'dtp=0.01']]
Ranged conditions
# show all timesteps when a primary star has a specific temperature range
timestep_conditions = [['10000 < teff_1 < 20000', 'dtp=0.0']]
Multiple conditions (AND)
# show all timesteps when a primary star has a specific temperature range
# AND it is on its main sequence (``kstar_1<=1``)
# any of these are equivalent and valid
timestep_conditions = [['10000 < teff_1 < 20000', 'kstar_1 <= 1', 'dtp=0.0']]
timestep_conditions = [['10000 < teff_1 < 20000 and kstar_1 <= 1', 'dtp=0.0']]
timestep_conditions = [['10000 < teff_1 < 20000 & kstar_1 <= 1', 'dtp=0.0']]
Multiple separate conditions (OR)
# use 1 Myr timesteps either when the primary is a WD
# OR when it's a black hole
timestep_conditions = [
['10 <= kstar_1 <= 12', 'dtp=1.0'],
['kstar_2 == 14', 'dtp=1.0']
]
Overlapping conditions
# use 1 Myr timesteps when the primary is a WD,
# but use 0.1 Myr timesteps when the secondary is a WD
timestep_conditions = [
['10 <= kstar_1 <= 12', 'dtp=1.0'],
['10 <= kstar_2 <= 12', 'dtp=0.1']
]
# if *both* stars are WDs, the last condition always takes
# precedence so the resolution is set to 0.1 Myr
Worked examples with COSMIC¶
Below we demonstrate three scenarios, setting dtp only during mass transfer, setting dtp to the same resolution for all of the evolution except for after the system merges or is disrupted, and finally an example of setting dtp only during the HMXB stage of the evolution.
All steps during mass transfer¶
First, print all time steps during mass transfer
In [1]: from cosmic.sample.initialbinarytable import InitialBinaryTable
In [2]: from cosmic.evolve import Evolve
In [3]: single_binary = InitialBinaryTable.InitialBinaries(
...: m1=7, m2=5, porb=3000,
...: ecc=0.42, tphysf=13700.0, kstar1=1, kstar2=1, metallicity=0.02
...: )
...:
In [4]: SSEDict = {'stellar_engine': 'sse'}
In [5]: BSEDict = {
...: "pts1": 0.001, "pts2": 0.01, "pts3": 0.02, "zsun": 0.014, "windflag": 3,
...: "eddlimflag": 0, "neta": 0.5, "bwind": 0.0, "hewind": 0.5, "beta": 0.125,
...: "xi": 0.5, "acc2": 1.5, "LBV_flag": 1, "alpha1": [1.0, 1.0],
...: "lambdaf": 0.0, "ce2stageflag": 0, "ceflag": 0, "cekickflag": 2,
...: "cemergeflag": 1, "cehestarflag": 0, "qcflag": 5,
...: "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],
...: "kickflag": 5, "sigma": 265.0, "bhflag": 1, "bhsigmafrac": 1.0,
...: "sigmadiv": -20.0, "ecsn": 2.25, "ecsn_mlow": 1.6, "aic": 1, "ussn": 1,
...: "polar_kick_angle": 90.0,
...: "natal_kick_array": [[-100.0, -100.0, -100.0, -100.0, 0.0], [-100.0, -100.0, -100.0, -100.0, 0.0]],
...: "mm_mu_ns": 400.0, "mm_mu_bh": 200.0, "remnantflag": 4,
...: "fryer_mass_limit": 0, "mxns": 3.0, "fryer_fmix": 1.0,
...: "fryer_mcrit_nsbh": 5.75, "rembar_massloss": 0.5, "wd_mass_lim": 1,
...: "maltsev_mode": 0, "maltsev_fallback": 0.5, "maltsev_pf_prob": 0.1,
...: "pisn": -2, "ppi_co_shift": 0.0, "ppi_extra_ml": 0.0, "bhspinflag": 0,
...: "bhspinmag": 0.0, "grflag": 1, "eddfac": 10, "gamma": -2, "don_lim": -1,
...: "acc_lim": [-1, -1], "smt_periastron_check": 0, "tflag": 1, "ST_tide": 1,
...: "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],
...: "ifflag": 1, "wdflag": 1, "epsnov": 0.001, "bdecayfac": 1,
...: "bconst": 3000, "ck": 1000, "rejuv_fac": 1.0, "rejuvflag": 0,
...: "bhms_coll_flag": 0, "htpmb": 1, "ST_cr": 1, "rtmsflag": 0
...: }
...:
# Define the condition for the time step
In [6]: timestep_conditions = [['RRLO_1>=1', 'dtp=0.0'], ['RRLO_2>=1', 'dtp=0.0']]
In [7]: bpp, bcm, initC, kick_info = Evolve.evolve(
...: initialbinarytable=single_binary, BSEDict=BSEDict, SSEDict=SSEDict,
...: timestep_conditions=timestep_conditions
...: )
...:
In [8]: print(bcm[['tphys', 'kstar_1', 'kstar_2', 'mass_1', 'mass_2', 'RRLO_1', 'RRLO_2']])
tphys kstar_1 kstar_2 mass_1 mass_2 RRLO_1 RRLO_2
0 0.000000 1 1 7.000000 5.000000 0.007112 0.006828
0 55.150466 6 1 6.656391 5.009692 1.001003 0.006622
0 105.213676 11 2 1.003568 5.003370 0.000552 1.000792
0 105.213676 15 5 0.000000 5.008469 -1.000000 0.000100
0 13700.000000 15 11 0.000000 1.018756 -1.000000 0.000100
Condition on binary state¶
Second, pick a certain resolution for the bcm array until the system merges or is disrupted and then only print the final state
In [9]: timestep_conditions = [['binstate=0', 'dtp=1.0']]
In [10]: bpp, bcm, initC, kick_info = Evolve.evolve(
....: initialbinarytable=single_binary,
....: BSEDict=BSEDict, SSEDict=SSEDict,
....: timestep_conditions=timestep_conditions
....: )
....:
In [11]: print(bcm[['tphys', 'kstar_1', 'kstar_2', 'mass_1', 'mass_2', 'bin_state']])
tphys kstar_1 kstar_2 mass_1 mass_2 bin_state
0 0.0 1 1 7.000000 5.000000 0
0 1.0 1 1 6.999499 4.999955 0
0 2.0 1 1 6.998990 4.999909 0
0 3.0 1 1 6.998471 4.999863 0
0 4.0 1 1 6.997943 4.999817 0
.. ... ... ... ... ... ...
0 103.0 11 1 1.003565 5.003677 0
0 104.0 11 1 1.003565 5.003677 0
0 105.0 11 2 1.003566 5.003527 0
0 106.0 15 11 0.000000 1.018755 1
0 13700.0 15 11 0.000000 1.018755 1
[108 rows x 6 columns]
Condition on evolutionary state¶
Finally, we show how to print a fine resolution only during the HMXB stage of the evolution.
In [12]: single_binary = InitialBinaryTable.InitialBinaries(
....: m1=80, m2=40, porb=350,
....: ecc=0.5, tphysf=13700.0,
....: kstar1=1, kstar2=1, metallicity=0.002
....: )
....:
In [13]: timestep_conditions = [['kstar_1 == 14', 'kstar_2 < 10','dtp=0.1'], ['kstar_2 == 14', 'kstar_1 < 10','dtp=0.1']]
In [14]: bpp, bcm, initC, kick_info = Evolve.evolve(
....: initialbinarytable=single_binary,
....: BSEDict=BSEDict, SSEDict=SSEDict,
....: timestep_conditions=timestep_conditions
....: )
....:
In [15]: print(bcm[['tphys', 'kstar_1', 'kstar_2', 'mass_1', 'mass_2', 'bin_state']])
tphys kstar_1 kstar_2 mass_1 mass_2 bin_state
0 0.00000 1 1 80.000000 40.000000 0
0 4.25508 14 1 23.192518 40.490392 0
0 4.35508 14 1 23.192542 40.410193 0
0 4.45508 14 1 23.192569 40.328570 0
0 4.55508 14 1 23.192598 40.245923 0
0 4.65508 14 1 23.192631 40.162792 0
0 4.75508 14 1 23.192667 40.079893 0
0 4.85508 14 1 23.192707 39.998143 0
0 4.95508 14 1 23.192752 39.918693 0
0 5.05508 14 1 23.192801 39.842924 0
0 5.15508 14 1 23.192855 39.772429 0
0 5.25508 14 1 23.192914 39.707754 0
0 5.35508 14 7 23.501820 14.515722 0
0 5.45508 14 7 23.501827 14.201333 0
0 5.55508 14 7 23.501834 13.870113 0
0 5.65508 14 7 23.501843 13.526101 0
0 5.75508 14 7 23.501851 13.173271 0
0 5.85508 14 14 23.501856 10.537572 2
0 13700.00000 14 14 23.501856 10.537572 2