Cosmology and Astrophysics from Relaxed Galaxy Clusters II

Cosmological Constraints

Introduction

Massive, dynamically relaxed clusters of galaxies represent a minority of the cluster population, and are a precious resource for probing both the astrophysics of clusters and cosmology. They provide precise and minimally biased mass estimates to studies of scaling relations and the growth of cosmic structure, and enable complementary contraints on cosmic expansion. Additionally, they represent the most natural targets for studying thermodynamic features of the intracluster medium with minimal systematics from deprojection or non-equilibrium processes.

In this paper, we employ a sample of relaxed clusters to obtain constraints on cosmology and the expansion of the Universe. This work is enabled by the prediction, from hydrodynamical simulations of cluster formation, that the ratio of gas mass to total mass (fgas) in rich clusters is approximately independent of redshift and mass, and has small cluster-to-cluster variability. X-ray mass estimates of these clusters (assuming hydrostatic equilibrium) are expected to be minimally biased by non-thermal pressure, and can now be directly calibrated using robust gravitational lensing measurements. By measuring gas and total masses for a sample of relaxed clusters, and employing simulation results for the expected depletion of gas in clusters compared to the Universe as a whole, we can constrain the cosmic ratio of baryon to total matter density. Combining this result with external constraints on the present-day expansion rate of the Universe (H0) and baryon density, tight and robust constraints can be placed on the density of matter in the Universe (Ωm). At the same time, the wide range in redshift covered by the fgas data set provides additional information about the history of cosmic expansion. These data provide corroboration of the finding (from type Ia supernovae) that the cosmic expansion is accelerating, and place competitive constraints on models of the dark energy which is thought to cause this acceleration.

Data and likelihood code

There are three ways to get the code. If you use any of these in your own work, you should cite this paper (and others in the series, if applicable).

  1. As source plus data files. The download clusters.zip contains
    • Data files.
    • A stand-alone code library for evaluating the likelihood for a given cosmological model. Requires the GNU Scientific Library. The code is written in C++, but C and Fortran interfaces are also included.
    • Additional files and code to use the data with COSMOMC. See instructions.
    • (Partially complete interfaces for Python and the cosmological code MontePython also exist, and can be obtained from the author by request.)
    For questions/problems related to these files, contact the authors.
  2. As of June 2015 (version 1.2), cosmosis includes this code under the module named "cluster_fgas".
  3. As of July 2014, the COSMOMC git has a clusters_fgas branch which includes the fgas code and data. You can email Antony Lewis, the developer of COSMOMC, for access. This is probably the easiest way to get started using our results with COSMOMC. The instructions below apply to an earlier version of COSMOMC, but will remain here just in case they're useful.

Important note: For the time being, we are witholding our gravitational lensing data, pending the acceptance of a paper specifically detailing the calibration of X-ray and lensing masses (Applegate et al.). The code is designed to work even when these additional data files are absent. The mass calibration result can be largely captured by a 0.9 ± 0.9 Gaussian prior on the ratio of lensing to X-ray masses, although this fails to account for a mild degeneracy with cosmological parameters. In COSMOMC, this is accomplished by adding the line prior[cl_cal] = 0.9 0.09 to the ini file, and fixing the remaining cl_ parameters (as done in the current download).

Reference for these results

Cosmology and Astrophysics from Relaxed Galaxy Clusters II: Cosmological Constraints.
A. B. Mantz, S. W. Allen, R. G. Morris, D. A. Rapetti, D. E. Applegate, P. L. Kelly, A. von der Linden, and R. W. Schmidt. MNRAS 440:2077, 2014. ADS, arXiv, BibTeX.

See also paper I: Sample Selection 

Instructions for using the source distribution with COSMOMC

