Energy Density

The function energy() calculates the components of the energy density. The energy density is calculated for all fields regardless of the value of noutput_flds. The output is in two files, one for the components of energy density and one for monitoring energy conservation.

The components are output to a file called energy_ext. The number of columns in this file varies with the number of fields and with the number of terms in the potential. Specifically, the first column contains the time, the next nflds columns contain the kinetic energy ${1 \over 2} \left({\partial f_{pr} \over \partial t_{pr}}\right)^2$, the next nflds columns contain the gradient energy ${1 \over 2 a^2} \vert \nabla f_{pr}\vert^2$, and the remaining columns contain the potential terms. The kinetic and gradient energies are calculated by energy(), whereas the potential energy is calculated by a model-specific function called potential_energy() in model.h. The output is in program units, meaning

$$\rho_{pr} \equiv {A^2 \over B^2} a^{-2s+2r} \rho$$ (5.41)

where

$$\rho = \left<{1 \over 2} \dot{f}^2 + {1 \over 2 a^2} \vert \nabla f\vert^2 + V\right>$$ (5.42)

is the physical energy density and $A$, $B$, $r$, and $s$ are the rescaling parameters defined in equation (6.2). Angle brackets denote averages over the grid. Plugging in the rescalings gives

$$\rho_{pr} = \left< {1 \over 2} f_{pr}'^2 - r {a' \over a} f_... ... a^{-2(s+1)} \vert \nabla_{pr} f_{pr}\vert^2 + V_{pr} \right>.$$ (5.43)

The second file generated by energy() is called conservation_ext. The first column contains the time and the second column contains a quantity used to monitor energy conservation. For the case of no expansion (Minkowski space) this quantity is simply the ratio of total energy density to energy density at the beginning of the run. Of course this quantity should remain close to $1$ throughout the run. In an expanding universe the situation is a little more complicated because energy density is not conserved, but rather decreases in a way determined by the expansion rate and the equation of state of the fields. This redshifting of energy is described by the continuity equation

$$\dot{\rho} + 3 {\dot{a} \over a} (\rho + 3 p) = 0.$$ (5.44)

Equation (5.44) can be derived from the Friedmann equations, so if these are being solved exactly then the continuity equation will be obeyed as well. So a simple way to check energy conservation in an expanding universe is to use one of the Friedmann equations, equations (6.19) and (6.20) to evolve the scale factor and then check if the other is being satisfied as well. In practice the program uses a combination of these two equations so either one can be used as a check of energy conservation. We use equation (6.20), i.e.

$$\left({\dot{a} \over a}\right)^2 = {8 \pi \over 3} \rho.$$ (5.45)

So for an expanding universe the second column of the conservation file contains the ratio of $\left({\dot{a} \over a}\right)^2$ to ${8 \pi \over 3} \rho$. To express this in terms of program variables we use equations (5.2) and (5.41) to give

$$H^2 / \left({8 \pi \over 3} \rho\right) = {3 \over 8 \pi} A^2 a^{2r-2} a'^2 {1 \over \rho_{pr}}.$$ (5.46)

As in the case of no expansion this quantity should remain close to $1$ throughout the run. We have generally found that the deviation of this ratio from one is comparable to the lack of energy conservation for runs done with the same models without expansion. We have also checked that using the second Friedmann equation, equation (6.19), gives essentially the same results as this method.

In the case of fixed power-law expansion the conservation file is not created.