Gravitational theory of cosmology, galaxies and galaxy clusters

Abstract

A modified gravitational theory explains early universe and late time cosmology, galaxy and galaxy cluster dynamics. The modified gravity (MOG) theory extends general relativity (GR) by three extra degrees of freedom: a scalar field G, enhancing the strength of the Newtonian gravitational constant $G_N$, a gravitational, spin 1 vector graviton field $\phi _\mu$, and the effective mass $\mu$ of the ultralight spin 1 graviton. For $t < t_\mathrm{rec}$, where $t_\mathrm{rec}$ denotes the time of recombination and re-ionization, the density of the vector graviton $\rho _\phi > \rho _b$, where $\rho _b$ is the density of baryons, while for $t > t_\mathrm{rec}$ we have $\rho _b > \rho _\phi$. The matter density is parameterized by $\Omega _M=\Omega _b+\Omega _\phi +\Omega _r$ where $\Omega _r=\Omega _\gamma +\Omega _\nu$. For the cosmological parameter values obtained by the Planck Collaboration, the CMB acoustical oscillation power spectrum, polarization and lensing data can be fitted as in the $\Lambda$CDM model. When the baryon density $\rho _b$ dominates the late time universe, MOG explains galaxy rotation curves, the dynamics of galaxy clusters, galaxy lensing and the galaxy clusters matter power spectrum without dominant dark matter.

1 Introduction

Dark matter was introduced to explain the stable dynamics of galaxies and galaxy clusters. General relativity (GR) with only ordinary baryon matter cannot explain the present accumulation of astrophysical and cosmological data without dark matter. However, dark matter has not been observed in laboratory experiments [1]. Therefore, it is important to consider a modified gravitational theory. The observed acceleration of the universe has also complicated the situation by needing a dark energy, either in the form of the cosmological constant vacuum energy or as a modification of GR.

The fully relativistic and covariant modified gravity (MOG) theory Scalar–Tensor–Vector–Gravity (STVG)  [2,3,4,5,6,7,8,9,10,11,12,13,14,15,16, 23,24,25,26] has been successfully applied to explain the rotation curves of galaxies and the dynamics of galaxy clusters. It is claimed in Ref. [14] that mass density profiles for galaxies near the centre of the clusters Abell 1689 and Abell 1835 are in conflict with the requirement of no dark matter in the fitting of MOG to the cluster data. The observational data in Fig. 1, Ref. [15], go below $r\lesssim 10$ kpc. For $r\lesssim 30$ kpc the observational data for the normalized acceleration in Fig. 1 Ref. [15] are well inside a galaxy interior. One should question the observational resolution at small r near the centre of galaxy clusters. While the observational data for $r\lesssim 10$ kpc may be correct, the dark matter models, such as the three fitted to the acceleration data of Abell 1689, and the MOG model must include in the fitting in Fig. 1 in Ref. [15] a modeling of the internal dynamics of galaxies for $r\lesssim 20$ kpc. The dark matter and MOG density models used to fit the Abell 1689 data in Fig. 1 in Ref. [15] did not account for the fitting of galaxies. Therefore, it is this deficiency and not the need for dark matter or the failure of MOG to fit the data for $r\lesssim 20$ kpc that is in conflict. A similar problem was encountered in the fitting of a large set of galaxy clusters in Ref. [10].

In addition to the enhanced gravitational scalar field strength G, the theory also has a repulsive gravitational field generated by a massive, spin 1 Proca vector graviton (VG) field $\phi _{\mu }$. The mass of the VG is determined by a scalar field $\mu$. The value of $\mu$ that fits the galaxy rotation curves and the cluster dynamics is $\mu \sim 0.01 - 0.04\,\mathrm{kpc}^{-1}$, corresponding to $\mu ^{-1}\sim 25 - 100$ kpc and a mass $m_\phi \sim ~ 10^{-28}\,\mathrm{eV}$. The mass of the VG in the early universe is $m_\phi \sim 10^{-22}$ eV. The massless spin 2 graviton associated with $g_{\mu \nu }$ and the massive VG are difficult and probably impossible to detect [19, 20]. This will explain the so far lack of success of detecting a dark matter particle.

