Monday, April 11, 2016

Enhancing the Resolution of N-body Merger Trees

We've just submitted a paper to the arXiv which describes a technique for enhancing the resolution of merger trees extracted from N-body simulations. The goal of this technique is to ensure that semi-analytic models applied to merger trees from cosmological N-body simulations give converged results.

This isn't as simple as it first appears - even if an N-body simulation has resolution sufficient to resolve the halos which you expect to host galaxies of interest with, say, 100 particles, that doesn't mean it has sufficient resolution to resolve the progenitors of that halo. Since structure formation in CDM is inherently hierarchical any halo gets built through the merging of smaller systems.

Importantly for semi-analytic models galaxies which form in those smaller, progenitor halos can influence the physical properties of the galaxy forming in the final halo. For example, galaxies in progenitor halos can produce metals which contaminate the halo of the main galaxy altering cooling rates, or can simply lock up baryons into stars making them unavailable for star formation in the main galaxy. If those progenitor halos are not resolved then this potentially important activity lower down the hierarchy is missed, and the main galaxy properties can be computed incorrectly.

To remedy this problem, we propose augmenting the resolution of N-body merger trees by grafting in higher resolution branches - generated using a Monte Carlo procedure (the PCH algorithm from Parkinson, Cole & Helly 2008) - which match the existing halo masses of the tree but provide a statistically valid sample of lower mass progenitor halos.

The process is illustrated by this schematic:

Panel a: A simplified diagram of a merger tree extracted from an N-body simulation. Circles represent halos (with radius indicating mass). Time increases in the direction shown by the arrow, and halos are located at quantized redshifts, labelled z0...z2 , and shown by horizontal, dot-dot-dashed lines. Dashed lines connect halos to their direct progenitors. One halo at z0 (highlighted by the yellow color and red outline) is selected for augmentation. This halo’s progenitors at z1 are also highlighted in yellow. Panel b: A trial tree (shown in blue) is generated using the PCH algorithm and compared to the halos of interest in the original tree. In this case, the match between trial tree and original tree halo masses is sufficiently close
to be deemed acceptable. Panel c: The accepted trial tree is grafted into the original tree. Note that the augmented node and its progenitors from the original tree are retained (so their masses are unchanged), but now with the structure of the trial tree grafted between them. Panel d: Where the trial tree has halos at z1 which did not match any halo in the existing tree (and are below the resolution of the original tree by construction), a new tree is grown from each such halo and attached. This process is repeated for each of the green halos in the original tree.

How well does it work? Very well! Check the paper for tests of the convergence of the technique and its ability to recover the statistics of N-body merger tree conditional mass functions. We also explore how well it leads to convergence in galaxy properties. Here's an example where we augment tree from the Millennium Simulation to match the resolution of the Millennium-II Simulation and then compare per-halo galaxy stellar mass functions:
Solid lines are for central galaxies, dashed are for satellite galaxies. Galaxies formed in unaugmented Millennium trees (red lines) differ significantly in their properties from those formed in Millennium-II trees (blue lines). These trees differ only in their resolutions, all other properties are identical. If we augment the Millennium trees to match the resolution of Millennium-II we get the green lines, which now agree much more closely with the results from the Millennium-II trees.

The augmenting code is available within Galacticus as a merger tree operator. A merger tree operator is applied to each tree before galaxy formation calculations begin and can make arbitrary modifications to the tree. In the case of augmenting we actually want to apply several operators. So, to augment you'd add something like this to your parameter file:

