Dynamic time resolution¶
COSMIC has the ability to set time resolution of the bcm array depending on the current state of the evolution.
In [1]: from cosmic.sample.initialbinarytable import InitialBinaryTable
In [2]: from cosmic.evolve import Evolve
General syntax¶
The timestep resolution for the bcm array can be set using timestep_conditions
parameter. The format of
this parameter is a list of 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=condition'
. The resolution
is in the form of a string with the format 'dtp=resolution'
.
The conditions are evaluated at each time step and if the condition is met, the resolution is set to the desired value. If multiple conditions are specified, a timestep is outputted if any of the conditions are met.
The resolution can be any positive float value given in units of Myr. If the resolution is set to 0.0, every time step will be outputted in the bcm array.
Examples¶
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 [3]: single_binary = InitialBinaryTable.InitialBinaries(m1=7.806106, m2=5.381412, porb=2858.942021,
...: ecc=0.601408, tphysf=13700.0, kstar1=1, kstar2=1, metallicity=0.02)
...:
In [4]: 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' : 1, 'zsun' : 0.014, 'bhms_coll_flag' : 0, 'don_lim' : -1, 'acc_lim' : -1, 'rtmsflag' : 0, 'wd_mass_lim': 1}
In [5]: # Define the condition for the time step
...: timestep_conditions = [['RRLO_1>=1', 'dtp=0.0'], ['RRLO_2>=1', 'dtp=0.0']]
...:
In [6]: bpp, bcm, initC, kick_info = Evolve.evolve(initialbinarytable=single_binary, BSEDict=BSEDict,
...: timestep_conditions=timestep_conditions)
...:
In [7]: 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.806106 5.381412 0.010953 0.010459
0 43.565604 5 1 7.346069 5.396914 1.000222 0.008730
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 [8]: timestep_conditions = [['binstate=0', 'dtp=1.0']]
In [9]: bpp, bcm, initC, kick_info = Evolve.evolve(initialbinarytable=single_binary, BSEDict=BSEDict,
...: timestep_conditions=timestep_conditions)
...:
In [10]: 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.806106 5.381412 0
0 1.0 1 1 7.805049 5.381334 0
0 2.0 1 1 7.803970 5.381256 0
0 3.0 1 1 7.802864 5.381177 0
0 4.0 1 1 7.801731 5.381097 0
0 5.0 1 1 7.800567 5.381017 0
0 6.0 1 1 7.799370 5.380936 0
0 7.0 1 1 7.798137 5.380854 0
0 8.0 1 1 7.796865 5.380771 0
0 9.0 1 1 7.795550 5.380687 0
0 10.0 1 1 7.794190 5.380602 0
0 11.0 1 1 7.792780 5.380516 0
0 12.0 1 1 7.791316 5.380430 0
0 13.0 1 1 7.789796 5.380342 0
0 14.0 1 1 7.788213 5.380253 0
0 15.0 1 1 7.786563 5.380163 0
0 16.0 1 1 7.784841 5.380073 0
0 17.0 1 1 7.783043 5.379980 0
0 18.0 1 1 7.781161 5.379887 0
0 19.0 1 1 7.779190 5.379793 0
0 20.0 1 1 7.777123 5.379697 0
0 21.0 1 1 7.774952 5.379600 0
0 22.0 1 1 7.772671 5.379501 0
0 23.0 1 1 7.770271 5.379402 0
0 24.0 1 1 7.767744 5.379301 0
0 25.0 1 1 7.765080 5.379198 0
0 26.0 1 1 7.762270 5.379094 0
0 27.0 1 1 7.759303 5.378988 0
0 28.0 1 1 7.756170 5.378881 0
0 29.0 1 1 7.752860 5.378772 0
0 30.0 1 1 7.749363 5.378661 0
0 31.0 1 1 7.745668 5.378549 0
0 32.0 1 1 7.741767 5.378434 0
0 33.0 1 1 7.737651 5.378318 0
0 34.0 1 1 7.733315 5.378200 0
0 35.0 1 1 7.728755 5.378080 0
0 36.0 1 1 7.723973 5.377958 0
0 37.0 1 1 7.718977 5.377834 0
0 38.0 1 1 7.713779 5.377708 0
0 39.0 4 1 7.651379 5.381606 0
0 40.0 4 1 7.568413 5.386139 0
0 41.0 4 1 7.515623 5.388395 0
0 42.0 4 1 7.468825 5.389992 0
0 43.0 4 1 7.423606 5.390657 0
Condition on evolutionary state¶
Finally, we show how to print a fine resolution only during the HMXB stage of the evolution.
In [11]: single_binary = InitialBinaryTable.InitialBinaries(m1=85.543645, m2=84.99784, porb=446.795757,
....: ecc=0.448872, tphysf=13700.0,
....: kstar1=1, kstar2=1, metallicity=0.002)
....:
In [12]: 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' : 1, 'zsun' : 0.014, 'bhms_coll_flag' : 0, 'don_lim' : -1, 'acc_lim' : -1, 'rtmsflag' : 0, 'wd_mass_lim': 1}
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,
....: 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.000000 1 1 85.543645 84.997840 0
0 4.071374 14 1 24.988585 106.891756 0
0 4.171374 14 1 24.988586 105.634889 0
0 4.271374 14 1 24.988586 104.335873 0
0 4.371374 14 1 24.988586 103.016130 0
0 4.471374 14 1 24.988586 101.707633 0
0 4.571374 14 1 24.988586 100.454871 0
0 4.671374 14 1 24.988586 99.314847 0
0 4.771374 14 1 24.988586 98.352198 0
0 4.871374 14 1 24.988590 89.040292 0