The difference between standard dark matter models and MOG is that dark matter models assume that GR is the correct theory of gravity and a dark matter particle such as WIMPS, axions and fuzzy dark matter [21] are postulated to belong to the standard particle model. MOG extends GR by the addition of new gravitational degrees of freedom. The cosmological observations attributed to undetected dark matter particles are explained as being of purely gravitational origin. The gravitational coupling of the VG to matter is universal with the gravitational charge $Q_g=\sqrt{8\pi \alpha G_N}M=\kappa M$, where $\alpha$ is a dimensionless parameter associated with the scalar coupling $G=G_N(1+\alpha )$, M is the mass of a body, $G_N$ is Newton’s gravitational constant. Moreover, $\kappa = \mathcal{O}(1/M_{PL})$ where $M_\mathrm{PL}= 1/\sqrt{8\pi G_N}=2.435\times 10^{18}$ GeV is the reduced Planck mass. As the universe expands to the present epoch, the VGs become ultra-relativistic and cannot form bound states such as galaxy and galaxy cluster halos. MOG becomes the dominant gravitational phase with the density of baryons $\rho _b$ exceeding the density $\rho _\phi$ of the VGs with the decreased mass $m_\phi \sim 10^{-28}$ eV. Our formulation of gravitation fits the cosmological data described by the $\Lambda$CDM model, as well as the present epoch galaxy and galaxy cluster data without dominant dark matter.

A complete description of the universe by a theory requires not only that the present day observations of galaxies and galaxy clusters require an explanation, but also that the wealth of cosmological observations is described as well by the theory. For MOG there are two possible approaches to an explanation of the cosmological data. One is to pursue a description without dark matter employing only baryonic matter and the baryonic density $\rho _b$. This approach has been explored with some success [25,26,27] but requires further exploration in the future. The second approach is to promote the gravitational VG (spin 1 vector graviton) degree of freedom as a dominant density $\rho _\phi$ in the early universe, so that it dominates the baryon density $\rho _b$ for the time before the reionization phase when stars and galaxies are first formed [2, 22]. After this time, there is a transition to modified gravity when $\rho _b > \rho _\phi$. This transition can be caused by a matter-energy phase transition. In the following, we present a further investigation of the second approach.

2 The Friedmann Equations

We base our cosmology on the homogeneous and isotropic Friedmann–Lemaitre–Robertson–Walker (FLRW) background metric:

$$\begin{aligned} ds^2=dt^2-a^2(t)\biggl [\frac{dr^2}{1-Kr^2}+r^2(d\theta ^2+\sin ^2\theta d\phi ^2)\biggr ], \end{aligned}$$

(1)

where $K=-1,0,1\,(\mathrm{length}^{-2})$ for open, flat and closed universes, respectively. We use the energy-momentum tensor of a perfect fluid in (A-17).

We have

$$\begin{aligned} \rho =\rho _b+\rho _\phi +\rho _G+\rho _\mu +\rho _r, \end{aligned}$$

(2)

where $\rho _b,\rho _\phi , \rho _G$ and $\rho _\mu$ denote the densities of baryon matter, the electrically neutral VG field $\phi _\mu$, the scalar G field, the effective mass $\mu$ of the VG and radiation, respectively. The radiation density $\rho _r=\rho _\gamma +\rho _\nu$, where $\rho _\gamma$ and $\rho _\nu$ denote the densities of photons and neutrinos, respectively.

In the following, we assume a spatially flat universe $K=0$ and $V_G=V_\mu =0$. Furthermore, we assume that the time dependence of the gravitational field strength is negligible, $\dot{G}\sim 0$, and we also assume that $\dot{\mu }\sim 0$. We obtain the approximate Friedmann equations:

$$\begin{aligned} \biggl (\frac{\dot{a}}{a}\biggr )^2=\frac{8\pi G\rho }{3}+\frac{\Lambda }{3}, \end{aligned}$$

(3)

$$\begin{aligned} \frac{\ddot{a}}{a}=-\frac{4\pi G}{3}(\rho +3p)+\frac{\Lambda }{3}. \end{aligned}$$

(4)

The energy conservation equation is

$$\begin{aligned} \dot{\rho }+3\frac{d\ln a}{dt}(\rho +p)=0. \end{aligned}$$

(5)

