From a8226f980c4392ecd74cf0f84a35332e34d890ac Mon Sep 17 00:00:00 2001 From: Sarah Date: Tue, 22 Sep 2020 16:59:57 +0100 Subject: [PATCH 1/8] Placeholder for halo circuloar velocity function --- skypy/halo/properties.py | 23 +++++++++++++++++++++++ 1 file changed, 23 insertions(+) create mode 100644 skypy/halo/properties.py diff --git a/skypy/halo/properties.py b/skypy/halo/properties.py new file mode 100644 index 000000000..77fc5cdc6 --- /dev/null +++ b/skypy/halo/properties.py @@ -0,0 +1,23 @@ +"""Halo properties module. + +This module provides methods to add simple properties to halos + +Models +====== + +.. autosummary:: + :nosignatures: + :toctree: ../api/ + + halo_circular_velocity + +""" + + + +def halo_circular_velocity(): + """stuff here + """ + + #halos.sort(mass, reverse=True) + \ No newline at end of file From c1a2f53ca166c003310e86edd935a5cb69b01d6b Mon Sep 17 00:00:00 2001 From: Sarah Date: Tue, 22 Sep 2020 17:31:42 +0100 Subject: [PATCH 2/8] drafted the circular velocity code pending testing --- skypy/halo/properties.py | 58 ++++++++++++++++++++++++++++++++++++---- 1 file changed, 53 insertions(+), 5 deletions(-) diff --git a/skypy/halo/properties.py b/skypy/halo/properties.py index 77fc5cdc6..954408fbc 100644 --- a/skypy/halo/properties.py +++ b/skypy/halo/properties.py @@ -13,11 +13,59 @@ """ +import numpy as np +from scipy import integrate +from functools import partial +from astropy import units +__all__ = [ + 'halo_circular_velocity', + ] -def halo_circular_velocity(): - """stuff here - """ +def halo_circular_velocity(M, Delta_v, redshift, cosmology): + """Halo circular velocity. + This function computes the halo circular velocity, setting it + equal to the virial velocity using equation (3) from [1]_. - #halos.sort(mass, reverse=True) - \ No newline at end of file + Parameters + ---------- + M : (nm,) array_like + Array for the halo mass, in units of solar mass. + We assume here that the halo mass is equal to the virial mass. + Delta_v : (nm,) array_like + The mean overdensity of the halo. + redshift : (nm,) array_like + The redshift of each halo. + cosmology : astropy.cosmology.Cosmology + Cosmology object providing methods for the evolution history of + omega_matter and omega_lambda with redshift. + + Returns + -------- + circular_velocity: (nm,) array_like + Halo circular velocity for a given mass array, cosmology and redshift, in + units of km s-1. + + Examples + --------- + >>> import numpy as np + >>> from skypy.halo import properties + + This example will compute the halo circular velocity, for a Planck15 cosmology at redshift 0. + + >>> from astropy.cosmology import Planck15 + >>> cosmo = Planck15 + >>> m = 10**np.arange(9.0, 12.0, 2) + >>> properties.halo_circular_velocity(m,!!!!!) + >>> mass.sheth_tormen_mass_function(m, k, Pk, D0, cosmo) + array([!!!, !!!]) + + References + ---------- + .. [1] Maller and Bullock 2004 MNRAS 355 694 DOI:10.1111/j.1365-2966.2004.08349.x + """ + + circular_velocity = 96.6 * (Delta_v * cosmology.omega_matter * cosmology.h **2) * + pow((1 + redshift)/0.3,0.5) * pow(M/1e11,1/3) + + return circular_velocity From 309223a4e2e3f4a7b1e86f4c5ebe0e28ccf86dd5 Mon Sep 17 00:00:00 2001 From: Sarah Date: Tue, 22 Sep 2020 17:54:34 +0100 Subject: [PATCH 3/8] managed to run this version and get some numbers out. noted some major issues to come back to at the end of the comments though --- skypy/halo/properties.py | 25 +++++++++++++++---------- 1 file changed, 15 insertions(+), 10 deletions(-) diff --git a/skypy/halo/properties.py b/skypy/halo/properties.py index 954408fbc..9ece1e57b 100644 --- a/skypy/halo/properties.py +++ b/skypy/halo/properties.py @@ -14,9 +14,6 @@ """ import numpy as np -from scipy import integrate -from functools import partial -from astropy import units __all__ = [ 'halo_circular_velocity', @@ -54,18 +51,26 @@ def halo_circular_velocity(M, Delta_v, redshift, cosmology): This example will compute the halo circular velocity, for a Planck15 cosmology at redshift 0. >>> from astropy.cosmology import Planck15 - >>> cosmo = Planck15 - >>> m = 10**np.arange(9.0, 12.0, 2) - >>> properties.halo_circular_velocity(m,!!!!!) - >>> mass.sheth_tormen_mass_function(m, k, Pk, D0, cosmo) - array([!!!, !!!]) + >>> cosmology = Planck15 + >>> M = 10**np.arange(9.0, 12.0, 2) + >>> Delta_v = np.arange(1.0, 1.1, 0.1) + >>> redshift = np.arange(0.3, 1, 0.5) + >>> properties.halo_circular_velocity(M, Delta_v, redshift, cosmology) + References ---------- .. [1] Maller and Bullock 2004 MNRAS 355 694 DOI:10.1111/j.1365-2966.2004.08349.x + + Notes of things to ask / think about: + * where to get Delta_v from - do we have a module for this already, or compute it in here too? + * fix h hack / check h in Maller & Bullock is h100 + * should we call it the circular_velocity or the virial_velocity? + * ditto M - should we call it M_v or pretend virial mass = halo mass throughout? + * should we output the virial radius here as well or is this already done elsewhere? """ - circular_velocity = 96.6 * (Delta_v * cosmology.omega_matter * cosmology.h **2) * - pow((1 + redshift)/0.3,0.5) * pow(M/1e11,1/3) + h = cosmology.H0 / 100 + circular_velocity = 96.6 * (Delta_v * cosmology.Om0 * h **2) * pow((1 + redshift)/0.3,0.5) * pow(M/1e11,1/3) return circular_velocity From 9a46b248e0eaac0db20d857e3e96ac6ea6cbf8d0 Mon Sep 17 00:00:00 2001 From: Sarah Bridle Date: Tue, 29 Sep 2020 14:27:18 +0100 Subject: [PATCH 4/8] Update skypy/halo/properties.py Co-authored-by: Lucia-Fonseca <55983939+Lucia-Fonseca@users.noreply.github.com> --- skypy/halo/properties.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/skypy/halo/properties.py b/skypy/halo/properties.py index 9ece1e57b..28cb15a01 100644 --- a/skypy/halo/properties.py +++ b/skypy/halo/properties.py @@ -22,7 +22,7 @@ def halo_circular_velocity(M, Delta_v, redshift, cosmology): """Halo circular velocity. This function computes the halo circular velocity, setting it - equal to the virial velocity using equation (3) from [1]_. + equal to the virial velocity using equation (3) and footnote 2 from [1]_. Parameters ---------- From 23dfe37966e02a0f80da5369a3d1571ef5c671ae Mon Sep 17 00:00:00 2001 From: Sarah Bridle Date: Tue, 29 Sep 2020 14:27:41 +0100 Subject: [PATCH 5/8] Update skypy/halo/properties.py Co-authored-by: Lucia-Fonseca <55983939+Lucia-Fonseca@users.noreply.github.com> --- skypy/halo/properties.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/skypy/halo/properties.py b/skypy/halo/properties.py index 28cb15a01..14a212199 100644 --- a/skypy/halo/properties.py +++ b/skypy/halo/properties.py @@ -27,8 +27,7 @@ def halo_circular_velocity(M, Delta_v, redshift, cosmology): Parameters ---------- M : (nm,) array_like - Array for the halo mass, in units of solar mass. - We assume here that the halo mass is equal to the virial mass. + Array for the virial mass, in units of solar mass. Delta_v : (nm,) array_like The mean overdensity of the halo. redshift : (nm,) array_like From d2aa0083ae78371e5ec97c7f9a36a42e0bc51f6d Mon Sep 17 00:00:00 2001 From: Sarah Date: Tue, 29 Sep 2020 14:54:39 +0100 Subject: [PATCH 6/8] Responded to comments to tidy up names etc --- skypy/halo/properties.py | 18 ++++++------------ 1 file changed, 6 insertions(+), 12 deletions(-) diff --git a/skypy/halo/properties.py b/skypy/halo/properties.py index 14a212199..3f50cea73 100644 --- a/skypy/halo/properties.py +++ b/skypy/halo/properties.py @@ -19,14 +19,14 @@ 'halo_circular_velocity', ] -def halo_circular_velocity(M, Delta_v, redshift, cosmology): +def halo_circular_velocity(halo_virial_mass, Delta_v, redshift, cosmology): """Halo circular velocity. This function computes the halo circular velocity, setting it equal to the virial velocity using equation (3) and footnote 2 from [1]_. Parameters ---------- - M : (nm,) array_like + halo_virial_mass : (nm,) array_like Array for the virial mass, in units of solar mass. Delta_v : (nm,) array_like The mean overdensity of the halo. @@ -60,16 +60,10 @@ def halo_circular_velocity(M, Delta_v, redshift, cosmology): References ---------- .. [1] Maller and Bullock 2004 MNRAS 355 694 DOI:10.1111/j.1365-2966.2004.08349.x - - Notes of things to ask / think about: - * where to get Delta_v from - do we have a module for this already, or compute it in here too? - * fix h hack / check h in Maller & Bullock is h100 - * should we call it the circular_velocity or the virial_velocity? - * ditto M - should we call it M_v or pretend virial mass = halo mass throughout? - * should we output the virial radius here as well or is this already done elsewhere? """ - h = cosmology.H0 / 100 - circular_velocity = 96.6 * (Delta_v * cosmology.Om0 * h **2) * pow((1 + redshift)/0.3,0.5) * pow(M/1e11,1/3) + virial_velocity = 96.6 * Unit('km s-1') * \ + (Delta_v * cosmology.Om0 * cosmology.h**2) * pow((1 + redshift)/0.3,0.5) * \ + pow(halo_virial_mass/(mass.to_value('Msun') / 1e11),1/3) - return circular_velocity + return virial_velocity From 685ca0ec5aa1fd38be89a733368e137782680776 Mon Sep 17 00:00:00 2001 From: Sarah Date: Tue, 29 Sep 2020 14:56:51 +0100 Subject: [PATCH 7/8] fixed over indented block pointed out by Lucia --- skypy/halo/properties.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/skypy/halo/properties.py b/skypy/halo/properties.py index 3f50cea73..50d960f81 100644 --- a/skypy/halo/properties.py +++ b/skypy/halo/properties.py @@ -28,11 +28,11 @@ def halo_circular_velocity(halo_virial_mass, Delta_v, redshift, cosmology): ---------- halo_virial_mass : (nm,) array_like Array for the virial mass, in units of solar mass. - Delta_v : (nm,) array_like - The mean overdensity of the halo. - redshift : (nm,) array_like - The redshift of each halo. - cosmology : astropy.cosmology.Cosmology + Delta_v : (nm,) array_like + The mean overdensity of the halo. + redshift : (nm,) array_like + The redshift of each halo. + cosmology : astropy.cosmology.Cosmology Cosmology object providing methods for the evolution history of omega_matter and omega_lambda with redshift. From 8e4149a4fa7cf76e49a3d2b020b2b41fbe314076 Mon Sep 17 00:00:00 2001 From: Sarah Date: Fri, 2 Oct 2020 23:46:10 +0100 Subject: [PATCH 8/8] Fixed errors in equation and fixed units - but distracted by reloading news headlines so please check! --- skypy/halo/properties.py | 14 ++++++++------ 1 file changed, 8 insertions(+), 6 deletions(-) diff --git a/skypy/halo/properties.py b/skypy/halo/properties.py index 50d960f81..e4fa6bd98 100644 --- a/skypy/halo/properties.py +++ b/skypy/halo/properties.py @@ -14,6 +14,7 @@ """ import numpy as np +from astropy.units import Unit __all__ = [ 'halo_circular_velocity', @@ -51,19 +52,20 @@ def halo_circular_velocity(halo_virial_mass, Delta_v, redshift, cosmology): >>> from astropy.cosmology import Planck15 >>> cosmology = Planck15 - >>> M = 10**np.arange(9.0, 12.0, 2) - >>> Delta_v = np.arange(1.0, 1.1, 0.1) + >>> halo_virial_mass = 10**np.arange(9.0, 12.0, 2) * Unit.Msun + >>> Delta_v = np.arange(201.0, 201.1, 0.1) >>> redshift = np.arange(0.3, 1, 0.5) - >>> properties.halo_circular_velocity(M, Delta_v, redshift, cosmology) + >>> properties.halo_circular_velocity(halo_virial_mass, Delta_v, redshift, cosmology) References ---------- - .. [1] Maller and Bullock 2004 MNRAS 355 694 DOI:10.1111/j.1365-2966.2004.08349.x + .. [1] Barnes and Haehnelt 2010 equation 3 https://arxiv.org/pdf/1403.1873.pdf """ virial_velocity = 96.6 * Unit('km s-1') * \ - (Delta_v * cosmology.Om0 * cosmology.h**2) * pow((1 + redshift)/0.3,0.5) * \ - pow(halo_virial_mass/(mass.to_value('Msun') / 1e11),1/3) + np.power(Delta_v * cosmology.Om0 * cosmology.h**2 / 24.4, 1.0/6.0) * \ + np.sqrt((1 + redshift) / 3.3) * \ + np.power(halo_virial_mass / (1.0e11 * Unit('Msun')), 1.0/3.0) return virial_velocity