Apart from compiling the package above, the installation consists of pasting a few short blocks of code into the COSMOMC source. These instructions work as of the December 2013 COSMOMC version and are not guaranteed in perpetuity. If you run into trouble or find errors in these instructions, feel free to contact me.

  1. Go to your COSMOMC directory.
  2. Unzip the download and run make in the fgas/src directory to compile the code. You may need to customize the Makefile.
  3. Move or link fgas/cosmomc/clusters.ini to the COSMOMC batch1 directory (optional, for convenience).
  4. Move or link fgas/cosmomc/clusters.f90 to the COSMOMC source directory (required).
  5. In the COSMOMC source/DataLikelihoods.f90, add these lines in the obvious places.
    use clusters
    call ClustersLikelihood_Add(DataLikelihoods, Ini)
      
  6. In the COSMOMC source/Makefile, make these additions as indicated:
    ### at the top is fine
    # clusters module
    CLUSTERS ?= ../fgas
    ifneq ($(CLUSTERS),)
    CLUSTERSO = clusters.o
    CLUSTERSF = clusters.f90
    endif
    
    ### add $(CLUSTERSO) to the end of the DATAMODULES definition
    
    ### can go anywhere after the OBJFILES and LINKFLAGS definitions
    ### note that you may need to specify the GSL path using -L in the 2nd statement
    ifneq ($(CLUSTERS),)
    OBJFILES += $(CLUSTERS)/src/fwrapper.o
    LINKFLAGS += -lgsl -lgslcblas
    endif
    
    ### can go anywhere after the LINKFLAGS definition
    ifneq ($(CLUSTERS),)
    LINKFLAGS += -lstdc++ -L$(CLUSTERS)/src -lclusters
    endif
      
  7. Add the following to modules.f90 in the camb subdirectory of COSMOMC. Just after the usual AngularDiameterDistance definition is a logical place.
        function AngularDiameterDistance2(z1, z2) ! z1 < z2
          real(dl) AngularDiameterDistance2
          real(dl), intent(in) :: z1, z2
          AngularDiameterDistance2 = CP%r/(1+z2)*rofchi(ComovingRadialDistance(z2)/CP%r - ComovingRadialDistance(z1)/CP%r)
        end function AngularDiameterDistance2
      
  8. You should now be able to make COSMOMC. If you've previously compiled, do a make clean first to force CAMB to recompile completely as well.

Below is an example of a minimal top-level .ini file using the clusters module. Note that, due to a legacy "use_clusters" keyword that is no longer used in COSMOMC as distributed but nevertheless has a default value, the use_clusters=T instruction must be in the top-level file.

DEFAULT(batch1/common_batch1.ini)

use_clusters = T
DEFAULT(batch1/clusters.ini)

file_root = chains/test