For the matter dominated universe with $p=0$ we have $\rho _b\propto 1/a^3$.

We introduce comoving coordinates:

$$\begin{aligned} \mathbf{x}=\frac{a_0\mathbf{r}}{a(t)}. \end{aligned}$$

(6)

The Fourier expanded density, pressure and gravitational potential in terms of plane waves with comoving wavenumber $\mathbf{k}$ are given by

$$\begin{aligned} \delta \rho ({\mathbf{r},t})=\frac{1}{(2\pi )^{3/2}}\bar{\rho }\int d^3k\delta _\mathbf{k}(t)\exp (i(a_0/a)\mathbf{k}\cdot \mathbf{r}), \end{aligned}$$

(7)

$$\begin{aligned} \delta p({\mathbf{r},t})=\frac{1}{(2\pi )^{3/2}}\int d^3k\delta p_\mathbf{k}(t)\exp (i(a_0/a)\mathbf{k}\cdot \mathbf{r}), \end{aligned}$$

(8)

$$\begin{aligned} \Phi ({\mathbf{r},t})=\frac{1}{(2\pi )^{3/2}}\int d^3k\Phi _\mathbf{k}(t)\exp (i(a_0/a)\mathbf{k}\cdot \mathbf{r}), \end{aligned}$$

(9)

where $\delta =\delta \rho /{\bar{\rho }}$ is the relative density perturbation contrast. We have for the gravitational potential:

$$\begin{aligned} k^2\Phi _\mathbf{k}=-4\pi G\biggl (\frac{a}{a_0}\biggr )^2\bar{\rho }\delta _\mathbf{k}, \end{aligned}$$

(10)

where the MOG potential for a given density $\rho ({\mathbf {r}})$ is

$$\begin{aligned} \Phi (\mathbf {r}) = - G_N \int d^3{\mathbf {r}}'\frac{\rho (\mathbf {r}')}{|\mathbf {r}-\mathbf {r}'|}\Big [1+\alpha -\alpha e^{-\mu |\mathbf {r}-\mathbf {r}'|}\Big ]. \end{aligned}$$

(11)

The MOG acceleration equation is given by

$$\begin{aligned} a_\mathrm{MOG}({\mathbf {r}})= & {} -G_N\int d^3{\mathbf {r}}'\frac{\rho ({\mathbf {r}}')({\mathbf {r}}-{\mathbf {r}}')}{|{\mathbf {r}}-{\mathbf {r}}'|^3} [1+\alpha \nonumber \\&-\alpha \exp (-\mu |{\mathbf {r}}-{\mathbf {r}}'|(1+\mu |{\mathbf {r}}-{\mathbf {r}}'|))]. \end{aligned}$$

(12)

The density of VGs in the early universe is given by

$$\begin{aligned} \Omega _\phi \sim 0.12\biggl (\frac{\kappa }{10^{18}\,\mathrm{GeV}}\biggr )^2\biggl (\frac{m_\phi }{10^{-22}\,\mathrm{eV}}\biggr )^{1/2}. \end{aligned}$$

(13)

The energy density $\rho _{\phi }$ scales as $1/a^3$ just like matter in the matter dominated universe. The mean velocity of the VGs is given by

$$\begin{aligned} {\bar{v}}\sim \frac{2\pi \hbar }{\lambda m_\phi }\sim 10\,\mathrm{kms^{-1}}\biggl (\frac{10^{-22}\,\mathrm{eV}}{m_\phi }\biggr ), \end{aligned}$$

(14)

where the VG de Broglie wave length $\lambda /2\pi \sim \hbar /m_\phi {\bar{v}}$.

3 MOG structure growth and angular acoustical power spectrum

We assume that when $\rho _{\phi }$ dominates in the early universe, $\alpha < 1$ and $G\sim G_N$. The particle densities are expressed as the ratios $\Omega _x=8\pi G_N\rho _x/3H^2$. In particular, we have for the baryon, VGs and the cosmological constant $\Lambda$:

$$\begin{aligned} \Omega _b=\frac{8\pi G_N\rho _b}{3H^2},\quad \Omega _\phi =\frac{8\pi G_N\rho _\phi }{3H^2},\quad \Omega _{\Lambda }=\frac{\Lambda }{3H^2}. \end{aligned}$$