<mergerTreeOperatorMethod value="sequence">
  <mergerTreeOperatorMethod value="regridTimes">
    <expansionFactorEnd value="1.0" />
    <expansionFactorStart value="0.047" />
    <regridCount value="60" />
    <snapshotRedshifts value="19.915688 18.243723 16.724525 15.343073 14.085914 12.940780
                  11.896569 10.943864 10.073462  9.277915  8.549912  7.883204
                  7.272188  6.711586  6.196833  5.723864  5.288833  4.888449
                  4.519556  4.179469  3.865683  3.575905  3.308098  3.060419
                  2.831182  2.618862  2.422044  2.239486  2.070027  1.912633
                  1.766336  1.630271  1.503636  1.385718  1.275846  1.173417
                  1.077875  0.988708  0.905463  0.827699  0.755036  0.687109
                  0.623590  0.564177  0.508591  0.456577  0.407899  0.362340
                  0.319703  0.279802  0.242469  0.207549  0.174898  0.144383
                  0.115883  0.089288  0.064493  0.041403  0.019933  0.000000" />
    <snapshotSpacing value="list" />
      <mergerTreeOperatorMethod value="pruneByMass">
    <massThreshold value="7.08e10" />
    <preservePrimaryProgenitor value="true" />
      <mergerTreeOperatorMethod value="augment">
    <attemptsMaximum value="10000" />
    <massCutOff value="7.08e10" />
    <massCutOffAttemptsMaximum value="50" />
    <massCutOffScaleFactor value="0.0" />
    <mergerTreeBuilderMethod value="cole2000">
          <accretionLimit value="0.1" />
          <branchIntervalStep value="true" />
          <mergeProbability value="0.1" />
          <mergerTreeMassResolutionMethod value="fixed">
            <massResolution value="2.5e9" />
          <randomSeedsFixed value="false" />
    <rescaleMaximum value="20" />
    <retryMaximum value="50" />
    <snapshotRedshifts value="19.915688 18.243723 16.724525 15.343073 14.085914 12.940780
                  11.896569 10.943864 10.073462  9.277915  8.549912  7.883204
                  7.272188  6.711586  6.196833  5.723864  5.288833  4.888449
                  4.519556  4.179469  3.865683  3.575905  3.308098  3.060419
                  2.831182  2.618862  2.422044  2.239486  2.070027  1.912633
                  1.766336  1.630271  1.503636  1.385718  1.275846  1.173417
                  1.077875  0.988708  0.905463  0.827699  0.755036  0.687109
                  0.623590  0.564177  0.508591  0.456577  0.407899  0.362340
                  0.319703  0.279802  0.242469  0.207549  0.174898  0.144383
                  0.115883  0.089288  0.064493  0.041403  0.019933  0.000000" />
    <toleranceScale value="0.15" />
      <mergerTreeOperatorMethod value="conditionalMF">
    <massRatioCount value="25" />
    <massRatioMaximum value="2.0" />
    <massRatioMinimum value="2.0e-6" />
    <outputGroupName value="treeStatisticsAfterInsertion" />
    <parentMassCount value="20" />
    <parentMassMaximum value="1.0e15" />
    <parentMassMinimum value="2.36e10" />
    <parentRedshifts value="0.0 0.0" />
    <primaryProgenitorDepth value="4" />
    <progenitorRedshifts value="0.508591 1.077875" />

We use a sequence operator, which simply applies the sequence of operators it contains in order to each tree. That sequence of operators does the following:
  1. The regridTimes operator ensures that the halos of the tree are located precisely at the given set of redshifts. It will snap halos to those redshifts if they differ slightly (e.g. due to lack of precision in the input data), and interpolate halos to snapshots where they are missing. This ensures that the tree is in a suitable state for augmenting.
  2. The pruneByMass operator prunes away branches below the given mass scale - we want to ensure that we don't attempt to augment low mass halos where the N-body-determined properties may be unreliable (due to limited numbers of particles).
  3. The augment operator does the work of augmenting the tree. We specify the mass at which the trees were cut off (by the pruneByMass operator), and give it a tree builder algorithm with which to construct new branches (this also specifies the mass resolution of those new branches).
  4. Finally, we use a conditionalMF operator to accumulate statistics on the tree (including the conditional mass function) and output those to the Galacticus output file.
This technique should allow converged galaxies properties to be predicted even in simulations which do not fully resolve the galaxy formation hierarchy.

Wednesday, May 13, 2015

A New Format For Galacticus Parameter Files

As part of the Galacticus v0.9.4 development, we've introduced a major change in the structure of Galacticus parameter files. This change provides three improvements:
  1. Allows for more compact XML;
  2. Makes associations between groups of parameters more clear (and prevents the need for extremely long parameter names!);
  3. Paves the way for more complete object composition within Galacticus defined directly through the XML parameter file.
As an example, here's the old format for parameter files showing setting of cosmological parameters:

  <!-- Cosmological parameters -->


In the new format, this same set of parameters is described as:

  <!-- Cosmological parameters and options -->
  <cosmologyParametersMethod value="simple">
    <HubbleConstant  value="70.2"   />
    <OmegaMatter     value="0.2725" />
    <OmegaDarkEnergy value="0.7275" />
    <OmegaBaryon     value="0.0455" />
    <temperatureCMB  value="2.72548"/>


The parameter name is now simply the name of each XML element, and the parameter value is stored as an attribute of the element. (You can also store the value in a <value> element within each parameter if necessary.) In this example, values of the "simple" cosmology parameters method are now specified as subparameters within the cosmologyParametersMethod element itself.

