Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Improve cluster velocity dispersion #101

Merged
merged 3 commits into from
Apr 18, 2024
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
38 changes: 21 additions & 17 deletions cogsworth/hydro/pop.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,8 +15,8 @@


class HydroPopulation(Population):
def __init__(self, star_particles, particle_size=1 * u.pc, virial_parameter=1.0, subset=None,
sampling_params={}, snapshot_type=None, **kwargs):
def __init__(self, star_particles, cluster_radius=3 * u.pc, cluster_mass=1e4 * u.Msun,
virial_parameter=1.0, subset=None, sampling_params={}, snapshot_type=None, **kwargs):
"""A population of stars sampled from a hydrodynamical zoom-in snapshot

Each star particle in the snapshot is converted to a (binary) stellar population, sampled with COSMIC,
Expand All @@ -28,8 +28,11 @@ def __init__(self, star_particles, particle_size=1 * u.pc, virial_parameter=1.0,
A table of star particles from a snapshot that contains the mass, metallicity, formation time,
position and velocity of each particle, as returned by
:func:`~cogsworth.hydro.rewind.rewind_to_formation`
particle_size : :class:`~astropy.units.Quantity`, optional
Size of gaussian for cluster for each star particle, by default 1*u.pc
cluster_radius : :class:`~astropy.units.Quantity`, optional
Size of gaussian for cluster for each star particle, by default 3 * u.pc
cluster_mass : :class:`~astropy.units.Quantity`, optional
Mass of cluster for each star particle for calculating velocity dispersion,
by default 1e4 * u.Msun
virial_parameter : `float`, optional
Virial parameter for each cluster, for setting velocity dispersions, by default 1.0
subset : `int` or `list` of `int`, optional
Expand All @@ -47,7 +50,8 @@ def __init__(self, star_particles, particle_size=1 * u.pc, virial_parameter=1.0,
Parameters to pass to :class:`~cogsworth.pop.Population`
"""
self.star_particles = star_particles
self.particle_size = particle_size
self.cluster_radius = cluster_radius
self.cluster_mass = cluster_mass
self.virial_parameter = virial_parameter
self._subset_inds = self.star_particles.index.values
if subset is not None and isinstance(subset, int):
Expand Down Expand Up @@ -117,7 +121,7 @@ def __getitem__(self, ind):
timestep_size=self.timestep_size, BSE_settings=self.BSE_settings,
sampling_params=self.sampling_params,
store_entire_orbits=self.store_entire_orbits,
virial_parameter=self.virial_parameter, particle_size=self.particle_size)
virial_parameter=self.virial_parameter, cluster_radius=self.cluster_radius)

new_pop.n_binaries = len(bin_nums)
new_pop.n_binaries_match = len(bin_nums)
Expand Down Expand Up @@ -228,31 +232,31 @@ def sample_initial_galaxy(self):
v_x = particles["v_x"].values * u.km / u.s
v_y = particles["v_y"].values * u.km / u.s
v_z = particles["v_z"].values * u.km / u.s

pos = np.random.normal([x.to(u.kpc).value, y.to(u.kpc).value, z.to(u.kpc).value],
self.particle_size.to(u.kpc).value / np.sqrt(3),
self.cluster_radius.to(u.kpc).value / np.sqrt(3),
size=(3, self.n_binaries_match)) * u.kpc

self._initial_galaxy = StarFormationHistory(self.n_binaries_match, immediately_sample=False)
self._initial_galaxy._x = pos[0]
self._initial_galaxy._y = pos[1]
self._initial_galaxy._z = pos[2]
self._initial_galaxy._tau = self._initial_binaries["tphysf"].values * u.Myr
self._initial_galaxy._Z = self._initial_binaries["metallicity"].values
self._initial_galaxy._which_comp = np.repeat("FIRE", len(self.initial_galaxy._tau))

v_R = (x * v_x + y * v_y) / (x**2 + y**2)**0.5
v_T = (x * v_y - y * v_x) / (x**2 + y**2)**0.5

vel_units = u.km / u.s
dispersion = dispersion_from_virial_parameter(self.virial_parameter,
self.particle_size,
particles["mass"].values * u.Msun)
self.cluster_radius, self.cluster_mass)
v_R, v_T, v_z = np.random.normal([v_R.to(vel_units).value,
v_T.to(vel_units).value,
v_z.to(vel_units).value],
dispersion.to(vel_units).value / np.sqrt(3),
size=(3, self.n_binaries_match)) * vel_units

self._initial_galaxy = StarFormationHistory(self.n_binaries_match, immediately_sample=False)
self._initial_galaxy._tau = self._initial_binaries["tphysf"].values * u.Myr
self._initial_galaxy._Z = self._initial_binaries["metallicity"].values
self._initial_galaxy._which_comp = np.repeat("FIRE", len(self.initial_galaxy._tau))

self._initial_galaxy._x = pos[0]
self._initial_galaxy._y = pos[1]
self._initial_galaxy._z = pos[2]
self._initial_galaxy._v_R = v_R
self._initial_galaxy._v_T = v_T
self._initial_galaxy._v_z = v_z
Loading