(15)

At the time of big-bang nucleosynthesis (BBN), we have $G\sim G_N$, guaranteeing that the production of elements at the time of BBN agrees with observations. We assume that at horizon entry until some time after decoupling, $\rho _\phi \gg \rho _b$, $\rho _\phi \gg \rho _G$, $\rho _\phi \gg \rho _\mu$, so that $\rho \sim \rho _b+\rho _\phi +\rho _r$ [2, 22].

The equation for the density perturbation for the non-relativistic single-component fluid is given by the Jeans equation:

$$\begin{aligned} \ddot{\delta }_\mathbf{k}+2H\dot{\delta }_\mathbf{k}+\bigg (\frac{c_s^2a_0^2k^2}{a^2}-4\pi G_N\bar{\rho }\biggr )\delta _\mathbf{k}=0, \end{aligned}$$

(16)

where $H={\dot{a}}/a$ and $\bar{\rho }$ denotes the mean density. Here, $c_s$ is the speed of sound:

$$\begin{aligned} c_s=\sqrt{\frac{dp}{d\rho }}, \end{aligned}$$

(17)

and perturbations with $\delta p/\delta \rho =(d\bar{p}/dt)/(d\bar{\rho }/dt)=c_s^2$ are called adiabatic perturbations and for these perturbations $p=p(\rho )$. The first term in the parenthesis of (16) is due to the adiabatic perturbation pressure contribution $\delta p=c_s^2\delta \rho$.

The nature of the solutions to (16) depends on the sign of the factor in parenthesis. Pressure attempts to resist compression, so when the pressure term dominates, we obtain an oscillatory solution comprising standing density (sound) waves. The second term is due to gravity and when this term dominates the perturbations grow. The Jeans wavenumber when the pressure and gravity terms are equal is given by

$$\begin{aligned} k_J=\frac{a}{a_0}\frac{\sqrt{4\pi G_N\bar{\rho }}}{c_s}, \end{aligned}$$

(18)

corresponding to the wavelength $\lambda _J=2\pi /k_J$ (Jeans length) and for non-relativistic matter $k_J\gg \mathcal{H}=(a/a_0)H$ and $c_s\ll 1$ where $\mathcal{H}$ is the comoving Hubble scale.

The Proca vector field $\phi _\mu$ is a neutral massive VG. Because the VGs are electrically neutral, they do not couple to photons and they can be treated as almost pressureless, so the pressure gradient term in the Jeans equation is absent and $c_s\sim 0$. The Jeans length is approximately zero and we have

$$\begin{aligned} \ddot{\delta }_\mathbf{k}+2H\dot{\delta }_\mathbf{k}-4\pi G_N\bar{\rho }\delta _\mathbf{k}=0. \end{aligned}$$

(19)

Fig. 1

The TT power spectrum obtained from Ref. [28], using the best-fit parameters of the Planck Collaboration [28]. The blue line shows the best-fit MOG spectrum. The lower panel shows the residuals with respect to the $\Lambda$CDM model and the MOG model

For the radiation dominant era $\rho _\phi < \rho _r$, but perturbations in radiation, $\delta \rho _r$, oscillate rapidly and are damped $\delta \rho _r\sim 0$, so they can be ignored compared to the $\delta \rho _\phi$ perturbations. The expansion law is that for a radiation dominated universe:

$$\begin{aligned} a(t)\propto t^{1/2},\quad H=\frac{1}{2t}=\sqrt{\frac{8\pi G_N}{3}\bar{\rho }}. \end{aligned}$$

(20)

Dividing (19) by $H^2$ we get

$$\begin{aligned} \frac{\ddot{\delta }_\mathbf{k}}{H^2}+2\frac{\dot{\delta }_\mathbf{k}}{H}=\frac{3}{2}\frac{\bar{\rho }_\phi }{\bar{\rho }}\delta _\mathbf{k}. \end{aligned}$$

(21)

The right-hand-side can be dropped because $\bar{\rho }_\phi < \bar{\rho }$. The solution to this equation is given by

$$\begin{aligned} \delta _\mathbf{k}=a+b\ln t, \end{aligned}$$

(22)