The following is optional, and implements the cosmological parametrization described in the paper, which is convenient for the fgas data.

  1. In source/params_CMB.f90, add the following code blocks in the indicated places:
    ! in the module definition
        Type, extends(CosmologyParameterization) :: clBackgroundParameterization
           logical :: flat = .false.
           real(mcp) :: H0_prior_mean = 0._mcp, H0_prior_std = 0._mcp
           real(mcp) :: ombh2_prior_mean = 0._mcp, ombh2_prior_std = 0._mcp
        contains
        procedure :: ParamArrayToTheoryParams => clBK_ParamArrayToTheoryParams
        procedure :: NonBaseParameterPriors => clBK_NonBaseParameterPriors
        procedure :: CalcDerivedParams => clBK_CalcDerivedParams
        procedure :: Initialize => clBK_Init
        end type clBackgroundParameterization
    
    ! at the top of SetTheoryParametrization
        Type(clBackgroundParameterization), pointer :: clBackgroundParam
    ! and in the following branching statement
        else if (paramtxt=='cluster_background') then
            allocate(clBackgroundParam)
            Parameterization => clBackgroundParam
            call clBackgroundParam%Initialize(Ini,Names)
    
    ! the rest can go at the end of the module
        !!! background parameterization with a little extra stuff (namely omb) for clusters
        subroutine clBK_Init(this, Ini, Names)
        Class(clBackgroundParameterization) :: this
        Type(TIniFile) :: Ini
        Type(TParamNames) :: Names
        character(LEN=Ini_max_string_len) prior
    
        call SetTheoryParameterNumbers(9,0)
        this%late_time_only = .true.
    
        ! special flag to require a flat universe, since omegak is not a free parameter
        this%flat = Ini_Read_Logical_File(Ini, 'enforce_flatness', this%flat)
    
        ! use the special H0 prior keys just because they already exist
        prior = Ini_Read_String_File(Ini, 'H0_prior', NotFoundFail=.false.)
        if (prior/='') then
            read(prior,*) this%H0_prior_mean, this%H0_prior_std
        end if
        ! also honor the omegabh2 prior, if there is one
        prior = Ini_Read_String_File(Ini, 'prior[omegabh2]', NotFoundFail=.false.)
        if (prior/='') then
            read(prior,*) this%ombh2_prior_mean, this%ombh2_prior_std
        end if
    
        call this%Init(Ini,Names, 'params_clbackground.paramnames')
    
        end subroutine clBK_Init
    
        function clBK_NonBaseParameterPriors(this,CMB)
        class(clBackgroundParameterization) :: this
        class(TTheoryParams) :: CMB
        real(mcp):: clBK_NonBaseParameterPriors
    
        select type (CMB)
        class is (CMBParams)
            clBK_NonBaseParameterPriors = logZero
            if (CMB%YHe < 0.0) return ! this catches ombh2 out of bounds with bbn_consistency (below)
            clBK_NonBaseParameterPriors = 0
            if (this%H0_prior_mean/=0._mcp) then
                clBK_NonBaseParameterPriors = clBK_NonBaseParameterPriors + ((CMB%H0 - this%H0_prior_mean)/this%H0_prior_std)**2/2
            end if
            if (this%ombh2_prior_mean/=0._mcp) then
                clBK_NonBaseParameterPriors = clBK_NonBaseParameterPriors + ((CMB%ombh2 - this%ombh2_prior_mean)/this%ombh2_prior_std)**2/2
            end if
        end select
        end function clBK_NonBaseParameterPriors
    
        subroutine clBK_ParamArrayToTheoryParams(this, Params, CMB)
        use bbn
        class(clBackgroundParameterization) :: this
        real(mcp) Params(:)
        class(TTheoryParams), target :: CMB
        real(mcp) fbaryon, omegam, h2
    
        select type (CMB)
        class is (CMBParams)
            fbaryon = Params(1)
            omegam = Params(2)
            CMB%H0 = 100.0 * Params(3)
            if (this%flat) then
               CMB%omv = 1.0 - omegam
            else
               CMB%omv = Params(4)
            end if
            CMB%omnuh2=Params(5)/neutrino_mass_fac*(3.046_mcp/3)**0.75_mcp
            CMB%w =    Params(6)
            CMB%wa =    Params(7)
            CMB%nnu =    Params(8)
    
            CMB%h=CMB%H0/100
            h2 = CMB%h**2
            CMB%omnu = CMB%omnuh2/h2
            CMB%omb = fbaryon * omegam
            CMB%ombh2 = CMB%omb*h2
            CMB%omc = (1.0 - fbaryon) * omegam - CMB%omnu
            CMB%omch2 = CMB%omc*h2
            CMB%zre=0
            CMB%tau=0
            CMB%omdmh2 = CMB%omch2+ CMB%omnuh2
            CMB%omdm = CMB%omdmh2/h2
            CMB%omk = 1 - CMB%omv - CMB%omb - CMB%omdm
            CMB%nufrac=CMB%omnuh2/CMB%omdmh2
            CMB%reserved=0
            CMB%fdm=0
            CMB%iso_cdm_correlated=0
            CMB%Alens=1
            if (bbn_consistency) then
               if (CMB%ombh2 > 0.005 .and. CMB%ombh2 < 0.040) then ! magic #s
                  CMB%YHe = yp_bbn(CMB%ombh2,CMB%nnu  - 3.046)
               else
                  CMB%Yhe = -1.0
               end if
            else
               CMB%YHe = Params(9)
            end if
        end select
        end subroutine clBK_ParamArrayToTheoryParams
    
    
        function clBK_CalcDerivedParams(this, P, Theory, derived) result (num_derived)
        class(clBackgroundParameterization) :: this
        Type(mc_real_pointer) :: derived
        class(TTheoryPredictions) :: Theory
        real(mcp) :: P(:)
        Type(CMBParams) CMB
        integer num_derived
    
        num_derived = 3
    
        allocate(Derived%P(num_derived))
    
        call this%ParamArrayToTheoryParams(P,CMB)
    
        derived%P(1) = CMB%omk
        derived%P(2) = CMB%omb
        derived%P(3) = CMB%Yhe
    
        end function clBK_CalcDerivedParams
     
  2. Put the following in params_clbackground.paramnames in your top-level COSMOMC directory:
    fbaryon       f_b # cosmic baryon fraction
    omegam        \Omega_m
    h             h
    omegal        \Omega_\Lambda
    mnu           \Sigma m_\nu              
    w             w         #equation of state parameter for scalar field dark energy today
    wa            w_a       #w_a variation
    nnu           N_{eff}   #effective number of neutrinos (only clearly defined for massless)
    yhe           Y_{He} # helium mass fraction
    omegak*       \Omega_K
    omegab*       \Omega_b
    yheused*      Y_{He} # used value (if set from bbn_consistency)
     
  3. Here is a generic .ini file with defaults for this parametrization:
    parameterization = cluster_background
    get_sigma8 = F
    enforce_flatness = T
    
    param[fbaryon] = 0.17 0.0 1.0 0.03 0.03
    param[omegam] = 0.27 0.0 1.0 0.04 0.04
    param[h] = 0.7 0.01 10.0 0.03 0.03
    param[omegal] = 0.73 0.0 2.0 0 0
    num_massive_neutrinos = 3
    param[mnu] = 0 0 0 0 0
    param[w] = -1 -5 0 0 0
    param[wa] = 0 0 0 0 0
    param[nnu] = 3.046 3.046 3.046 0 0
    param[yhe] = 0.24 0.24 0.24 0 0