yt.analysis_modules.halo_finding.rockstar.rockstar module

Operations to get Rockstar loaded up

class yt.analysis_modules.halo_finding.rockstar.rockstar.InlineRunner[source]

Bases: yt.utilities.parallel_tools.parallel_analysis_interface.ParallelAnalysisInterface

comm = None
get_dependencies(fields)
partition_index_2d(axis)
partition_index_3d(ds, padding=0.0, rank_ratio=1)
partition_index_3d_bisection_list()

Returns an array that is used to drive _partition_index_3d_bisection, below.

partition_region_3d(left_edge, right_edge, padding=0.0, rank_ratio=1)

Given a region, it subdivides it into smaller regions for parallel analysis.

run(handler, pool)[source]
setup_pool()[source]
class yt.analysis_modules.halo_finding.rockstar.rockstar.RockstarHaloFinder(ts, num_readers=1, num_writers=None, outbase='rockstar_halos', particle_type='all', force_res=None, total_particles=None, dm_only=False, particle_mass=None, min_halo_size=25)[source]

Bases: yt.utilities.parallel_tools.parallel_analysis_interface.ParallelAnalysisInterface

Spawns the Rockstar Halo finder, distributes dark matter particles and finds halos.

The halo finder requires dark matter particles of a fixed size. Rockstar has three main processes: reader, writer, and the server which coordinates reader/writer processes.

Parameters:
  • ts (DatasetSeries, Dataset) – This is the data source containing the DM particles. Because halo IDs may change from one snapshot to the next, the only way to keep a consistent halo ID across time is to feed Rockstar a set of snapshots, ie, via DatasetSeries.
  • num_readers (int) – The number of reader can be increased from the default of 1 in the event that a single snapshot is split among many files. This can help in cases where performance is IO-limited. Default is 1. If run inline, it is equal to the number of MPI threads.
  • num_writers (int) – The number of writers determines the number of processing threads as well as the number of threads writing output data. The default is set to comm.size-num_readers-1. If run inline, the default is equal to the number of MPI threads.
  • outbase (str) – This is where the out*list files that Rockstar makes should be placed. Default is ‘rockstar_halos’.
  • particle_type (str) – This is the “particle type” that can be found in the data. This can be a filtered particle or an inherent type.
  • force_res (float) – This parameter specifies the force resolution that Rockstar uses in units of Mpc/h. If no value is provided, this parameter is automatically set to the width of the smallest grid element in the simulation from the last data snapshot (i.e. the one where time has evolved the longest) in the time series: ds_last.index.get_smallest_dx().in_units("Mpc/h").
  • total_particles (int) – If supplied, this is a pre-calculated total number of particles present in the simulation. For example, this is useful when analyzing a series of snapshots where the number of dark matter particles should not change and this will save some disk access time. If left unspecified, it will be calculated automatically. Default: None.
  • particle_mass (float) – If supplied, use this as the particle mass supplied to rockstar. Otherwise, the smallest particle mass will be identified and calculated internally. This is useful for multi-dm-mass simulations. Note that this will only give sensible results for halos that are not “polluted” by lower resolution particles. Default: None.
Returns:

Return type:

None

Examples

To use the script below you must run it using MPI: mpirun -np 4 python run_rockstar.py

>>> import yt
>>> yt.enable_parallelism()
>>> from yt.analysis_modules.halo_finding.rockstar.api import \
...     RockstarHaloFinder
>>> # create a particle filter to remove star particles
>>> @yt.particle_filter("dark_matter", requires=["creation_time"])
... def _dm_filter(pfilter, data):
...     return data["creation_time"] <= 0.0
>>> def setup_ds(ds):
...     ds.add_particle_filter("dark_matter")
>>> es = yt.simulation("enzo_tiny_cosmology/32Mpc_32.enzo", "Enzo")
>>> es.get_time_series(setup_function=setup_ds, redshift_data=False)
>>> rh = RockstarHaloFinder(es, num_readers=1, num_writers=2,
...                         particle_type="dark_matter")
>>> rh.run()
comm = None
get_dependencies(fields)
partition_index_2d(axis)
partition_index_3d(ds, padding=0.0, rank_ratio=1)
partition_index_3d_bisection_list()

Returns an array that is used to drive _partition_index_3d_bisection, below.

partition_region_3d(left_edge, right_edge, padding=0.0, rank_ratio=1)

Given a region, it subdivides it into smaller regions for parallel analysis.

run(block_ratio=1, callbacks=None, restart=False)[source]
class yt.analysis_modules.halo_finding.rockstar.rockstar.StandardRunner(num_readers, num_writers)[source]

Bases: yt.utilities.parallel_tools.parallel_analysis_interface.ParallelAnalysisInterface

comm = None
get_dependencies(fields)
partition_index_2d(axis)
partition_index_3d(ds, padding=0.0, rank_ratio=1)
partition_index_3d_bisection_list()

Returns an array that is used to drive _partition_index_3d_bisection, below.

partition_region_3d(left_edge, right_edge, padding=0.0, rank_ratio=1)

Given a region, it subdivides it into smaller regions for parallel analysis.

run(handler, wg)[source]
setup_pool()[source]