where a and b are constants. We find that the perturbations due to the VGs grow at most logarithmically during the radiation dominated era of the universe. The increased expansion rate due to the presence of the smooth radiation component slows down the growth of perturbations.

Because the VGs are almost pressureless and the speed of sound $c_s\sim 0$, the VG Jeans length is approximately zero, so all scales are larger than the Jeans scale, and there is no oscillatory behavior for the VG particles. Instead, perturbations grow at all length scales. On the other hand, the $\rho _{b\gamma }$ perturbations oscillate before decoupling. The dark energy determined by the cosmological constant $\Lambda$ is a smooth component of energy and only becomes important at later times.

Due to decoupling the baryon Jeans length is comparable to the Hubble length, so adiabatic baryon perturbations will oscillate before decoupling. The baryon density perturbation begins to grow only after the decoupling time $t=t_{dec}$, because before decoupling the baryon-photon pressure prevents any growth. Without the VG density $\rho _\phi$ the baryon perturbations grow as $\delta _b\propto a\propto t^{2/3}$ after decoupling, and the baryon density perturbations $\delta \rho _b$ catch up to the VG density perturbations $\delta \rho _\phi$. However, the CMB anisotropy due to baryon density perturbations observed at $t=t_{dec}$ are too small $\sim 10^{-4}$ to produce a growth factor of 1100 needed to generate the presently observed large-scale structure. The VG density perturbations solve the problem, for they begin to grow earlier at horizon entry and by $t=t_{dec}$ they are significantly greater than the baryon density perturbations.

The baryon sound wave oscillations due to the baryon-photon pressure prior to the decoupling time produce acoustical peaks in the angular power spectrum, $\mathcal{D}_\ell =\ell (\ell +1)C_\ell /2\pi$. Because the form of the MOG Friedmann equation (3) is the same as in the ${\Lambda }CDM$ model, we can match the $\Lambda$CDM calculation of the CMB angular power spectrum. The calculation of the CMB power spectrum is displayed in Fig.1. The calculation of the power spectrum in the ${\Lambda }CDM$ model is duplicated in MOG, using the Planck 2018 Collaboration best-fit parameter values [28]: $\Omega _bh^2=0.0224\pm 0.0001$, $\Omega _\phi h^2=0.120\pm 0.001$, $\Omega _\Lambda =0.680\pm 0.013$, $n_s=0.965\pm 0.004$, $\sigma _8=0.811\pm 0.006$, $H_0=67.4\pm 0.5\,\mathrm{km}\,\mathrm{sec}^{-1}\,\mathrm{Mpc}^{-1}$, together with the remaining parameters in the fitting process.

4 Late time universe

As the universe expands beyond the time of decoupling, the gravitational attraction between the baryons becomes enhanced as the parameter $\alpha$ increases in size and $G=G_\infty =G_N(1+\alpha )$. The density $\rho _\phi$ and the increasing enhancement of the size of G deepen the baryon gravitational potential well. Eventually, as the large scale structures form $\rho _\phi < \rho _b$ and the baryon dominated MOG dynamics takes over. Once the combined $\rho _b$ and $\rho _\phi$ density perturbations have grown sufficiently to produce stars and galaxies at about 400 million years after the big bang, the MOG non-relativistic dynamics for baryons takes over to determine the final evolution and dynamics of galaxies and galaxy clusters in the present epoch. The effective mass $\mu$ undergoes a decrease, and from the best-fit values from galaxy rotation curves and cluster dynamics [9, 10, 12], $\mu =0.01 - 0.04\,\mathrm{kpc}^{-1}$, corresponding to $m_\phi \sim 10^{-28}\,\mathrm{eV}$.

From (13) for a VG mass $m_\phi \sim 10^{-28}$ eV, the VG density in the present universe is $\Omega _\phi \le 10^{-3}$. From (14) the mean velocity of the VGs is now of order the speed of light c and the VGs are ultra-relativistic particles, so galaxy and galaxy cluster halos of VGs cannot form as stable bound systems. The de Broglie wavelength, in natural units, $\lambda _\phi \sim 1/m_\phi$, of the ultra-light VGs is of the size of galaxies. For the predictions of galaxy rotation curves and the galaxy cluster dynamics the best-fit values of $\alpha$ are determined by Eq. (A-22).

