A Cosmology Calculator for the World Wide Web

There are many definitions of distance in cosmology, and it is often frustrating to calculate the kind of distance needed for a given situation. In particular, press releases (W. Clavin 2006, unpublished)² very often use the "light‐travel time distance," D_ltt = c(t₀ - t_em), even though this distance has many undesirable characteristics. I wrote the cosmology calculator in order to be able to calculate D_ltt along with the observable angular and luminosity distances D_A and D_L, respectively. It also computes the proper radial distance, which is the only distance that is compatible with the common definition of distance for moving objects: the spatial separation at a common time. The common time used is right now, or t₀. This is called D_now below. The term D_now, also known as the comoving radial distance, is not measurable from a single position, but can in principle be measured by a family of comoving observers. It is the one distance that satisfies the Hubble law exactly— dD/dt = H(t)D—without any corrections. Other distances are defined using quantities taken from more than one time and generally only satisfy the Hubble law to first order.

The following results have been given many times before. Hogg (1999) gives equations that are equivalent to the ones in this paper, but for completeness in documenting the cosmology calculator, the equations are spelled out below.

The metric is given by

where S(x) is sinh (x) if Ω_tot<1, and sin (x) for Ω_tot>1. The radius of curvature is given by R₀ = (c/H₀)/(|1 - Ω_tot|)^1/2, since I use the normalization that a = 1 at the current time, so the redshift is given by 1 + z = 1/a. The past light cone is given by cdt = a(t)R₀dψ, so the comoving radial distance is

and of course the light‐travel time distance is given by

The exact dynamics of the universe can be captured using the energy equation from the Newtonian approximation to cosmology, which gives

General relativity modifies the acceleration equation by including the pressure as a source of gravity, but this just cancels the variation of the enclosed mass M(<a) caused by the pressure, so the energy equation is exact. Thus, we can write a ˙ as H₀ X, with

where q₀ is the deceleration parameter. Let us define a quantity Z that is the comoving distance divided by the Hubble radius c/H₀. Then

To first order, this agrees with the redshift, so Z = z + .... This definition gives D_now = (cZ/H₀) and D_ltt = (c/H₀) ∫_{1/(1 + z)}¹da/ X. Thus, the angular size distance—defined as D_A = R/ θ, where R is the transverse proper length of an object that subtends an angle θ on the sky—is given by

We can define a function J(x), given by

Then

Note that the luminosity distance D_L is defined to make the inverse square law work perfectly for bolometric fluxes, so that F_bol = L/(4πD²_L) for an object of luminosity L.

The enclosed volume calculation requires the integral of either sin²ψdψ or sinh²ψdψ. The closed case becomes

Note that ψ = (|1 - Ω_tot|)^1/2Z(z). The open case gives

The ratio of V(<z) to the Euclidean formula (4π/3)(R₀ψ)³ is given by

in both the open and closed cases.

The cosmology calculator is implemented as a Web page (CosmoCalc)³ that has a large number of JavaScript definitions in the header, followed by immediately executed JavaScript that writes a frame set to the current document. The frame set calls for the input form (CCform)⁴ and the output page (CCout).⁵ If JavaScript is not enabled, or if there is an error in the JavaScript, then the body of CosmoCalc is displayed. This body is just an error message saying that JavaScript must be enabled.

I have received several requests for the code used in the cosmology calculator; but since the code is in JavaScript, it is in the HTML files in ASCII form. It is easy to save the page source to get the code, and it is easy to modify the code using any text editor. If your modifications introduce an error, then you will see the error message saying that JavaScript must be enabled. This only means that you must find the error in your modified CosmoCalc. I have had to do this dozens of times, so do not be discouraged.

Even if you do not intend to modify the code, downloading the three HTML files will let you run the calculator locally when not connected to the Internet.

The numerical evaluation of the integrals for Z and t is done using the midpoint rule with 1000 panels. While this is not a very efficient use of CPU time, it was very easy to code. And with JavaScript, the CPU time is consumed on the user's computer instead of on the server. The functions being integrated go smoothly to zero as a goes to zero, as long as Ω_r>0.

Another hidden aspect of the cosmology calculator is that it automatically sets the radiation density Ω_rh² to the value appropriate for T₀ = 2.72528 K and three massless neutrino species, Ω_rh² = 4.165 × 10^-5. Here h = H₀/(100 km s^-1 Mpc^-1), and this factor includes a small (<1%) boost in the neutrino density due to a slight transfer to e⁺,e^- annihilation energy into neutrinos (Hannestad & Madsen 1995). Therefore, if you want to verify some simple cases, such as the empty universe, you should use a large value of the Hubble constant, which reduces the relative importance of this radiation term. For example, the "Open" button with H₀ = 97.78 and Ω_M = 0 gives an age of the universe of 9.934 Gyr, which is 0.7% from the expected 10 Gyr. But using H₀ = 977.8 gives an age of 999 Myr, which is only 0.1% from the expected 1 Gyr.

For very early times, the age of the universe is adjusted to allow for the existence of e⁺,e^-, μ⁺,μ^-, etc., using the terms g_* and g_*S for the standard model of particle physics given by Kolb & Turner (1990).

