Subhalo time evolution

The history-based approach of HBT-HERONS means that the evolution of subhaloes is tracked simultaneously with their identification in each simulation output. No additional algorithm is therefore required to link subhaloes forward in time, so the information required to build merger trees is already contained within the catalogues.

Since HBT-HERONS analyses a simulation from early to late times, the information used to follow the evolution of subhaloes is primarily represented through the notion of descendant subhaloes. There are two types of fate a given subhalo can have at any given output, and the dataset used to find the descendant of a subhalo varies depending on which of the two it experienced:

It remains self-bound. The subhalo with the same TrackId, which is unique and time-persistent for each subhalo, is its descendant.
It merges with another subhalo. Depending on the way in which this happened, i.e. through disruption or because their cores overlapped in phase-space, the DescendantTrackId or SinkTrackId values identify which subhalo it merged with.

We provide examples on how to use these datasets to follow the complete evolution of a subhalo of interest, and to identify subhaloes that merged with it.

Evolution of a single subhalo

Tip

If you are analysing a simulation with many subhaloes and time outputs, we strongly recommend using the toolbox/catalogue_cleanup/SortCatalogues.py script. This script merges all files generated by different MPI ranks into a single one, and sorts subhaloes (and their properties) according to their TrackId. Doing the sorting helps speed up subsequent file reads since the index location of a subhalo is known beforehand, as it is simply its TrackId.

The TrackId of a subhalo is a unique identifier that persists in time throughout the simulation. Thus, following the evolution of a subhalo from the first output when it is found in the simulation (SnapshotOfBirth) until the last output when it is resolved as self-bound (SnapshotOfDeath) only requires knowing its TrackId. In fact, the subhalo can still be tracked as an orphan subhalo after its "death", but only a subset of properties are computed in that case.

Code example

Here we show how to follow the mass evolution of the most massive subhalo identified in the last output of a HBT-HERONS analysis.

After running toolbox/catalogue_cleanup/SortCatalogues.pyWithout running toolbox/catalogue_cleanup/SortCatalogues.py

import h5py
from glob import glob

# Get ordered list to the sorted subhalo catalogues. Create a dictionary to access its paths by output number.
catalogue_paths = sorted(glob("<SORTED_CATALOGUE_BASE_PATH>/OrderedSubSnap_*.hdf5"))
catalogue_paths = dict([(int(path[-8:-8+3]),path) for path in catalogue_paths])
max_output_number = list(catalogue_paths)[-1]

# Get the TrackId of the most massive subhalo at the last available output,
# when it was first identified and when it disrupted/merged.
with h5py.File(catalogue_paths[max_output_number]) as catalogue:
    TrackId_to_follow = catalogue["Subhalos"]["Mbound"][()].argmax()
    output_start = catalogue["Subhalos"]["SnapshotOfBirth"][TrackId_to_follow]
    output_end   = catalogue["Subhalos"]["SnapshotOfDeath"][TrackId_to_follow]

# If output_end is equal to -1, that means it is still resolved at the time when the output was saved.
output_end = output_end if output_end != -1 else max_output_number

# Create an array to hold values we are interested in tracking (number of bound particles)
number_resolved_outputs = list(catalogue_paths).index(output_end) \
                        - list(catalogue_paths).index(output_start) + 1
Mbound_evolution   = - np.ones(number_resolved_outputs)
Snapshot_evolution = - np.ones(number_resolved_outputs)

# Iterate over catalogues to obtain Nbound value of the entry with the TrackId we want to follow.
entry_nr = 0
for output_number, path in catalogue_paths.items():
    if (output_number < output_start) | (output_number > output_end):
        continue

    with h5py.File(path) as catalogue:
        Mbound_evolution[entry_nr] = catalogue["Subhalos/Mbound"][TrackId_to_follow]
        Snapshot_evolution[entry_nr] = output_number
    entry_nr += 1