Fig. 2

After applying the appropriate window function, MOG (thick red line) shows agreement with the luminous red galaxy survey mass power spectrum data. The MOG fit to the data is comparable to the CDM prediction (thin blue line) [25, 26]

The calculation of the matter power spectrum for $\rho _b >\rho _\phi$ will have unit oscillations. However, the finite size of galaxy survey samples and the associated window function used to produce presently available power spectra mask any such oscillations. Applying a window function to the MOG prediction for the matter power spectrum smooths out the power spectrum curve. The enhanced size of $G=G_N(1+\alpha )$ with $\alpha \sim 19$ predicts the right shape for the galaxy matter power spectrum curve, resulting in a fit to the data [25, 26]. The GR prediction without dark matter and with $G=G_N$ cannot produce the correct magnitude or shape for the matter power spectrum. In Fig. 2, we show the predicted MOG matter power spectrum [25, 26].

In future galaxy surveys which utilize a large enough number of galaxies, with galaxies detected at sufficiently large redshift z, and with the use of a sufficiently narrow enough window function, it should be possible to detect any significant unit oscillations in the matter power spectrum, which can distinguish in the present universe between a dominant dark matter model and MOG without dark matter.

5 Conclusions

By assuming that the density of VGs associated with the massive and neutral vector field $\phi _\mu$ in MOG theory is the dominant density in the early universe with a VG mass, $m_\phi \sim 10^{-22}\,\mathrm{eV}$, the perturbations $\delta \rho _\phi$ satisfy the pressureless Jeans equation, allowing for an enhanced growth from the time of horizon entry and radiation-matter equality to produce large scale stellar and galaxy structure after the time of decoupling and recombination. The baryons obey the Jeans equation with baryon-photon pressure, so that the baryons oscillate between the time of horizon entry and decoupling, generating the baryon-pressure acoustical waves detected in the angular power spectrum in the CMB. The angular power spectrum calculation matches the $\Lambda$CDM model fit to the Planck collaboration 2018 data.

After the universe expands until the time of stellar and galaxy formation when $\rho _\phi < \rho _b$ and the mass of the VGs become ultra-relativistic with $m_\phi \sim 10^{-28}\,\mathrm{eV}$, the non-relativistic MOG Newtonian acceleration law, including the repulsive Yukawa interaction and the attractive gravitation due to the enhanced value of G, determines the rotation curves of galaxies and the dynamics of galaxy clusters.

The matter power spectrum in MOG can, with an appropriate window function, fit the galaxy matter power spectrum data. A critical test of MOG is whether significant baryon oscillations in the power spectrum begin to show as the number of observed luminous red galaxies increases and the size of the window function decreases. If the smooth ${\Lambda }CDM$ model fit to the matter power spectrum persists with a large enough increase in observed galaxies in galaxy redshift surveys, then this would rule out MOG in favor of the existence of dark matter halos in galaxies.

The ultra-light VG particles and their associated Proca field $\phi _\mu$ are an integral part of the MOG theory determined by the STVG action. The MOG theory can fit observational data from the solar system to galaxy and galaxy clusters, and the early universe cosmological data, and provides a gravitational description of the universe at both small and large distance scales.

Appendix A

The MOG action is given by

$$\begin{aligned} S=S_G+S_\phi +S_S+S_M, \end{aligned}$$

(A-1)

where $S_M$ is the matter action and

$$\begin{aligned} S_G&=\frac{1}{16\pi }\int d^4x\sqrt{-g}\left[ \frac{1}{G}(R+2\Lambda )\right] , \end{aligned}$$

(A-2)

$$\begin{aligned} S_\phi&=\int d^4x\sqrt{-g}\left[ -\frac{1}{4}B^{\mu \nu }B_{\mu \nu }+\frac{1}{2}\mu ^2\phi ^\mu \phi _\mu \right] , \end{aligned}$$

(A-3)

$$\begin{aligned} S_S&=\int d^4x\sqrt{-g}\left[ \frac{1}{G^3}\left( \frac{1}{2}g^{\mu \nu }\partial _\mu G\partial _\nu G-V_G\right) \right. \nonumber \\&\quad \left. +\frac{1}{\mu ^2G}\left( \frac{1}{2}g^{\mu \nu }\partial _\mu \mu \partial _\nu \mu -V_\mu \right) \right] , \end{aligned}$$