If you have old-format parameter files, you can convert them to the new format using the migration script:

scripts/aux/ oldParametersFile.xml newParametersFile.xml

Internally, this new format will be used to allow detailed construction of objects (e.g. the ability to construct multiple cosmological models each with different sets of cosmology parameters), and, for example, to construct nested sets of filters which can be applied to merger tree construction (e.g. build a tree → measure the tree information content → prune the tree → measure the information content again).

Thursday, December 18, 2014

Galacticus v0.9.3 released

We're happy to announce the release of Galacticus v0.9.3. As usual, you can get Galacticus direct from our repo at BitBucket, just use:

hg clone -u v0.9.3

to grab v0.9.3, or you can download a source tarball. v0.9.3 becomes the stable version (supported and will receive bug fixes), v0.9.2 becomes deprecated (no longer supported and will not receive bug fixes) and development moves to v0.9.4.

One of the major internal changes to Galacticus in v0.9.3 has been to begin moving all of the algorithms to a fully object-oriented framework. This mostly makes for cleaner (and less!) code, but has some notable effects visible to users. In particular, there are now several cases where an algorithm functions by modifying another algorithm. A good example of this is the new warm dark matter halo concentrations algorithm. This works by taking a CDM halo concentration algorithm and modifying it. This functionality is possible because of the new object-oriented design.

Some of the other notable changes and improvements to Galacticus are described below.

New components

Several new component classes have been added to Galacticus, most of which allow the tracking of various useful statistics of galaxies and halos:
  • Galaxy dynamics - this component records a time series of properties related to bar instability in galaxy disks (instability timescale and adiabatic ratio);
  • Hot halo mass outflows - tracks the mass of gas in a galaxy's hot halo which arrived there directly via outflows;
  • Mass flows - tracks the cumulative mass of gas which has ever cooled onto a galaxy
  • Host halo history - tracks the maximum mass ever reached by a satellite galaxy's host halo;
  • Age statistics - tracks the stellar mass-weighted mean ages of disk and spheroid components of galaxies;

 New outputs

  • Output filters:
    • Filter on a combination of stellar mass and spheroid-to-total mass ratio;
    • The lightcone filter can now prune trees which lie entirely outside the lightcone so that they do not have to be processed at all;
      • This filter also now handles non-small angle geometries correctly.
  • Physical properties:
    • The λR parameter of Cappellari et al. (2007);
    • Stellar half-mass radii; 
    • Satellite status (i.e. if satellite is an orphan);

New physical models

New tools and functionality

  • A script is now provided to help migrate Galacticus parameter files between versions - it intelligently updates your parameter file to handle changes in parameter names and definitions between versions;
  • The script used to launch batches of Galacticus models has been re-written and modularized, making it easier to launch models on different systems (PBS, SLURM, etc.);
  • Reading of merger trees from file is now split into data import (reading the data from file) and data processing (manipulating the data into a form that Galacticus can use), allowing for import from multiple different merger tree file formats;
  • Updates of external codes:
  • Luminosity filters:
    • Arbitrary tophat filters can now be specified using topHat_Lmin_Lmax_R - this gets expanded into a top-hat filter  between wavelengths Lmin and Lmax with resolution R;
    • A filter redshift can now be set to "all" to cause the filter to be used for all available output redshifts;
  • Optimizations:
    • HDF5 output can now be more effectively buffered which results in much lower I/O overheads when outputting large numbers of redshifts;
    • Stellar population luminosities integrated under filters can be stored to file for rapid re-use;
  • Help:
    • Lists of suitable components to support requested functionality are automatically reported whenever the selected component is insufficient;
    • When a requested method is unrecognized, a list of available methods is printed.

Friday, December 12, 2014

Tracking Satellite Halo Orbits

The final major new addition to v0.9.3 of Galacticus has just been added. This is the subhalo evolution code written by Anthony Pullen which was used in our recent work on the subhalo population in warm dark matter halos.

The new model is similar to previous works by Taylor & Babul (2001), Benson et al. (2002), and Zentner et al. (2005). It tracks the orbit of each subhalo directly, accounting for the gravitational potential of the host halo, dynamical friction, tidal mass loss, and tidal shocking, but ignoring subhalo-subhalo interactions to minimize the computational cost.

The result is that Galacticus can predict the radial and mass distributions of subhalos directly, without the need to extract these from N-body simulations. An example of the distribution of halos is shown below (image credit: Anthony Pullen):