A slightly modified version of the cosmology calculator (DlttCalc) is also posted online.⁶ The input form asks for the light‐travel time in Gyr, instead of redshift. The redshift is found by evaluating the integral for light‐travel time in steps of −0.001 in a, starting from a = 1 until the input value is exceeded, and then interpolating to get a and thus z. Once z is known, the calculations proceed as in CosmoCalc.

Since the cosmology calculator was first written, there have been two developments that have changed the kinds of models that people want to run. I have created a new version of the cosmology calculator (ACC) that has more options and have made it available on the Web.⁷ One development is the discovery of neutrino oscillations, indicating that the assumption of massless neutrinos is not correct. The neutrino temperature is (4/11)^1/3T₀ = 1.95 K, and a typical momentum for a thermal neutrino is ≈ 3kT/c. This corresponds to the rest energy of a neutrino with mc² ≈ 0.0005 eV. Since the neutrinos thermalized while still relativistic, their distribution is that of a relativistic Fermi‐Dirac particle, so the neutrino density is given by

The number of spin states for a neutrino is one, but to allow for antineutrinos, one should use g_s = 2.

The integral over neutrino momentum has to be evaluated for each step of the integration over a, so it needs to be done efficiently, even when the work is done on the user's computer. In the low‐temperature limit when kT_ν≪m_νc², the integral over x evaluates to m_ν(3/4)ζ(3)Γ(3), while in the high‐temperature limit kT_ν≫m_νc², the integral evaluates to (kT_ν/c²)(7/8)ζ(4)Γ(4). Both limits can be evaluated correctly by approximating the integral using

with the single integration knot at x₁ = [(7/8)ζ(4)Γ(4)]/[(3/4)ζ(3)Γ(3)] = 3.151. This approximation has a maximum error of <3%. However, a better approximation is

with

and β = 1.842, which has a maximum error <0.3%. The mass at which neutrinos are semirelativistic at the current epoch is

The final result is that the effective neutrino density can be written as

The normalization of 93.14 eV is from Mangano et al. (2005) and is 1.05% higher than the nominal, due to residual coupling of annihilation energy to the neutrinos. In the relativistic limit, the density is 1.53% higher. Increasing T_ν by 0.48% over the nominal (4/11)^1/3T₀ when computing m_rel allows for the relativistic limit, and this boost is included in equation (17).

This Ω_ν(z) is just a function of z that gives Ω_ν at z = 0 and should not be confused with the actual Ω_ν at z‐values other than zero. In ACC, the neutrino density at z = 0 is subtracted from the input Ω_M, giving separate Ω_CM = Ω_M - Ω_ν(z = 0) for the cold dark matter (CDM) and baryons, and Ω_ν for the neutrinos. Of course, the neutrinos are not included in the radiation, so Ω_rh² = 2.778 × 10^-5(T₀/2.72528)⁴.

For a hierarchical neutrino mass pattern with masses ≈0.001, 0.009, and 0.049 eV, the change in distances introduced by neutrino masses is negligible. At redshifts up to 5 in the Wilkinson Microwave Anisotropy Probe concordance flat ΛCDM cosmology with Ω_M = 0.27 and H₀ = 71, the changes are less than 0.01%. This mass pattern is the default when ACC is loaded. But even the more massive nearly degenerate neutrino mass patterns, such as ≈0.33, 0.33, and 0.33 eV, have a minimal effect on the distances and times.

The second change in the paradigm is the introduction of the equation of state w(z) = P/ρc² as a parameter in the model. I have implemented w(z) as w = w₀ + 2w^'(1 - a), following Linder (2003), who added the factor of 2 normalization to the simple form used by Chevallier & Polarski (2001). This functional form behaves well at high redshift, and it allows an analytic solution for the dark energy density as a function of z. This solution is

The defaults are w = -1 and w^' = 0, and in that case, ACC reduces to the CosmoCalc case.

Unlike the neutrino masses, changes in w have substantial effects on distances and ages. Changing w to −0.7 instead of the default −1 changes the age of the universe by 6% and the luminosity distance at z = 1 by 7% when Ω_M = 0.27, Ω_DE = 0.73, and H₀ = 71 are left unchanged. But for w = -0.7, the model that is consistent with both supernovae (Riess et al. 2004, "gold+silver") and the cosmic microwave background (Page et al. 2003, peak positions) changes to an open model with Ω_M = 0.19, Ω_DE = 0.74, and H₀ = 84, and this gives an age change of 17% and a luminosity distance change of 21% at z = 1. However, if Ω_M and Ω_DE are allowed to vary as free parameters, the observable supernova signal in D_L is reduced to only 0.3% (Wright 2005).

Finally, ACC allows one to input T₀, allowing for easier tests of simple cases. With these changes, the X‐function is

The rest of the calculation of distances is unchanged.

These cosmology calculators are suitable for interactive use, providing fairly quick answers for single cases. Users who wish to use the code for large‐scale calculations should translate it into a compiled language and change the quadrature formula. In particular, if fitting to data sets with redshifts and distance, the data should be sorted by redshift, and the distance integrals should be evaluated for all objects in a single pass through the sorted data. Plotting figures is an obvious case in which redshifts are computed in order: Figure 1 shows the distance measures discussed in this paper as a function of z for three different models. The formulae presented here were translated into Postscript for this figure, resulting in a 4 kbyte file.

Zoom In Zoom Out Reset image size