(A-4)

where $B_{\mu \nu }=\partial _\mu \phi _\nu -\partial _\nu \phi _\mu$ and $V_G$ and $V_\mu$ are potentials. Note that we choose units such that $c=1$, and use the metric signature $[+,-,-,-]$. The Ricci tensor is

$$\begin{aligned} R_{\mu \nu }=\partial _\lambda \Gamma ^\lambda _{\mu \nu }-\partial _\nu \Gamma ^\lambda _{\mu \lambda } +\Gamma ^\lambda _{\mu \nu }\Gamma ^\sigma _{\lambda \sigma }-\Gamma ^\sigma _{\mu \lambda }\Gamma ^\lambda _{\nu \sigma }. \end{aligned}$$

(A-5)

The matter stress-energy tensor is obtained by varying the matter action $S_M$ with respect to the metric:

$$\begin{aligned} T^{\mu \nu }_M=-2(-g)^{-1/2}\delta S_M/\delta g_{\mu \nu }\,. \end{aligned}$$

(A-6)

Varying $S_\phi$ and $S_S$ with respect to the metric yields

$$\begin{aligned} T^{\mu \nu }_\phi= & {} -2(-g)^{-1/2}[\delta S_\phi /\delta g_{\mu \nu }], \end{aligned}$$

(A-7)

$$\begin{aligned} T^{\mu \nu }_S= & {} -2(-g)^{-1/2}[\delta S_S/\delta g_{\mu \nu }]. \end{aligned}$$

(A-8)

Combining these gives the total stress-energy tensor

$$\begin{aligned} T^{\mu \nu }=T^{\mu \nu }_\mathrm{M}+T^{\mu \nu }_\phi +T^{\mu \nu }_S. \end{aligned}$$

(A-9)

We have

$$\begin{aligned} T_{\phi \mu \nu }=-\frac{1}{4\pi }\bigg [B^\alpha _\mu B_{\nu \alpha }-g_{\mu \nu }\bigg (B^{\rho \alpha }B_{\rho \alpha }+\frac{1}{2}\mu ^2\phi ^\mu \phi _\mu \bigg )\bigg ].\nonumber \\ \end{aligned}$$

(A-10)

The MOG field equations are given by

$$\begin{aligned}&G_{\mu \nu }-\Lambda g_{\mu \nu }+Q_{\mu \nu }=8\pi GT_{\mu \nu }\,, \end{aligned}$$

(A-11)

$$\begin{aligned}&\frac{1}{\sqrt{-g}}\partial _\mu \biggl (\sqrt{-g}B^{\mu \nu }\biggr )+\mu ^2\phi ^\nu =J^\nu , \end{aligned}$$

(A-12)

$$\begin{aligned}&\partial _\sigma B_{\mu \nu }+\partial _\mu B_{\nu \sigma }+\partial _\nu B_{\sigma \mu }=0. \end{aligned}$$

(A-13)

Here, $G_{\mu \nu }=R_{\mu \nu }-\frac{1}{2} g_{\mu \nu } R$ is the Einstein tensor and

$$\begin{aligned} Q_{\mu \nu }&=\frac{2}{G^2}(\partial ^\alpha G \partial _\alpha G\,g_{\mu \nu } - \partial _\mu G\partial _\nu G)\nonumber \\&\quad - \frac{1}{G}(\Box G\,g_{\mu \nu } - \nabla _\mu \partial _\nu G) \end{aligned}$$

(A-14)

is a term resulting from the the presence of second derivatives of $g_{\mu \nu }$ in R in $S_G$. We also have the field equations:

$$\begin{aligned} \Box G=K,\quad \Box \mu =L, \end{aligned}$$

(A-15)

where $\nabla _\mu$ is the covariant derivative with respect to the metric $g_{\mu \nu }$, $\Box =\nabla ^\mu \nabla _\mu$, $K=K(G,\mu ,\phi _\mu )$ and $L=L(G,\mu ,\phi _\mu )$.

Combining the Bianchi identities, $\nabla _\nu G^{\mu \nu }=0$, with the field equations (A-12) yields the conservation law