import sys
sys.path.append("<HBT-HERONS_PATH>/toolbox")
from HBTReader import HBTReader

# The reader will parse the base folder and parameter file to identify number
# of outputs and if any are missing.
catalogue = HBTReader("<HBT-HERONS_CATALOGUE_BASE_PATH>")

# Get the TrackId of the most massive subhalo. The reader loads the latest
# available snapshot by default and all subhalo properties.
subhaloes = catalogue.LoadSubhalos()
TrackId_to_follow = subhaloes["TrackId"][subhaloes["Mbound"].argmax()]

# Get its bound mass evolution, which returns by default all properties and
# the associated scale factors and snashot output numbers.
subhalo_evolution  = catalogue.GetTrackEvolution(TrackId_to_follow)
Mbound_evolution   = subhalo_evolution["Mbound"]
Snapshot_evolution = subhalo_evolution["Snapshot"]

We now have an array that contains the number of bound particles for that subhalo across time. We can plot it using the code below:

import matplotlib.pyplot as plt

fig, ax1 = plt.subplots(1)
ax1.plot(Snapshot_evolution, Mbound_evolution, 'k-')
ax1.set_xlabel('Output Number')
ax1.set_ylabel('Total mass of bound particles [HBT-HERONS units]')
ax1.set_yscale('log')
plt.show()

Identifying subhalo mergers

The catalogues also provide sufficient infomation to identify which subhaloes merged together, helping to establish links between different evolutionary branches through so-called secondary progenitors. Obtaining the secondary progenitors of a given subhalo is more involved that following its main branch. One complication is that subhaloes disappear from the HBT-HERONS catalogues in two different ways:

Subhalo disruption as a consequence of subhaloes no longer being self-bound.
Subhalo sinking when the subhalo is still self-bound but its core coalesces in phase-space with that of another subhalo.

HBT-HERONS provides different information depending on which of the two processes led to the removal of a subhalo from the simulation. Commonly used merger trees do not provide this two-category classification, as following the entire process of subhalo sinking is difficult and is inadecuately followed in traditional subhalo finders.

Subhalo disruption

Disruption occurs when the subhalo is no longer considered to be self-bound. There are two conditions within HBT-HERONS used to determine whether a subhalo is self-bound:

The total number of bound particles has to be equal or greater than MinNumPartOfSub.
The total number of bound tracer particles has to be equal or greater than MinNumTracerPartOfSub.

If either of these conditions is not satisfied, the subhalo is considered as disrupted (Nbound = 0) and is subsequently tracked as an orphan subhalo.

To identify the descendant of a disrupted subhalo, HBT-HERONS uses the NumTracersForDescendants most bound tracer particles when the subhalo was last self-bound. The descendant subhalo (DescendantTrackId) is the one that is self-bound and contains the majority of the aforementioned particles.

Code example

Here we show how to identify all subhaloes that disrupted and hence merged with the same subhalo of the first example.

After running toolbox/catalogue_cleanup/SortCatalogues.pyWithout running toolbox/catalogue_cleanup/SortCatalogues.py

import h5py
from glob import glob

# Get ordered list to the sorted subhalo catalogues. Create a dictionary to access its paths by output number.
catalogue_paths = sorted(glob("<SORTED_CATALOGUE_BASE_PATH>/OrderedSubSnap_*.hdf5"))
catalogue_paths = dict([(int(path[-8:-8+3]),path) for path in catalogue_paths])
max_output_number = list(catalogue_paths)[-1]

# Get the TrackId of the most massive subhalo in the last available output, and
# all progenitor TrackIds.
with h5py.File(catalogue_paths[max_output_number]) as catalogue:
    TrackId_to_follow = catalogue["Subhalos"]["Mbound"][()].argmax()
    disrupted_progenitors  = np.where(catalogue["Subhalos/DescendantTrackId"][()] == TrackId_to_follow)[0]