Each circle represents a subhalo within a Milky Way-mass dark matter halo at z=0, with the radius indicating the current tidal radius of the subhalo.

Wednesday, December 3, 2014

Modeling Reionization

As the result of work by Caltech SURF student, Daniel McAndrew, Galacticus now has a fully self-consistent calculation of reionization (and the ionization and thermal state of the intergalactic medium in general). This model uses the spectrum of photons emitted by model galaxies and AGN to compute the evolution of ionization states of hydrogen and helium in the intergalactic medium (IGM), along with the IGM temperature. The thermal evolution of the IGM is then used to compute the filtering mass, and this affects later galaxy evolution via the gas accretion rate into halos (Naoz & Barkana, 2007).

The result is an internally self-consistent model which can predict when reionization occurs. There are some simplifications (the biggest among them probably being the escape fraction of ionizing photons from galaxies, and the clumping factor in the IGM) but we hope that future work will address some of those.

To activate this calculation, include the following in your parameters file:

  <!-- IGM evolver -->
  <!-- Background radiation -->

  <!-- Halo accretion options -->

The first block of parameters switches Galacticus to using an internal calculation for the state of the IGM, instructs it to solve for IGM properties as a function of time, and specifies that IGM properties should be updated 10 times per decade of cosmic time. Specifically, at each of these time intervals, solving of galaxy evolution is halted and the IGM evolved up to this time using the currently computed photoionizing background spectrum.

The second block of parameters activates an internal calculation of cosmic background radiation, in which the background is computed from the emissivities of model galaxies and AGN. The number of points at which to tabulate the background per decade of wavelength and cosmic time are specified.

Finally, the third block of parameters tells Galacticus to use the Naoz & Barkana (2007) prescription for computing gas accretion into halos from the IGM. This prescription uses the filtering mass to determine accretion rates, and will take the filtering mass from the internal IGM evolution calculation.

Once completed, data on the IGM and background radiation are written to the Galacticus output file in the igmProperties and backgroundRadiation groups respectively.

Sunday, November 2, 2014

Computing Emission Line Properties

The latest update to Galacticus v0.9.3 allows calculation of emission line luminosities for galaxies. Much of the work to make this happen was carried out by Pomona College student Gabe Currier over the summer of 2013. The calculation is based on the methodology of Panuzzo et al. (2003). Briefly, HII region models are constructed using Cloudy for a variety of gas densities and metallicities, and HI, HeI, and OII ionizing luminosities. Emission line luminosities are then computed by interpolating in these tables based on the instantaneous properties of model galaxies.

To compute emission line luminosities it is therefore necessary to run Galacticus including the following rest-frame luminosity filters at each redshift of interest: Lyc, HeliumContinuum, OxygenContinuum. This causes the ionizing luminosity for each three species to be computed and output.

Emission line luminosities can then be found using Galacticus's data extraction modules, by simply importing the Galacticus::EmissionLines module. Emission line luminosites (in units of Solar luminosities) can then be accessed as named properties using names such as  


which will return the Hα luminosity at z=0. Equivalent properties for disk and spheroid are provided (simply replace "total" with "disk" or "spheroid").

Currently available lines are:
  • balmerAlpha6563
  • balmerBeta4861
  • oxygenII3726
  • oxygenII3729
  • oxygenIII4959
  • oxygenIII5007
  • nitrogenII6584
  • sulfurII6731
  • sulfurII6716

Additionally, if a line is requested as


then the line luminosity is computed under the provided filter, and the luminosity is returned in units of maggies for easy conversion to AB magnitudes.

Sunday, September 14, 2014

Migrating parameter files between Galacticus versions

The names and allowed values of parameters often change between versions of Galacticus - due to changes in naming conventions, attempts to unify parameter names, and changes in the structure of the code.

To permit easy and error-free migration between versions a script is provided to translate parameter files from earlier to later versions - preserving the layout of your parameter file, including comments. To migrate a parameter file simply use:

scripts/aux/ parameters.xml newParameters.xml

By default, this script will translate from the previous to the current version of Galacticus. If your parameter file contains a version element then this will be used to determine which version of Galacticus the parameter file was constructed for. The migration script will then migrate the parameter file through all intermediate versions to bring it into compliance with the current version. You can also specify input and output versions directly:

scripts/aux/ parameters.xml newParameters.xml --inputVersion 0.9.0 --outputVersion 0.9.3

will convert parameters.xml from version 0.9.0 syntax to version 0.9.3 syntax.