$$\begin{aligned} \nabla _\nu T^{\mu \nu }+\frac{1}{G}\partial _\nu G\,T^{\mu \nu } - \frac{1}{8\pi G}\nabla _\nu Q^{\mu \nu }=0 \,. \end{aligned}$$

(A-16)

It is a key premise of MOG that all baryonic matter possesses, in proportion to its mass M, positive gravitational charge: $Q_g=\kappa \,M$. This charge serves as the source of the vector field $\phi ^\mu$. Moreover, $\kappa =\sqrt{8\pi (G-G_N)}=\sqrt{8\pi \alpha G_N}$, where $G_N$ is Newton’s gravitational constant and $\alpha =(G-G_N)/G_N\ge 0$. Variation of $S_M$ with respect to the vector field $\phi ^\mu$, then yields the MOG 4-current $J_\mu =-(-g)^{-1/2}\delta S_M/\delta \phi ^\mu$.

For the case of a perfect fluid:

$$\begin{aligned} T^{\mu \nu }=(\rho +p)u^\mu u^\nu -pg^{\mu \nu }, \end{aligned}$$

(A-17)

where $\rho$ and p are the matter density and pressure, respectively, and $u^\mu$ is the 4-velocity of an element of the fluid. We obtain from (A-17) and $u^\mu u_\mu =1$ the 4-current:

$$\begin{aligned} J_\mu =\kappa T_{M\mu \nu }u^\nu =\kappa \rho u_\mu . \end{aligned}$$

(A-18)

It is shown in [23] that, with the assumption $\nabla _\mu J^\mu =0$, (A-16) reduces to

$$\begin{aligned} \nabla _\nu T_M^{\mu \nu }=B_\nu ^{~\mu } J^\nu \,. \end{aligned}$$

(A-19)

The field equations (A-12) and (A-13) are to be considered as effective field equations for the $B_{\mu \nu }$ field. The scalar field G is in effect a running coupling constant dependent on the binding energy or mass of a body. This also applies to the scalar field $\alpha$. The field equations for $B_{\mu \nu }$ do not solve to give unique universal values of $\alpha$ or $\mu$. The solutions for G or $\alpha$ depend on the renormalization group value of G analogous to the behavior of coupling constants in the standard model of particle physics.

The equation of motion of a particle is given by

$$\begin{aligned} \frac{d^2x^\mu }{ds^2}+{\Gamma _{\alpha \beta }}^\mu u^\alpha u^\beta =\frac{q_g}{m}{B^\mu }_\alpha u^\alpha , \end{aligned}$$

(A-20)

where $q_g=\kappa m$ and m is the particle mass. We have $q_g/m=\kappa$ and the equation of motion satisfies the equivalence principle. Due to the Lorentz-type gravitational force, massive particles do not free fall on geodesics. Photons have $m_\gamma =0$ and satisfy the null geodesic equation:

$$\begin{aligned} k^\mu \nabla _\mu k^\nu = 0, \end{aligned}$$

(A-21)

where $k^\mu$ is the 4-momentum vector of photons and $k^2=g^{\mu \nu }k_\mu k_\nu =0$. Gravitational radiation follows the same null geodesics as photons [17, 18], and gravitational waves travel with the speed of light.

For the low energy behavior of $G=G_N(1+\alpha )$ and $\mu$, we can adopt approximately constant values for G and $\mu$: $\partial _\mu G\approx 0$ and $\partial _\mu \mu \approx 0$, where $\alpha$ depends on the density and mass of a body. A phenomenological formula for $\alpha$ is:

$$\begin{aligned} \alpha =\alpha _{\inf }\frac{M}{(\sqrt{M}+E)^2}, \end{aligned}$$

(A-22)

where $\alpha _{inf}\sim \mathcal{O}(10)$ and $E=2.5\times 10^4\,M_{\odot }$. For $r\ll \mu ^{-1}\sim 25 -100$ kpc the MOG acceleration (12) reduces to the Newtonian acceleration ${\mathbf {a}}_\mathrm{Newt}$. This is consistent with $\alpha \sim 10^{-9}$ obtained from (A-22) for $M=1\, M_{\odot }$, guaranteeing that the solar system observational data is satisfied.