import sys
sys.path.append("<HBT-HERONS_PATH>/toolbox")
from HBTReader import HBTReader

# The reader will parse the base folder and parameter file to identify number
# of outputs and if any are missing.
catalogue = HBTReader("<HBT-HERONS_CATALOGUE_BASE_PATH>")

# Get the TrackId of the most massive subhalo. The reader loads the latest
# available snapshot by default and all subhalo properties.
subhaloes = catalogue.LoadSubhalos()
TrackId_to_follow = subhaloes["TrackId"][subhaloes["Mbound"].argmax()]
disrupted_progenitors = subhaloes["TrackId"][subhaloes["DescendantTrackId"] == TrackId_to_follow]

The disrupted_progenitors array contains a TrackId for each subhalo that disrupted and whose core ended bound to the example subhalo (TrackId_to_follow). If you want to find every subhalo that is connected to the example subhalo through disruption, e.g. to build its full merger tree, you will need to repeat the above code for each of the disrupted_progenitors until a complete list is built. The evolution of each subhalo prior to disrupting can be followed in the same way as for an individual subhalo.

Subhalo sinking

The sinking of a subhalo refers to process of its core becoming indistinguishable in phase-space from the core of another subhalo, whilst both remain self-bound. The coalescence occurs as dynamical friction makes the core of the least massive of the pair to "sink" towards the core of more massive one. As the efficiency of dynamical friction increases when the masses of both subhaloes becomes comparable, sinking only happens for non-negligible mass ratios.

Formally, HBT-HERONS checks the phase-space offset between the cores of subhaloes that have a hierarchical connection. If it identifies them to be within a threshold distance, HBT-HERONS flags the least massive of the two as being sunk and removes it from the simulation. The TrackId of the most massive of the two is assigned as the SinkTrackId of the subhalo that sunk.

Code example

Here we show how to identify all subhaloes that directly sunk with the same subhalo of the first example.

After running toolbox/catalogue_cleanup/SortCatalogues.pyWithout running toolbox/catalogue_cleanup/SortCatalogues.py

import h5py
from glob import glob

# Get ordered list to the sorted subhalo catalogues. Create a dictionary to access its paths by output number.
catalogue_paths = sorted(glob("<SORTED_CATALOGUE_BASE_PATH>/OrderedSubSnap_*.hdf5"))
catalogue_paths = dict([(int(path[-8:-8+3]),path) for path in catalogue_paths])
max_output_number = list(catalogue_paths)[-1]

# Get the TrackId of the most massive subhalo in the last available output, and
# all progenitor TrackIds.
with h5py.File(catalogue_paths[max_output_number]) as catalogue:
    TrackId_to_follow = catalogue["Subhalos"]["Mbound"][()].argmax()
    sink_progenitors  = np.where(catalogue["Subhalos/SinkTrackId"][()] == TrackId_to_follow)[0]

import sys
sys.path.append("<HBT-HERONS_PATH>/toolbox")
from HBTReader import HBTReader

# The reader will parse the base folder and parameter file to identify number
# of outputs and if any are missing.
catalogue = HBTReader("<HBT-HERONS_CATALOGUE_BASE_PATH>")

# Get the TrackId of the most massive subhalo. The reader loads the latest
# available snapshot by default and all subhalo properties.
subhaloes = catalogue.LoadSubhalos()
TrackId_to_follow = subhaloes["TrackId"][subhaloes["Mbound"].argmax()]
sink_progenitors = subhaloes["TrackId"][subhaloes["SinkTrackId"] == TrackId_to_follow]

The sink_progenitors array contains a TrackId for each subhalo that sunk with the example subhalo (TrackId_to_follow). If you want to find every subhalo that is connected to the example subhalo through sinking, e.g. to build its full merger tree, you will need to repeat the above code for each of the sink_progenitors until a complete list is built. The evolution of each subhalo prior to sinking can be followed in the same way as for an individual subhalo.