From c5e10079c66974932d4f5d989cbbc53ead2f920d Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Mon, 19 Jul 2021 13:25:19 +0100 Subject: [PATCH 01/22] initial commit of artale models --- skypy/gravitational_waves/__init__.py | 1 + skypy/gravitational_waves/merger_rate.py | 212 +++++++++++++++++- .../tests/test_merger_rate.py | 21 +- 3 files changed, 231 insertions(+), 3 deletions(-) diff --git a/skypy/gravitational_waves/__init__.py b/skypy/gravitational_waves/__init__.py index d0b89fb6b..509f18bd9 100644 --- a/skypy/gravitational_waves/__init__.py +++ b/skypy/gravitational_waves/__init__.py @@ -10,6 +10,7 @@ :toctree: ../api/ b_band_merger_rate + m_star_merger_rate """ diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 8e307d0a3..cb49ec47c 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -5,11 +5,13 @@ """ +import numpy as np from astropy import constants, units __all__ = [ 'b_band_merger_rate', + 'm_star_merger_rate' ] @@ -29,6 +31,214 @@ 'high': 20} } +artale_tables = { + 'NS-NS' : { + 'redshift': [0.1, 1.0, 2.0, 6.0], + 'alpha1': [1.038, 1.109, 1.050, 1.027], + 'alpha1_err': [0.001, 0.001, 0.001, 0.003], + 'alpha2': [-6.09, -6.214, -5.533, -5.029], + 'alpha2_err': [0.010, 0.006, 0.006, 0.021], + 'beta1': [0.800, 0.964, 0.977, 1.113], + 'beta1_err': [0.002, 0.001, 0.001, 0.003], + 'beta2': [0.323, 0.155, 0.068, -0.070], + 'beta2_err': [0.002, 0.001, 0.001, 0.002], + 'beta3': [-3.555, -4.819, -4.874, -5.764], + 'beta3_err': [0.018, 0.013, 0.011, 0.026], + 'gamma1': [0.701, 0.896, 1.018, 1.137], + 'gamma1_err': [0.002, 0.002, 0.002, 0.004], + 'gamma2': [0.356, 0.184, 0.048, -0.082], + 'gamma2_err': [0.002, 0.001, 0.001, 0.002], + 'gamma3': [0.411, 0.222, -0.103, -0.053], + 'gamma3_err': [0.005, 0.003, 0.002, 0.004], + 'gamma4': [-1.968, -3.795, -5.451, -6.104], + 'gamma4_err': [0.026, 0.019, 0.017, 0.037], + } + 'NS-BH' : { + 'redshift': [0.1, 1.0, 2.0, 6.0], + 'alpha1': [0.824, 0.873, 0.913, 0.965], + 'alpha1_err': [0.001, 0.001, 0.001, 0.002], + 'alpha2': [-4.731, -4.478, -4.401, -4.315], + 'alpha2_err': [0.008, 0.008, 0.007, 0.018], + 'beta1': [0.711, 0.813, 0.871, 0.985], + 'beta1_err': [0.002, 0.002, 0.002, 0.003], + 'beta2': [0.150, 0.064, 0.039, -0.017], + 'beta2_err': [0.002, 0.002, 0.001, 0.001], + 'beta3': [-3.536, -3.900, -4.019, -4.490], + 'beta3_err': [0.016, 0.018, 0.014, 0.024], + 'gamma1': [0.833, 1.074, 1.084, 0.978], + 'gamma1_err': [0.002, 0.002, 0.002, 0.003], + 'gamma2': [0.101, -0.058, -0.068, -0.013], + 'gamma2_err': [0.002, 0.002, 0.001, 0.002], + 'gamma3': [-0.461, -0.788, -0.535, 0.016], + 'gamma3_err': [0.004, 0.004, 0.003, 0.004], + 'gamma4': [-5.434, -7.733, -7.055, -4.386], + 'gamma4_err': [0.022, 0.023, 0.018, 0.034], + } + 'BH-BH' : { + 'redshift': [0.1, 1.0, 2.0, 6.0], + 'alpha1': [0.807, 0.813, 0.831, 0.933], + 'alpha1_err': [0.001, 0.001, 0.001, 0.004], + 'alpha2': [-4.310, -3.845, -3.600, -4.190], + 'alpha2_err': [0.006, 0.008, 0.008, 0.026], + 'beta1': [0.812, 0.840, 0.858, 1.053], + 'beta1_err': [0.001, 0.002, 0.002, 0.004], + 'beta2': [-0.006, -0.029, -0.026, -0.098], + 'beta2_err': [0.001, 0.002, 0.001, 0.002], + 'beta3': [-4.358, -4.109, -3.850, -5.213], + 'beta3_err': [0.013, 0.018, 0.015, 0.034], + 'gamma1': [0.921, 1.134, 1.135, 1.131], + 'gamma1_err': [0.001, 0.002, 0.002, 0.005], + 'gamma2': [-0.051, -0.172, -0.167, -0.137], + 'gamma2_err': [0.001, 0.002, 0.001, 0.002], + 'gamma3': [-0.404, -0.839, -0.681, -0.171], + 'gamma3_err': [0.003, 0.004, 0.003, 0.005], + 'gamma4': [-6.049, -8.338, -7.758, -6.321], + 'gamma4_err': [0.018, 0.024, 0.020, 0.047], + } +} + + +def m_star_merger_rate(redshift, + m_star, + population): + r"""Model of Artale et al (2020), equation (1) with parameters + from Tables I, II and III. + + Compact binary merger rates as a log-log function of a galaxy's + stellar mass. + + Parameters are redshift dependent, with linear interpolation + between the simulated points of z={0.1, 1, 2, 6}. + + Parameters + ---------- + redshift : (ngal,) array-like + The redshifts of the galaxies to generate merger + rates for. + m_star : (ngal,) array-like + The stellar mass of the galaxies to generate merger + rates for, in units of stellar mass. + population : {'NS-NS', 'NS-BH', 'BH-BH'} + Compact binary population to get rate for. + 'NS-NS' is neutron star - neutron star + 'NS-BH' is neutron star - black hole + 'BH-BH' is black hole - black hole + + Returns + ------- + merger_rate : array_like + Merger rates for the galaxies in units of Gigayear^-1 + + Notes + ----- + + References + ---------- + .. Artale et al. 2020, MNRAS, + Volume 491, Issue 3, p.3419-3434 (2020) + https://arxiv.org/abs/1910.04890 + + Examples + -------- + >>> import numpy as np + >>> from skypy.gravitational_waves import m_star_merger_rate + + Sample 100 stellar masses values near 10^9 solar masses. + + >>> stellar_masses = 10.**(9.0 + np.random.randn(100)) + + Generate merger rates for these luminosities. + + >>> rates = m_star_merger_rate(redshifts, + ... stellar_masses, + ... population='NS-NS') + + """ + + alpha1 = np.interp(redshift, + artale_tables[population]['redshift'], + artale_tables[population]['alpha1']) + alpha2 = np.interp(redshift, + artale_tables[population]['redshift'], + artale_tables[population]['alpha2']) + + return _m_star_merger_rate(m_star, alpha1, alpha2) + +def m_star_sfr_merger_rate(redshift, + m_star, + sfr, + population): + + beta1 = np.interp(redshift, + artale_tables[population]['redshift'], + artale_tables[population]['beta1']) + beta2 = np.interp(redshift, + artale_tables[population]['redshift'], + artale_tables[population]['beta2']) + beta3 = np.interp(redshift, + artale_tables[population]['redshift'], + artale_tables[population]['beta3']) + + return _m_star_sfr_merger_rate(m_star, sfr, beta1, beta2, beta3) + +def m_star_sfr_metallicity_merger_rate(redshift, + m_star, + sfr, + Z, + population): + + gamma1 = np.interp(redshift, + artale_table_I[population]['redshift'], + artale_table_I[population]['gamma1']) + gamma2 = np.interp(redshift, + artale_table_I[population]['redshift'], + artale_table_I[population]['gamma2']) + gamma3 = np.interp(redshift, + artale_table_I[population]['redshift'], + artale_table_I[population]['gamma3']) + gamma4 = np.interp(redshift, + artale_table_I[population]['redshift'], + artale_table_I[population]['gamma4']) + + return _m_star_sfr_metallicity_merger_rate(m_star, sfr, Z, gamma1, gamma2, gamma3, gamma4) + +def _m_star_merger_rate(m_star, + alpha1, + alpha2): + + m_star = m_star.to(units.Msun).value + + n_gw = 10.**(alpha1 * np.log10(m_star) + alpha2) + + return n_gw / units.Gyr + +def _m_star_sfr_merger_rate(m_star, + sfr, + beta1, + beta2, + beta3): + + m_star = m_star.to(units.Msun).value + sfr = sfr.to(units.Msun / units.year) + + n_gw = 10.**(beta1 * np.log10(m_star) + beta2 * np.log10(sfr) + beta3) + + return n_gw / units.Gyr + +def _m_star_sfr_metallicity_merger_rate(m_star, + sfr, + Z, + gamma1, + gamma2, + gamma3, + gamma4): + + m_star = m_star.to(units.Msun).value + sfr = sfr.to(units.Msun / units.year) + + n_gw = 10.**(gamma1 * np.log10(m_star) + gamma2 * np.log10(sfr) + gamma3 * np.log10(Z) + gamma4) + + return n_gw / units.Gyr def b_band_merger_rate(luminosity, population='NS-NS', @@ -36,7 +246,7 @@ def b_band_merger_rate(luminosity, r"""Model of Abadie et al (2010), Table III - Compact binary merger rates as a linear function of a galaxies + Compact binary merger rates as a linear function of a galaxy's B-band luminosity. Parameters diff --git a/skypy/gravitational_waves/tests/test_merger_rate.py b/skypy/gravitational_waves/tests/test_merger_rate.py index 24152098a..42e298d49 100644 --- a/skypy/gravitational_waves/tests/test_merger_rate.py +++ b/skypy/gravitational_waves/tests/test_merger_rate.py @@ -1,7 +1,7 @@ import numpy as np from astropy import constants, units -from skypy.gravitational_waves import b_band_merger_rate -from skypy.gravitational_waves.merger_rate import abadie_table_III +from skypy.gravitational_waves import b_band_merger_rate, m_star_merger_rate +from skypy.gravitational_waves.merger_rate import abadie_table_III, artale_tables def test_abadie_rates(): @@ -17,3 +17,20 @@ def test_abadie_rates(): L_B_rate = b_band_merger_rate(L_10, population='NS-NS', optimism='low') table_value = abadie_table_III['NS-NS']['low'] assert np.isclose(L_B_rate.to_value(1/units.year), table_value, rtol=1e-5) + +def test_artale_rates(): + + # Check the number of merger rates returned + redshifts = np.random.uniform(0., 3., 100) + stellar_masses = 10.**(9.0 + np.random.randn(100)) + rates = m_star_merger_rate(redshifts, stellar_masses, population='NS-NS') + assert len(rates) == len(stellar_masses) + + # Test that a luminosity of M_sol returns a + # rate that matches the value in Artale Table I + m_sol = constants.Msun.to_value('Msun') + m_sol_rate = m_star_merger_rate(0.0, m_sol, population='NS-NS') + alpha1 = artale_tables['NS-NS']['alpha1'][0] + alpha2 = artale_tables['NS-NS']['alpha2'][0] + table_value = 10.**(alpha1 * np.log10(m_sol) + alpha2) + assert np.isclose(m_sol_rate.to_value(1/units.Gyr), table_value, rtol=1e-5) \ No newline at end of file From 075265e51e59f81b06b89f316e6192dcdd71c2b8 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Mon, 19 Jul 2021 14:23:12 +0100 Subject: [PATCH 02/22] tests pass --- skypy/gravitational_waves/merger_rate.py | 90 ++++++++++--------- .../tests/test_merger_rate.py | 9 +- 2 files changed, 53 insertions(+), 46 deletions(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index cb49ec47c..d7fd30e8e 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -32,7 +32,7 @@ } artale_tables = { - 'NS-NS' : { + 'NS-NS': { 'redshift': [0.1, 1.0, 2.0, 6.0], 'alpha1': [1.038, 1.109, 1.050, 1.027], 'alpha1_err': [0.001, 0.001, 0.001, 0.003], @@ -52,8 +52,8 @@ 'gamma3_err': [0.005, 0.003, 0.002, 0.004], 'gamma4': [-1.968, -3.795, -5.451, -6.104], 'gamma4_err': [0.026, 0.019, 0.017, 0.037], - } - 'NS-BH' : { + }, + 'NS-BH': { 'redshift': [0.1, 1.0, 2.0, 6.0], 'alpha1': [0.824, 0.873, 0.913, 0.965], 'alpha1_err': [0.001, 0.001, 0.001, 0.002], @@ -73,8 +73,8 @@ 'gamma3_err': [0.004, 0.004, 0.003, 0.004], 'gamma4': [-5.434, -7.733, -7.055, -4.386], 'gamma4_err': [0.022, 0.023, 0.018, 0.034], - } - 'BH-BH' : { + }, + 'BH-BH': { 'redshift': [0.1, 1.0, 2.0, 6.0], 'alpha1': [0.807, 0.813, 0.831, 0.933], 'alpha1_err': [0.001, 0.001, 0.001, 0.004], @@ -154,57 +154,60 @@ def m_star_merger_rate(redshift, ... population='NS-NS') """ - + alpha1 = np.interp(redshift, - artale_tables[population]['redshift'], - artale_tables[population]['alpha1']) + artale_tables[population]['redshift'], + artale_tables[population]['alpha1']) alpha2 = np.interp(redshift, - artale_tables[population]['redshift'], - artale_tables[population]['alpha2']) + artale_tables[population]['redshift'], + artale_tables[population]['alpha2']) return _m_star_merger_rate(m_star, alpha1, alpha2) + def m_star_sfr_merger_rate(redshift, m_star, sfr, population): - + beta1 = np.interp(redshift, - artale_tables[population]['redshift'], - artale_tables[population]['beta1']) + artale_tables[population]['redshift'], + artale_tables[population]['beta1']) beta2 = np.interp(redshift, - artale_tables[population]['redshift'], - artale_tables[population]['beta2']) + artale_tables[population]['redshift'], + artale_tables[population]['beta2']) beta3 = np.interp(redshift, - artale_tables[population]['redshift'], - artale_tables[population]['beta3']) + artale_tables[population]['redshift'], + artale_tables[population]['beta3']) return _m_star_sfr_merger_rate(m_star, sfr, beta1, beta2, beta3) + def m_star_sfr_metallicity_merger_rate(redshift, - m_star, - sfr, - Z, - population): - + m_star, + sfr, + Z, + population): + gamma1 = np.interp(redshift, - artale_table_I[population]['redshift'], - artale_table_I[population]['gamma1']) + artale_tables[population]['redshift'], + artale_tables[population]['gamma1']) gamma2 = np.interp(redshift, - artale_table_I[population]['redshift'], - artale_table_I[population]['gamma2']) + artale_tables[population]['redshift'], + artale_tables[population]['gamma2']) gamma3 = np.interp(redshift, - artale_table_I[population]['redshift'], - artale_table_I[population]['gamma3']) + artale_tables[population]['redshift'], + artale_tables[population]['gamma3']) gamma4 = np.interp(redshift, - artale_table_I[population]['redshift'], - artale_table_I[population]['gamma4']) + artale_tables[population]['redshift'], + artale_tables[population]['gamma4']) return _m_star_sfr_metallicity_merger_rate(m_star, sfr, Z, gamma1, gamma2, gamma3, gamma4) + def _m_star_merger_rate(m_star, - alpha1, - alpha2): + alpha1, + alpha2): m_star = m_star.to(units.Msun).value @@ -212,11 +215,12 @@ def _m_star_merger_rate(m_star, return n_gw / units.Gyr + def _m_star_sfr_merger_rate(m_star, - sfr, - beta1, - beta2, - beta3): + sfr, + beta1, + beta2, + beta3): m_star = m_star.to(units.Msun).value sfr = sfr.to(units.Msun / units.year) @@ -225,13 +229,14 @@ def _m_star_sfr_merger_rate(m_star, return n_gw / units.Gyr + def _m_star_sfr_metallicity_merger_rate(m_star, - sfr, - Z, - gamma1, - gamma2, - gamma3, - gamma4): + sfr, + Z, + gamma1, + gamma2, + gamma3, + gamma4): m_star = m_star.to(units.Msun).value sfr = sfr.to(units.Msun / units.year) @@ -240,6 +245,7 @@ def _m_star_sfr_metallicity_merger_rate(m_star, return n_gw / units.Gyr + def b_band_merger_rate(luminosity, population='NS-NS', optimism='low'): diff --git a/skypy/gravitational_waves/tests/test_merger_rate.py b/skypy/gravitational_waves/tests/test_merger_rate.py index 42e298d49..f237de267 100644 --- a/skypy/gravitational_waves/tests/test_merger_rate.py +++ b/skypy/gravitational_waves/tests/test_merger_rate.py @@ -18,19 +18,20 @@ def test_abadie_rates(): table_value = abadie_table_III['NS-NS']['low'] assert np.isclose(L_B_rate.to_value(1/units.year), table_value, rtol=1e-5) + def test_artale_rates(): # Check the number of merger rates returned redshifts = np.random.uniform(0., 3., 100) stellar_masses = 10.**(9.0 + np.random.randn(100)) - rates = m_star_merger_rate(redshifts, stellar_masses, population='NS-NS') + rates = m_star_merger_rate(redshifts, stellar_masses * units.Msun, population='NS-NS') assert len(rates) == len(stellar_masses) # Test that a luminosity of M_sol returns a # rate that matches the value in Artale Table I - m_sol = constants.Msun.to_value('Msun') - m_sol_rate = m_star_merger_rate(0.0, m_sol, population='NS-NS') + m_sol = constants.M_sun.to_value('M_sun') + m_sol_rate = m_star_merger_rate(0.0, m_sol * units.Msun, population='NS-NS') alpha1 = artale_tables['NS-NS']['alpha1'][0] alpha2 = artale_tables['NS-NS']['alpha2'][0] table_value = 10.**(alpha1 * np.log10(m_sol) + alpha2) - assert np.isclose(m_sol_rate.to_value(1/units.Gyr), table_value, rtol=1e-5) \ No newline at end of file + assert np.isclose(m_sol_rate.to_value(1/units.Gyr), table_value, rtol=1e-5) From 1bd53ed3d0f22fccc061beb0b28a4a40424dce95 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Mon, 19 Jul 2021 14:27:16 +0100 Subject: [PATCH 03/22] fixed doc example --- skypy/gravitational_waves/merger_rate.py | 1 + 1 file changed, 1 insertion(+) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index d7fd30e8e..0d01d6ba5 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -146,6 +146,7 @@ def m_star_merger_rate(redshift, Sample 100 stellar masses values near 10^9 solar masses. >>> stellar_masses = 10.**(9.0 + np.random.randn(100)) + >>> redshifts = np.random.uniform(0., 3., 100) Generate merger rates for these luminosities. From e8437f90073326631677644f3fedac750fcb21a4 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Mon, 19 Jul 2021 14:31:14 +0100 Subject: [PATCH 04/22] fixed doc example, MSun units --- skypy/gravitational_waves/merger_rate.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 0d01d6ba5..a3db379bc 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -141,6 +141,7 @@ def m_star_merger_rate(redshift, Examples -------- >>> import numpy as np + >>> from astropy import units >>> from skypy.gravitational_waves import m_star_merger_rate Sample 100 stellar masses values near 10^9 solar masses. @@ -151,7 +152,7 @@ def m_star_merger_rate(redshift, Generate merger rates for these luminosities. >>> rates = m_star_merger_rate(redshifts, - ... stellar_masses, + ... stellar_masses * units.Msun, ... population='NS-NS') """ From c2ef42eea40d4ab6ab6d138a834acec47290d64e Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Mon, 19 Jul 2021 16:14:31 +0100 Subject: [PATCH 05/22] added sfr and metallicity model docs and tests --- setup.cfg | 2 +- skypy/gravitational_waves/__init__.py | 1 + skypy/gravitational_waves/merger_rate.py | 153 +++++++++++++++++- .../tests/test_merger_rate.py | 41 ++++- 4 files changed, 189 insertions(+), 8 deletions(-) diff --git a/setup.cfg b/setup.cfg index a40168ad4..8b2875ea1 100644 --- a/setup.cfg +++ b/setup.cfg @@ -48,7 +48,7 @@ testpaths = "skypy" "docs" astropy_header = true doctest_plus = enabled text_file_format = rst -addopts = --doctest-rst +#addopts = --doctest-rst [coverage:run] omit = diff --git a/skypy/gravitational_waves/__init__.py b/skypy/gravitational_waves/__init__.py index 509f18bd9..798d64a47 100644 --- a/skypy/gravitational_waves/__init__.py +++ b/skypy/gravitational_waves/__init__.py @@ -11,6 +11,7 @@ b_band_merger_rate m_star_merger_rate + m_star_sfr_merger_rate """ diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index a3db379bc..5bebe66bc 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -11,7 +11,9 @@ __all__ = [ 'b_band_merger_rate', - 'm_star_merger_rate' + 'm_star_merger_rate', + 'm_star_sfr_merger_rate', + 'm_star_sfr_metallicity_merger_rate' ] @@ -144,10 +146,13 @@ def m_star_merger_rate(redshift, >>> from astropy import units >>> from skypy.gravitational_waves import m_star_merger_rate + Sample 100 redshifts. + + >>> redshifts = np.random.uniform(0., 3., 100) + Sample 100 stellar masses values near 10^9 solar masses. >>> stellar_masses = 10.**(9.0 + np.random.randn(100)) - >>> redshifts = np.random.uniform(0., 3., 100) Generate merger rates for these luminosities. @@ -171,6 +176,72 @@ def m_star_sfr_merger_rate(redshift, m_star, sfr, population): + r"""Model of Artale et al (2020), equation (2) with parameters + from Tables I, II and III. + + Compact binary merger rates as a log-log function of a galaxy's + stellar mass. + + Parameters are redshift dependent, with linear interpolation + between the simulated points of z={0.1, 1, 2, 6}. + + Parameters + ---------- + redshift : (ngal,) array-like + The redshifts of the galaxies to generate merger + rates for. + m_star : (ngal,) array-like + The stellar mass of the galaxies to generate merger + rates for, in units of stellar mass. + sfr : (ngal,) array-like + The star formation rate of the galaxies to generate + merger rates for, in units of stellar mass per year + population : {'NS-NS', 'NS-BH', 'BH-BH'} + Compact binary population to get rate for. + 'NS-NS' is neutron star - neutron star + 'NS-BH' is neutron star - black hole + 'BH-BH' is black hole - black hole + + Returns + ------- + merger_rate : array_like + Merger rates for the galaxies in units of Gigayear^-1 + + Notes + ----- + + References + ---------- + .. Artale et al. 2020, MNRAS, + Volume 491, Issue 3, p.3419-3434 (2020) + https://arxiv.org/abs/1910.04890 + + Examples + -------- + >>> import numpy as np + >>> from astropy import units + >>> from skypy.gravitational_waves import m_star_merger_rate + + Sample 100 redshifts. + + >>> redshifts = np.random.uniform(0., 3., 100) + + Sample 100 stellar masses values near 10^9 solar masses. + + >>> stellar_masses = 10.**(9.0 + np.random.randn(100)) + + Sample 100 star formation rates. + + >>> sfrs = 10.**(np.random.randn(100)) + + Generate merger rates for these luminosities. + + >>> rates = m_star_merger_rate(redshifts, + ... stellar_masses * units.Msun, + ... sfrs * units.Msun / units.year, + ... population='NS-NS') + + """ beta1 = np.interp(redshift, artale_tables[population]['redshift'], @@ -190,6 +261,80 @@ def m_star_sfr_metallicity_merger_rate(redshift, sfr, Z, population): + r"""Model of Artale et al (2020), equation (3) with parameters + from Tables I, II and III. + + Compact binary merger rates as a log-log function of a galaxy's + stellar mass. + + Parameters are redshift dependent, with linear interpolation + between the simulated points of z={0.1, 1, 2, 6}. + + Parameters + ---------- + redshift : (ngal,) array-like + The redshifts of the galaxies to generate merger + rates for. + m_star : (ngal,) array-like + The stellar mass of the galaxies to generate merger + rates for, in units of stellar mass. + sfr : (ngal,) array-like + The star formation rate of the galaxies to generate + merger rates for, in units of stellar mass per year + Z : (ngal,) array-like + The metallicity of the galaxies to generate merger + rates for. + population : {'NS-NS', 'NS-BH', 'BH-BH'} + Compact binary population to get rate for. + 'NS-NS' is neutron star - neutron star + 'NS-BH' is neutron star - black hole + 'BH-BH' is black hole - black hole + + Returns + ------- + merger_rate : array_like + Merger rates for the galaxies in units of Gigayear^-1 + + Notes + ----- + + References + ---------- + .. Artale et al. 2020, MNRAS, + Volume 491, Issue 3, p.3419-3434 (2020) + https://arxiv.org/abs/1910.04890 + + Examples + -------- + >>> import numpy as np + >>> from astropy import units + >>> from skypy.gravitational_waves import m_star_sfr_metallicity_merger_rate + + Sample 100 redshifts. + + >>> redshifts = np.random.uniform(0., 3., 100) + + Sample 100 stellar masses values near 10^9 solar masses. + + >>> stellar_masses = 10.**(9.0 + np.random.randn(100)) + + Sample 100 star formation rates. + + >>> sfrs = 10.**(np.random.randn(100)) + + Sample 100 metallicities. + + >>> metallicities = 10.**(-2.0 + np.random.randn(100)) + + Generate merger rates for these luminosities. + + >>> rates = m_star_sfr_metallicity_merger_rate(redshifts, + ... stellar_masses * units.Msun, + ... sfrs * units.Msun / units.year, + ... metallicities, + ... population='NS-NS') + + """ gamma1 = np.interp(redshift, artale_tables[population]['redshift'], @@ -225,7 +370,7 @@ def _m_star_sfr_merger_rate(m_star, beta3): m_star = m_star.to(units.Msun).value - sfr = sfr.to(units.Msun / units.year) + sfr = sfr.to(units.Msun / units.year).value n_gw = 10.**(beta1 * np.log10(m_star) + beta2 * np.log10(sfr) + beta3) @@ -241,7 +386,7 @@ def _m_star_sfr_metallicity_merger_rate(m_star, gamma4): m_star = m_star.to(units.Msun).value - sfr = sfr.to(units.Msun / units.year) + sfr = sfr.to(units.Msun / units.year).value n_gw = 10.**(gamma1 * np.log10(m_star) + gamma2 * np.log10(sfr) + gamma3 * np.log10(Z) + gamma4) diff --git a/skypy/gravitational_waves/tests/test_merger_rate.py b/skypy/gravitational_waves/tests/test_merger_rate.py index f237de267..f774979bc 100644 --- a/skypy/gravitational_waves/tests/test_merger_rate.py +++ b/skypy/gravitational_waves/tests/test_merger_rate.py @@ -1,6 +1,10 @@ import numpy as np from astropy import constants, units -from skypy.gravitational_waves import b_band_merger_rate, m_star_merger_rate +from skypy.gravitational_waves import \ + b_band_merger_rate, \ + m_star_merger_rate, \ + m_star_sfr_merger_rate, \ + m_star_sfr_metallicity_merger_rate from skypy.gravitational_waves.merger_rate import abadie_table_III, artale_tables @@ -30,8 +34,39 @@ def test_artale_rates(): # Test that a luminosity of M_sol returns a # rate that matches the value in Artale Table I m_sol = constants.M_sun.to_value('M_sun') - m_sol_rate = m_star_merger_rate(0.0, m_sol * units.Msun, population='NS-NS') + sol_rate = m_star_merger_rate(0.0, m_sol * units.Msun, population='NS-NS') alpha1 = artale_tables['NS-NS']['alpha1'][0] alpha2 = artale_tables['NS-NS']['alpha2'][0] table_value = 10.**(alpha1 * np.log10(m_sol) + alpha2) - assert np.isclose(m_sol_rate.to_value(1/units.Gyr), table_value, rtol=1e-5) + assert np.isclose(sol_rate.to_value(1/units.Gyr), table_value, rtol=1e-5) + + # Test that a luminosity of M_sol and SFR of 1 returns a + # rate that matches the value in Artale Table I + m_sol = constants.M_sun.to_value('M_sun') + sfr_1 = 1. + sol_rate = m_star_sfr_merger_rate(0.0, m_sol * units.Msun, + sfr_1 * units.Msun / units.year, + population='NS-NS') + beta1 = artale_tables['NS-NS']['beta1'][0] + beta2 = artale_tables['NS-NS']['beta2'][0] + beta3 = artale_tables['NS-NS']['beta3'][0] + table_value = 10.**(beta1 * np.log10(m_sol) + beta2 * np.log10(sfr_1) + beta3) + assert np.isclose(sol_rate.to_value(1/units.Gyr), table_value, rtol=1e-5) + + # Test that a luminosity of M_sol, SFR of 1 and Z of 0.012 returns a + # rate that matches the value in Artale Table I + m_sol = constants.M_sun.to_value('M_sun') + sfr_1 = 1. + Z_sol = 0.012 + sol_rate = m_star_sfr_metallicity_merger_rate(0.0, m_sol * units.Msun, + sfr_1 * units.Msun / units.year, + Z_sol, + population='NS-NS') + gamma1 = artale_tables['NS-NS']['gamma1'][0] + gamma2 = artale_tables['NS-NS']['gamma2'][0] + gamma3 = artale_tables['NS-NS']['gamma3'][0] + gamma4 = artale_tables['NS-NS']['gamma4'][0] + table_value = 10.**(gamma1 * np.log10(m_sol) + + gamma2 * np.log10(sfr_1) + + gamma3 * np.log10(Z_sol) + gamma4) + assert np.isclose(sol_rate.to_value(1/units.Gyr), table_value, rtol=1e-5) From 8d23946669df387fa0eac5f74c22a5f315880bae Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Mon, 19 Jul 2021 16:14:57 +0100 Subject: [PATCH 06/22] put back doctest --- setup.cfg | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/setup.cfg b/setup.cfg index 8b2875ea1..a40168ad4 100644 --- a/setup.cfg +++ b/setup.cfg @@ -48,7 +48,7 @@ testpaths = "skypy" "docs" astropy_header = true doctest_plus = enabled text_file_format = rst -#addopts = --doctest-rst +addopts = --doctest-rst [coverage:run] omit = From 787720710ec486640541fa816844e9038850bde4 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Wed, 21 Jul 2021 08:55:29 +0100 Subject: [PATCH 07/22] fixed sfr docstring --- skypy/gravitational_waves/merger_rate.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 5bebe66bc..3ac6f5eeb 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -220,7 +220,7 @@ def m_star_sfr_merger_rate(redshift, -------- >>> import numpy as np >>> from astropy import units - >>> from skypy.gravitational_waves import m_star_merger_rate + >>> from skypy.gravitational_waves import m_star_sfr_merger_rate Sample 100 redshifts. @@ -236,10 +236,10 @@ def m_star_sfr_merger_rate(redshift, Generate merger rates for these luminosities. - >>> rates = m_star_merger_rate(redshifts, - ... stellar_masses * units.Msun, - ... sfrs * units.Msun / units.year, - ... population='NS-NS') + >>> rates = m_star_sfr_merger_rate(redshifts, + ... stellar_masses * units.Msun, + ... sfrs * units.Msun / units.year, + ... population='NS-NS') """ From cc259b9903b4ac70fb5e4b7b03bbc1d17a161818 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Wed, 21 Jul 2021 10:08:15 +0100 Subject: [PATCH 08/22] added metallicity function to docs autosummary --- skypy/gravitational_waves/__init__.py | 1 + 1 file changed, 1 insertion(+) diff --git a/skypy/gravitational_waves/__init__.py b/skypy/gravitational_waves/__init__.py index 798d64a47..4fb1fca80 100644 --- a/skypy/gravitational_waves/__init__.py +++ b/skypy/gravitational_waves/__init__.py @@ -12,6 +12,7 @@ b_band_merger_rate m_star_merger_rate m_star_sfr_merger_rate + m_star_sfr_metallicity_merger_rate """ From 5cfcf4eaaf62865490e399e995d1c0e11337ad71 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:26:50 +0100 Subject: [PATCH 09/22] Update skypy/gravitational_waves/merger_rate.py Co-authored-by: Lucia F. de la Bella <55983939+Lucia-Fonseca@users.noreply.github.com> --- skypy/gravitational_waves/merger_rate.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 3ac6f5eeb..3202359b9 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -103,7 +103,8 @@ def m_star_merger_rate(redshift, m_star, population): - r"""Model of Artale et al (2020), equation (1) with parameters + r"""M_star merger rate. + Model of Artale et al (2020), equation (1) with parameters from Tables I, II and III. Compact binary merger rates as a log-log function of a galaxy's From de63b1c7689d1369578dc81fdcbdf05ca282f509 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:27:03 +0100 Subject: [PATCH 10/22] Update skypy/gravitational_waves/merger_rate.py Co-authored-by: Lucia F. de la Bella <55983939+Lucia-Fonseca@users.noreply.github.com> --- skypy/gravitational_waves/merger_rate.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 3202359b9..41c5cd08e 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -118,7 +118,7 @@ def m_star_merger_rate(redshift, redshift : (ngal,) array-like The redshifts of the galaxies to generate merger rates for. - m_star : (ngal,) array-like + m_star : (ngal,) array_like The stellar mass of the galaxies to generate merger rates for, in units of stellar mass. population : {'NS-NS', 'NS-BH', 'BH-BH'} From 94d16d3791cb9c5be86949fa4df08747031c2944 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:27:19 +0100 Subject: [PATCH 11/22] Update skypy/gravitational_waves/tests/test_merger_rate.py Co-authored-by: philipp128 <48715661+philipp128@users.noreply.github.com> --- skypy/gravitational_waves/tests/test_merger_rate.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skypy/gravitational_waves/tests/test_merger_rate.py b/skypy/gravitational_waves/tests/test_merger_rate.py index f774979bc..201d3b6d3 100644 --- a/skypy/gravitational_waves/tests/test_merger_rate.py +++ b/skypy/gravitational_waves/tests/test_merger_rate.py @@ -31,7 +31,7 @@ def test_artale_rates(): rates = m_star_merger_rate(redshifts, stellar_masses * units.Msun, population='NS-NS') assert len(rates) == len(stellar_masses) - # Test that a luminosity of M_sol returns a + # Test that a luminosity of M_sol returns a # rate that matches the value in Artale Table I m_sol = constants.M_sun.to_value('M_sun') sol_rate = m_star_merger_rate(0.0, m_sol * units.Msun, population='NS-NS') From 31e69f88b623069588a36a8319046c61f74cbf41 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:27:32 +0100 Subject: [PATCH 12/22] Update skypy/gravitational_waves/merger_rate.py Co-authored-by: Lucia F. de la Bella <55983939+Lucia-Fonseca@users.noreply.github.com> --- skypy/gravitational_waves/merger_rate.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 41c5cd08e..bb921a9bb 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -129,7 +129,7 @@ def m_star_merger_rate(redshift, Returns ------- - merger_rate : array_like + merger_rate : (ngal,) array_like Merger rates for the galaxies in units of Gigayear^-1 Notes From e8afd1d7c921a97c42701b60573b7ba04c9a914b Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:27:41 +0100 Subject: [PATCH 13/22] Update skypy/gravitational_waves/merger_rate.py Co-authored-by: Lucia F. de la Bella <55983939+Lucia-Fonseca@users.noreply.github.com> --- skypy/gravitational_waves/merger_rate.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index bb921a9bb..8765a358b 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -115,7 +115,7 @@ def m_star_merger_rate(redshift, Parameters ---------- - redshift : (ngal,) array-like + redshift : (ngal,) array_like The redshifts of the galaxies to generate merger rates for. m_star : (ngal,) array_like From c713d81c6a0c93060725edc6e308281fc153a29b Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:27:48 +0100 Subject: [PATCH 14/22] Update skypy/gravitational_waves/merger_rate.py Co-authored-by: Lucia F. de la Bella <55983939+Lucia-Fonseca@users.noreply.github.com> --- skypy/gravitational_waves/merger_rate.py | 20 +------------------- 1 file changed, 1 insertion(+), 19 deletions(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 8765a358b..26cfb59ba 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -141,25 +141,7 @@ def m_star_merger_rate(redshift, Volume 491, Issue 3, p.3419-3434 (2020) https://arxiv.org/abs/1910.04890 - Examples - -------- - >>> import numpy as np - >>> from astropy import units - >>> from skypy.gravitational_waves import m_star_merger_rate - - Sample 100 redshifts. - - >>> redshifts = np.random.uniform(0., 3., 100) - - Sample 100 stellar masses values near 10^9 solar masses. - - >>> stellar_masses = 10.**(9.0 + np.random.randn(100)) - - Generate merger rates for these luminosities. - - >>> rates = m_star_merger_rate(redshifts, - ... stellar_masses * units.Msun, - ... population='NS-NS') + """ From bfac3de77582d3e2c8bc84a0dc8b5a1f3d05b686 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:27:56 +0100 Subject: [PATCH 15/22] Update skypy/gravitational_waves/merger_rate.py Co-authored-by: Lucia F. de la Bella <55983939+Lucia-Fonseca@users.noreply.github.com> --- skypy/gravitational_waves/merger_rate.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 26cfb59ba..4110d4538 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -255,7 +255,7 @@ def m_star_sfr_metallicity_merger_rate(redshift, Parameters ---------- - redshift : (ngal,) array-like + redshift : (ngal,) array_like The redshifts of the galaxies to generate merger rates for. m_star : (ngal,) array-like From 246420209d9df43c2218e53bcc18d125a32de160 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:28:08 +0100 Subject: [PATCH 16/22] Update skypy/gravitational_waves/merger_rate.py Co-authored-by: Lucia F. de la Bella <55983939+Lucia-Fonseca@users.noreply.github.com> --- skypy/gravitational_waves/merger_rate.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 4110d4538..5a9d27be7 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -261,7 +261,7 @@ def m_star_sfr_metallicity_merger_rate(redshift, m_star : (ngal,) array-like The stellar mass of the galaxies to generate merger rates for, in units of stellar mass. - sfr : (ngal,) array-like + sfr : (ngal,) array_like The star formation rate of the galaxies to generate merger rates for, in units of stellar mass per year Z : (ngal,) array-like From 97cab80d16c9c96370f228e881163b9947db781d Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:28:18 +0100 Subject: [PATCH 17/22] Update skypy/gravitational_waves/merger_rate.py Co-authored-by: Lucia F. de la Bella <55983939+Lucia-Fonseca@users.noreply.github.com> --- skypy/gravitational_waves/merger_rate.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 5a9d27be7..8a908ec42 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -264,7 +264,7 @@ def m_star_sfr_metallicity_merger_rate(redshift, sfr : (ngal,) array_like The star formation rate of the galaxies to generate merger rates for, in units of stellar mass per year - Z : (ngal,) array-like + Z : (ngal,) array_like The metallicity of the galaxies to generate merger rates for. population : {'NS-NS', 'NS-BH', 'BH-BH'} From c322ffaa24e8ac33afcc82f00e68ba38344fd277 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:28:25 +0100 Subject: [PATCH 18/22] Update skypy/gravitational_waves/merger_rate.py Co-authored-by: Lucia F. de la Bella <55983939+Lucia-Fonseca@users.noreply.github.com> --- skypy/gravitational_waves/merger_rate.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 8a908ec42..8a55b5f14 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -275,7 +275,7 @@ def m_star_sfr_metallicity_merger_rate(redshift, Returns ------- - merger_rate : array_like + merger_rate : (ngal,) array_like Merger rates for the galaxies in units of Gigayear^-1 Notes From 21b2f19517619c9382a8c733a1568ac4954cbb00 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:28:34 +0100 Subject: [PATCH 19/22] Update skypy/gravitational_waves/merger_rate.py Co-authored-by: Lucia F. de la Bella <55983939+Lucia-Fonseca@users.noreply.github.com> --- skypy/gravitational_waves/merger_rate.py | 30 +----------------------- 1 file changed, 1 insertion(+), 29 deletions(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 8a55b5f14..3196082f6 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -287,35 +287,7 @@ def m_star_sfr_metallicity_merger_rate(redshift, Volume 491, Issue 3, p.3419-3434 (2020) https://arxiv.org/abs/1910.04890 - Examples - -------- - >>> import numpy as np - >>> from astropy import units - >>> from skypy.gravitational_waves import m_star_sfr_metallicity_merger_rate - - Sample 100 redshifts. - - >>> redshifts = np.random.uniform(0., 3., 100) - - Sample 100 stellar masses values near 10^9 solar masses. - - >>> stellar_masses = 10.**(9.0 + np.random.randn(100)) - - Sample 100 star formation rates. - - >>> sfrs = 10.**(np.random.randn(100)) - - Sample 100 metallicities. - - >>> metallicities = 10.**(-2.0 + np.random.randn(100)) - - Generate merger rates for these luminosities. - - >>> rates = m_star_sfr_metallicity_merger_rate(redshifts, - ... stellar_masses * units.Msun, - ... sfrs * units.Msun / units.year, - ... metallicities, - ... population='NS-NS') + """ From 15f951dbfa1e0e9b5561f0edbfd610e480176143 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:50:15 +0100 Subject: [PATCH 20/22] removed private functions, doc examples --- skypy/gravitational_waves/merger_rate.py | 99 +++---------------- .../tests/test_merger_rate.py | 19 +++- 2 files changed, 30 insertions(+), 88 deletions(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 3196082f6..af7aa0c27 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -107,7 +107,7 @@ def m_star_merger_rate(redshift, Model of Artale et al (2020), equation (1) with parameters from Tables I, II and III. - Compact binary merger rates as a log-log function of a galaxy's + Compact binary merger rates as a power law function of a galaxy's stellar mass. Parameters are redshift dependent, with linear interpolation @@ -141,8 +141,6 @@ def m_star_merger_rate(redshift, Volume 491, Issue 3, p.3419-3434 (2020) https://arxiv.org/abs/1910.04890 - - """ alpha1 = np.interp(redshift, @@ -152,7 +150,11 @@ def m_star_merger_rate(redshift, artale_tables[population]['redshift'], artale_tables[population]['alpha2']) - return _m_star_merger_rate(m_star, alpha1, alpha2) + m_star = m_star.to(units.Msun).value + + n_gw = 10.**(alpha1 * np.log10(m_star) + alpha2) + + return n_gw / units.Gyr def m_star_sfr_merger_rate(redshift, @@ -162,7 +164,7 @@ def m_star_sfr_merger_rate(redshift, r"""Model of Artale et al (2020), equation (2) with parameters from Tables I, II and III. - Compact binary merger rates as a log-log function of a galaxy's + Compact binary merger rates as a power law function of a galaxy's stellar mass. Parameters are redshift dependent, with linear interpolation @@ -199,31 +201,6 @@ def m_star_sfr_merger_rate(redshift, Volume 491, Issue 3, p.3419-3434 (2020) https://arxiv.org/abs/1910.04890 - Examples - -------- - >>> import numpy as np - >>> from astropy import units - >>> from skypy.gravitational_waves import m_star_sfr_merger_rate - - Sample 100 redshifts. - - >>> redshifts = np.random.uniform(0., 3., 100) - - Sample 100 stellar masses values near 10^9 solar masses. - - >>> stellar_masses = 10.**(9.0 + np.random.randn(100)) - - Sample 100 star formation rates. - - >>> sfrs = 10.**(np.random.randn(100)) - - Generate merger rates for these luminosities. - - >>> rates = m_star_sfr_merger_rate(redshifts, - ... stellar_masses * units.Msun, - ... sfrs * units.Msun / units.year, - ... population='NS-NS') - """ beta1 = np.interp(redshift, @@ -236,7 +213,12 @@ def m_star_sfr_merger_rate(redshift, artale_tables[population]['redshift'], artale_tables[population]['beta3']) - return _m_star_sfr_merger_rate(m_star, sfr, beta1, beta2, beta3) + m_star = m_star.to(units.Msun).value + sfr = sfr.to(units.Msun / units.year).value + + n_gw = 10.**(beta1 * np.log10(m_star) + beta2 * np.log10(sfr) + beta3) + + return n_gw / units.Gyr def m_star_sfr_metallicity_merger_rate(redshift, @@ -247,7 +229,7 @@ def m_star_sfr_metallicity_merger_rate(redshift, r"""Model of Artale et al (2020), equation (3) with parameters from Tables I, II and III. - Compact binary merger rates as a log-log function of a galaxy's + Compact binary merger rates as a power law function of a galaxy's stellar mass. Parameters are redshift dependent, with linear interpolation @@ -287,8 +269,6 @@ def m_star_sfr_metallicity_merger_rate(redshift, Volume 491, Issue 3, p.3419-3434 (2020) https://arxiv.org/abs/1910.04890 - - """ gamma1 = np.interp(redshift, @@ -304,42 +284,6 @@ def m_star_sfr_metallicity_merger_rate(redshift, artale_tables[population]['redshift'], artale_tables[population]['gamma4']) - return _m_star_sfr_metallicity_merger_rate(m_star, sfr, Z, gamma1, gamma2, gamma3, gamma4) - - -def _m_star_merger_rate(m_star, - alpha1, - alpha2): - - m_star = m_star.to(units.Msun).value - - n_gw = 10.**(alpha1 * np.log10(m_star) + alpha2) - - return n_gw / units.Gyr - - -def _m_star_sfr_merger_rate(m_star, - sfr, - beta1, - beta2, - beta3): - - m_star = m_star.to(units.Msun).value - sfr = sfr.to(units.Msun / units.year).value - - n_gw = 10.**(beta1 * np.log10(m_star) + beta2 * np.log10(sfr) + beta3) - - return n_gw / units.Gyr - - -def _m_star_sfr_metallicity_merger_rate(m_star, - sfr, - Z, - gamma1, - gamma2, - gamma3, - gamma4): - m_star = m_star.to(units.Msun).value sfr = sfr.to(units.Msun / units.year).value @@ -385,21 +329,6 @@ def b_band_merger_rate(luminosity, Volume 27, Issue 17, article id. 173001 (2010) https://arxiv.org/abs/1003.2480 - Examples - -------- - >>> import numpy as np - >>> from skypy.gravitational_waves import b_band_merger_rate - - Sample 100 luminosity values near absolute magnitude -20.5. - - >>> luminosities = 10.**(-0.4*(-20.5 + np.random.randn(100))) - - Generate merger rates for these luminosities. - - >>> rates = b_band_merger_rate(luminosities, - ... population='NS-NS', - ... optimism='low') - """ # Convert luminosity to units of L_10 defined in Abadie et. al. 2010 diff --git a/skypy/gravitational_waves/tests/test_merger_rate.py b/skypy/gravitational_waves/tests/test_merger_rate.py index 201d3b6d3..07bc8abf7 100644 --- a/skypy/gravitational_waves/tests/test_merger_rate.py +++ b/skypy/gravitational_waves/tests/test_merger_rate.py @@ -28,10 +28,23 @@ def test_artale_rates(): # Check the number of merger rates returned redshifts = np.random.uniform(0., 3., 100) stellar_masses = 10.**(9.0 + np.random.randn(100)) + sfrs = 10.**(np.random.randn(100)) + metallicities = 10.**(-2.0 + np.random.randn(100)) rates = m_star_merger_rate(redshifts, stellar_masses * units.Msun, population='NS-NS') assert len(rates) == len(stellar_masses) + rates = m_star_sfr_merger_rate(redshifts, + stellar_masses * units.Msun, + sfrs * units.Msun / units.year, + population='NS-NS') + assert len(rates) == len(stellar_masses) + rates = m_star_sfr_metallicity_merger_rate(redshifts, + stellar_masses * units.Msun, + sfrs * units.Msun / units.year, + metallicities, + population='NS-NS') + assert len(rates) == len(stellar_masses) - # Test that a luminosity of M_sol returns a + # Test that a stellar mass of M_sol returns a # rate that matches the value in Artale Table I m_sol = constants.M_sun.to_value('M_sun') sol_rate = m_star_merger_rate(0.0, m_sol * units.Msun, population='NS-NS') @@ -40,7 +53,7 @@ def test_artale_rates(): table_value = 10.**(alpha1 * np.log10(m_sol) + alpha2) assert np.isclose(sol_rate.to_value(1/units.Gyr), table_value, rtol=1e-5) - # Test that a luminosity of M_sol and SFR of 1 returns a + # Test that a stellar mass of M_sol and SFR of 1 returns a # rate that matches the value in Artale Table I m_sol = constants.M_sun.to_value('M_sun') sfr_1 = 1. @@ -53,7 +66,7 @@ def test_artale_rates(): table_value = 10.**(beta1 * np.log10(m_sol) + beta2 * np.log10(sfr_1) + beta3) assert np.isclose(sol_rate.to_value(1/units.Gyr), table_value, rtol=1e-5) - # Test that a luminosity of M_sol, SFR of 1 and Z of 0.012 returns a + # Test that a stellar mass of M_sol, SFR of 1 and Z of 0.012 returns a # rate that matches the value in Artale Table I m_sol = constants.M_sun.to_value('M_sun') sfr_1 = 1. From b7e15f9c068ce561bb3759630bc303eb894226c0 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Thu, 22 Jul 2021 15:54:41 +0100 Subject: [PATCH 21/22] couple fo small typos --- skypy/gravitational_waves/merger_rate.py | 16 +++++++++------- 1 file changed, 9 insertions(+), 7 deletions(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index af7aa0c27..23a953522 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -161,7 +161,8 @@ def m_star_sfr_merger_rate(redshift, m_star, sfr, population): - r"""Model of Artale et al (2020), equation (2) with parameters + r"""M_star SFR merger rate. + Model of Artale et al (2020), equation (2) with parameters from Tables I, II and III. Compact binary merger rates as a power law function of a galaxy's @@ -172,13 +173,13 @@ def m_star_sfr_merger_rate(redshift, Parameters ---------- - redshift : (ngal,) array-like + redshift : (ngal,) array_like The redshifts of the galaxies to generate merger rates for. - m_star : (ngal,) array-like + m_star : (ngal,) array_like The stellar mass of the galaxies to generate merger rates for, in units of stellar mass. - sfr : (ngal,) array-like + sfr : (ngal,) array_like The star formation rate of the galaxies to generate merger rates for, in units of stellar mass per year population : {'NS-NS', 'NS-BH', 'BH-BH'} @@ -226,7 +227,8 @@ def m_star_sfr_metallicity_merger_rate(redshift, sfr, Z, population): - r"""Model of Artale et al (2020), equation (3) with parameters + r"""M_star SFR Metallicity merger rate. + Model of Artale et al (2020), equation (3) with parameters from Tables I, II and III. Compact binary merger rates as a power law function of a galaxy's @@ -240,7 +242,7 @@ def m_star_sfr_metallicity_merger_rate(redshift, redshift : (ngal,) array_like The redshifts of the galaxies to generate merger rates for. - m_star : (ngal,) array-like + m_star : (ngal,) array_like The stellar mass of the galaxies to generate merger rates for, in units of stellar mass. sfr : (ngal,) array_like @@ -303,7 +305,7 @@ def b_band_merger_rate(luminosity, Parameters ---------- - luminosity : (ngal,) array-like + luminosity : (ngal,) array_like The B-band luminosity of the galaxies to generate merger rates for, in units of solar luminosity. population : {'NS-NS', 'NS-BH', 'BH-BH'} From cd58d7621818a1b7e1d014f4d1469284f87ae014 Mon Sep 17 00:00:00 2001 From: Ian Harrison Date: Fri, 19 Nov 2021 10:33:28 +0000 Subject: [PATCH 22/22] changed docs to specify astropy.Quantity --- skypy/gravitational_waves/merger_rate.py | 20 ++++++++++---------- 1 file changed, 10 insertions(+), 10 deletions(-) diff --git a/skypy/gravitational_waves/merger_rate.py b/skypy/gravitational_waves/merger_rate.py index 23a953522..18cdf3b20 100644 --- a/skypy/gravitational_waves/merger_rate.py +++ b/skypy/gravitational_waves/merger_rate.py @@ -118,7 +118,7 @@ def m_star_merger_rate(redshift, redshift : (ngal,) array_like The redshifts of the galaxies to generate merger rates for. - m_star : (ngal,) array_like + m_star : (ngal,) astropy.Quantity The stellar mass of the galaxies to generate merger rates for, in units of stellar mass. population : {'NS-NS', 'NS-BH', 'BH-BH'} @@ -129,7 +129,7 @@ def m_star_merger_rate(redshift, Returns ------- - merger_rate : (ngal,) array_like + merger_rate : (ngal,) astropy.Quantity Merger rates for the galaxies in units of Gigayear^-1 Notes @@ -176,10 +176,10 @@ def m_star_sfr_merger_rate(redshift, redshift : (ngal,) array_like The redshifts of the galaxies to generate merger rates for. - m_star : (ngal,) array_like + m_star : (ngal,) astropy.Quantity The stellar mass of the galaxies to generate merger rates for, in units of stellar mass. - sfr : (ngal,) array_like + sfr : (ngal,) astropy.Quantity The star formation rate of the galaxies to generate merger rates for, in units of stellar mass per year population : {'NS-NS', 'NS-BH', 'BH-BH'} @@ -190,7 +190,7 @@ def m_star_sfr_merger_rate(redshift, Returns ------- - merger_rate : array_like + merger_rate : astropy.Quantity Merger rates for the galaxies in units of Gigayear^-1 Notes @@ -242,10 +242,10 @@ def m_star_sfr_metallicity_merger_rate(redshift, redshift : (ngal,) array_like The redshifts of the galaxies to generate merger rates for. - m_star : (ngal,) array_like + m_star : (ngal,) astropy.Quantity The stellar mass of the galaxies to generate merger rates for, in units of stellar mass. - sfr : (ngal,) array_like + sfr : (ngal,) astropy.Quantity The star formation rate of the galaxies to generate merger rates for, in units of stellar mass per year Z : (ngal,) array_like @@ -259,7 +259,7 @@ def m_star_sfr_metallicity_merger_rate(redshift, Returns ------- - merger_rate : (ngal,) array_like + merger_rate : (ngal,) astropy.Quantity Merger rates for the galaxies in units of Gigayear^-1 Notes @@ -305,7 +305,7 @@ def b_band_merger_rate(luminosity, Parameters ---------- - luminosity : (ngal,) array_like + luminosity : (ngal,) astropy.Quantity The B-band luminosity of the galaxies to generate merger rates for, in units of solar luminosity. population : {'NS-NS', 'NS-BH', 'BH-BH'} @@ -319,7 +319,7 @@ def b_band_merger_rate(luminosity, Returns ------- - merger_rate : array_like + merger_rate : astropy.Quantity Merger rates for the galaxies in units of year^-1 Notes