Dipolar Dark Matter with Massive Bigravity

# Dipolar Dark Matter with Massive Bigravity

Luc Blanchet and Lavinia Heisenberg
###### Abstract

Massive gravity theories have been developed as viable IR modifications of gravity motivated by dark energy and the problem of the cosmological constant. On the other hand, modified gravity and modified dark matter theories were developed with the aim of solving the problems of standard cold dark matter at galactic scales. Here we propose to adapt the framework of ghost-free massive bigravity theories to reformulate the problem of dark matter at galactic scales. We investigate a promising alternative to dark matter called dipolar dark matter (DDM) in which two different species of dark matter are separately coupled to the two metrics of bigravity and are linked together by an internal vector field. We show that this model successfully reproduces the phenomenology of dark matter at galactic scales (i.e. MOND) as a result of a mechanism of gravitational polarisation. The model is safe in the gravitational sector, but because of the particular couplings of the matter fields and vector field to the metrics, a ghost in the decoupling limit is present in the dark matter sector. However, it might be possible to push the mass of the ghost beyond the strong coupling scale by an appropriate choice of the parameters of the model. Crucial questions to address in future work are the exact mass of the ghost, and the cosmological implications of the model.

NORDITA-2015-57

Dipolar Dark Matter

with Massive Bigravity

• Institut d’Astrophysique de Paris — UMR 7095 du CNRS, Université Pierre & Marie Curie, 98bis boulevard Arago, 75014 Paris, France

• Nordita, KTH Royal Institute of Technology and Stockholm University, Roslagstullsbacken 23, 10691 Stockholm, Sweden

• Department of Physics & The Oskar Klein Centre, AlbaNova University Centre, 10691 Stockholm, Sweden

E-mails: blanchet@iap.fr, laviniah@kth.se

## 1 Introduction

In the last century, cosmology has progressively developed from a philosophical to an empirical scientific discipline, witnessing high precision cosmological observations, which culminated with the standard model of cosmology, the -CDM model [1]. The standard model is based on General Relativity (GR) and is in great agreement with the abundance of the light elements in the big bang nucleosynthesis, the anisotropies of the cosmic microwave background (CMB), the baryon acoustic oscillations, lensing and the observed large scale structures. It notably relies on the presence of dark energy in form of a cosmological constant , giving rise to the accelerated expansion of the universe.

### 1.1 Massive gravity context

In the prevailing view, the cosmological constant corresponds to the constant vacuum energy density which receives large quantum corrections. Unfortunately, so far there is no successful mechanism to explain the observed unnatural tiny value of the cosmological constant, see e.g. [2, 3]. This difficulty has sparked a whole industry studying modifications of gravity in the infra-red (IR) invoking new dynamical degrees of freedom.

One tempting route is massive gravity, motivated by the possibility that the graviton has a mass. It was a challenge over forty years to construct a covariant non-linear theory for massive gravity. The foregathered remarkable amount of effort finally rose to the challenge [4, 5, 6, 7, 8, 9, 10, 11], to construct the potential interactions in a way that got rid of the Boulware–Deser (BD) ghost [12]. This theory has been further extended to more general models by adding additional degrees of freedom. A very important extension is massive bigravity theory [13].

Once the interactions in the gravitational sector were guaranteed to be ghost-free, a natural follow up question was how to couple this theory to the matter fields without spoiling the ghost-freedom. First attempts were already discussed in the original paper of bigravity [13]. If one couples the matter fields to both metrics simultaneously, this reintroduces the BD ghost [14, 15]. Furthermore, the one-loop quantum corrections detune the special potential interactions at an unacceptable low scale and hence this way of coupling is not a consistent one [15]. The safer way is to couple the matter fields to just one metric. In this way, the quantum corrections give rise to contributions in form of cosmological constants. One could also try to couple the matter field to the massless mode, which is unfortunately also not ghost free [16]. Another possible way of coupling the matter fields to both metrics is through a composite effective metric constructed out of both metrics [15, 17, 18, 19], which is unique in the sense that it is the only non-minimal matter coupling that maintains ghost-freedom in the decoupling limit [20, 21, 22]. Furthermore, the quantum corrections are guaranteed to maintain the nice potential structure. Other important consequence of this new matter coupling is the fact that it helps to evade the no-go result [23] for the flat Friedmann–Lemaître–Robertson–Walker (FLRW) background. A detailed perturbed ADM analysis revealed the existence of a BD ghost originating from an operator involving spatial derivatives [15, 17]. Therefore, the ghost will probably reappear for highly anisotropic solutions. Since the ghost remains absent up to the strong coupling scale, the matter coupling can be considered in an effective field theory sense at the very least till the strong coupling scale. The precise cut-off of the theory or mass of the ghost has still to be established.

The absence of the BD ghost is not only important at the classical level, but also at the quantum level. For the classical ghost-freedom the relative tuning of the potential interactions is the key point. Therefore one has to ensure that the quantum corrections do not detune the potential interactions. Concerning the decoupling limit, it is easy to show that the theory receives no quantum corrections via the non-renormalization theorem due to the antisymmetric structure of the interactions [24] (in fact the same antisymmetric structure of the Galileon interactions protect them from quantum corrections [25, 26, 27]). Beyond the decoupling limit, the quantum corrections of the matter loops maintain the potential interactions provided that the above criteria are fulfilled. Concerning the graviton loops, they do destroy the relative tuning of the potential interactions. Nevertheless, this detuning is harmless since the mass of the corresponding BD ghost is never below the cut-off scale of the theory [28]. The bimetric version of the theory shares the same property [29].

Massive gravity/bigravity theory has a rich phenomenology. The decoupling limit admits stable self-accelerating solutions [30], where the helicity-0 degree of freedom of the massive graviton plays the role of a condensate whose energy density sources self-acceleration. Unfortunately these solutions suffer from strong coupling issues due to the vanishing kinetic term of the vector modes [30, 31]. With the original massive gravity and the restriction of flat reference metric, one faces the no-go result, namely that there are no flat FLRW solution [23]. One can construct self-accelerating open FLRW solutions [32], which however have three instantaneous modes [33] and suffer from a nonlinear ghost instability [34]. The attempt to promote the reference metric to de Sitter [35, 36] also failed due to the presence of the Higuchi ghost [35]. In order to avoid these difficulties, one either gives up on the FLRW symmetries [23], or invokes new additional degrees of freedom [37, 13, 38].

Thanks to the freedom gained in the inclusion of the second kinetic term in the bimetric extension [13], there now exists many elaborate works concerning the cosmology of the bigravity theory, see e.g. [39, 40, 41]. In the case of minimally coupled matter fields and small graviton mass, the theory admits several interesting branches of solutions. Unfortunately, among them the self-accelerating branch is unstable due to the presence of three instantaneous modes, and a second branch of solutions admits an early time gradient instability [42]. Nevertheless, there exists attempts to overcome the gradient instability either by viable though finely tuned solutions in the case of a strongly interacting bimetric theory with  [43, 44], or by demanding that as was proposed in [45]. See also other recent works concerning the phenomenology of bimetric gravity [46, 47, 48, 49, 50].

### 1.2 Dark matter context

On a quite different but equally fascinating topic is the phenomenology of dark matter at galactic scales. The evidence for dark matter is through the measurement of the rotation curves of spiral galaxies which turn out to be approximately flat, contrary to the Newtonian prediction based on ordinary baryonic matter [51, 52]. The standard explanation is that the disk of galaxies is embedded into the quasi-spherical potential of a huge halo of dark matter, and that this dark matter is the same as the cold dark matter (CDM) which is evidenced at large cosmological scales notably with the fit of -CDM to the CMB anisotropies [53, 54]. Unfortunately this explanation faces severe challenges when compared to observations at galactic scales [55, 56]. There are predictions of the -CDM model that are not observed, like the phase-space correlation of galaxy satellites and the generic formation of dark matter cusps in the central regions of galaxies. Even worse, there are observations which are not predicted by -CDM, such as the tight correlation between the mass discrepancy (luminous vs. dynamical mass) which measures the presence of dark matter and the involved scale of acceleration, the famous baryonic Tully-Fisher (BTF) relation for spiral galaxies [57, 58, 59], and its equivalent for elliptical galaxies, the Faber-Jackson relation [60].

Instead of additional, non-visible mass, Milgrom [61, 62, 63] proposed an amendment to the Newtonian laws of motion in order to account for the phenomenology of dark matter in galaxies, dubbed MOND for MOdified Newtonian Dynamics. According to the MOND hypothesis the change has a relevant influence on the motion only for very small accelerations, as they occur at astronomical scales, below the critical value measured for instance from the BTF relation.111A striking observation is the numerical coincidence that  [55, 56]. A more elaborate version of MOND is the Bekenstein-Milgrom [64] modification of the Poisson equation of Newtonian gravity. All the challenges met by -CDM at galactic scales are then solved — sometimes with incredible success — by the MOND formula [55, 56]. Unfortunately, MOND faces difficulties in explaining the dark matter distribution at the larger scales of galaxy clusters [65, 66, 67, 68, 69]. It is also severely constrained by observations in the solar system [70, 71].

Reconciling -CDM at cosmological scales and MOND at galactic scales into a single relativistic theory is a great challenge. There has been extensive works on this. One approach is to modify gravity by invoking new gravitational degrees of freedom without the presence of dark matter [72, 73, 74, 75, 76, 77, 78, 79, 80, 81]. Primary examples are the tensor-vector-scalar (TeVeS) theory [73, 74] and generalized Einstein-Æther theories [75, 76]. Another approach is a new form of dark matter à la MOND, called dipolar dark matter (DDM). It is based on the dielectric analogy of MOND [82] — a remarkable property of MOND with possible far-reaching implications. Indeed the MOND equation represents exactly the gravitational analogue (in the non-relativistic limit) of the Gauss equation of electrostatics modified by polarisation effects in non-linear dielectric media. Some early realizations of this analogy were proposed in [83, 84, 85] and shown to also reproduce the cosmological model -CDM. The best way to interpret this property is by a mechanism of gravitational polarisation, involving two different species of particles coupled to different gravitational potentials. This was the motivation for pushing forward a more sophisticated model in the context of a bimetric extension of GR [86, 87]. In this model the two species of dark matter particles that couple to the two metrics are linked through an internal vector field. The model [87] is able to recover MOND at galactic scales and the cosmological model -CDM at large scales. Furthermore the post-Newtonian parameters (PPN) in the solar system are the same as those of GR. Unfortunately, even if this model delivers a very interesting phenomenology, it is plagued by ghost instabilities in the gravitational sector already at the linear order of perturbations and within the decoupling limit [88]. A more promising model was then proposed in [88], with the hope to cure these instabilities in the gravitational sector.

In the present work we will study in detail the theoretical and phenomenological consistency of this recently proposed model for dipolar dark matter in the context of bimetric gravity [88]. We will first work out the covariant field equations in Sec. 2 and analyze the linear field equations in Sec. 3. As next, we will investigate in Sec. 4 the mechanism of gravitational polarisation in the non relativistic approximation, and show how it successfully recovers the MOND phenomenology on galactic scales. Finally we will pay attention in Sec. 5 to the matter sector and investigate the interactions in the decoupling limit. The outcome of Sec. 5 is that a ghost instability is still present in the (dark) matter sector of the model. The paper ends with some short conclusions in Sec. 6.

## 2 The covariant theory

A new relativistic model for dipolar dark matter has been recently proposed in [88]. Basically this model is defined by combining the specific dark matter sector of a previous model [87] with gravity in the form of ghost-free bimetric extensions of GR [13]. The dark matter in this model consists of two types of particles respectively coupled to the two dynamical metrics and of bigravity. In addition, these two sectors are linked together via an additional internal field in the form of an Abelian vector field . This vector field is coupled to the effective composite metric of bigravity which is built out of the two metrics and , and which has been shown to be allowed in the matter Lagrangian in the effective field theory sense at least up to the strong coupling scale [15, 17]. The total matter-plus-gravity Lagrangian reads222The metric signature convention is . We adopt geometrical units with the gravitational constant, the Planck constant and the speed of light being unity, , unless specified otherwise.

 L =√−g(M2g2Rg−ρbar−ρg)+√−f(M2f2Rf−ρf) +√−geff[m2M2eff+Aμ(jμg−jμf)+λM2effW(X)]. (2.1)

Here and are the Ricci scalars of the two metrics, and the scalar energy densities of the ordinary matter modelled simply by pressureless baryons,333We have in mind that the baryons actually represent the full standard model of particle physics. and the two species of pressureless dark matter particles are denoted by , and respectively. In addition , stand for the mass currents of the dark matter, defined by

 jμg=√−g√−geffJμg% andjμf=√−f√−geffJμf, (2.2)

where and are the corresponding conserved dark matter currents associated with the respective metrics and , thus satisfying and . Here and denote the scalar densities and , are the four velocities normalized to and . We also introduced two constants and specifying the ratios between the charge and the inertial mass of the dark matter particles.

The metric independent matter degrees of freedom are the coordinate densities and , and the coordinate velocities and . The associated metric independent currents are and . They are conserved in the ordinary sense, and . They relate to the classical currents , or , by

 J∗μg=√−gJμg=√−geffjμgandJ∗μf=√−fJμf=√−geffjμf.

Notice that the coupling term in the action (2.1) is actually independent of any metric, neither nor nor . We shall study in detail in Sec. 5 the implications of the term for the decoupling limit of the theory.

The vector field is generated by the dark matter currents and plays the role of a “graviphoton”. The presence of this internal field is necessary to stabilize the dipolar medium and will yield the wanted mechanism of gravitational polarisation, as we shall see in Sec. 4. It has a non-canonical kinetic term , where is a function which is left unspecified at this stage, and

 X=−FμνFμν4λ. (2.3)

Here is a constant and the field strength is constructed with the effective composite metric given by (2.4) below, i.e. with not depending on the metric as usual.

The model (2.1) is defined by several constant parameters, the coupling constants , and , the mass of the graviton , the charge to mass ratios and of dark matter, the constant associated with the vector field, and the arbitrary constants and entering the effective metric (2.4). Parts of these constants, as well as the precise form of the function , will be determined in Sec. 4 when we demand that the model reproduces the required physics of dark matter at galactic scales. In particular will be related to the MOND acceleration scale , and the model will finally depend on and the mass of the graviton . Additionally, there will be still some remaining freedom in the parameters and , and in the coupling constants and .

In the model (2.1) the ghost-free potential interactions between the two metrics and take the particular form of the square root of the determinant of the effective metric [15, 17, 19]

 geffμν=α2gμν+2αβG% effμν+β2fμν, (2.4)

where and are arbitrary constants, and we have defined

 Geffμν=gμρXρν=fμρYρν, (2.5)

where the square root matrix is defined by , together with its inverse . Interchanging the two metrics and does not change the form of the metric (2.4)–(2.5) except for a redefinition of the parameters and . Notice that can be proved to be automatically symmetric [89, 16], and in [88] it was shown, that it corresponds to the composite metric considered in the previous model [87]. The square root of the determinant of is given in either forms respectively associated with the or metrics, by

 (2.6)

More explicitly, introducing the elementary symmetric polynomials and of the square root matrices or ,444As usual we denote the traces of rank-2 tensors as , , and so on. Recall in particular that , with and .

 e0(X) = 1, (2.7a) e1(X) = [X], (2.7b) e2(X) = 12([X]2−[X2]), (2.7c) e3(X) = (2.7d) e4(X) = 124([X]4−6[X]2[X2]+3[X2]2+8[X][X3]−6[X4]), (2.7e)

we can write, still in symmetric guise,

 √−geff=√−g4∑n=0α4−nβnen(X)=√−f4∑n=0αnβ4−nen(Y). (2.8)

In this form we see that (2.6) corresponds to the right form of the acceptable potential interactions between the metrics and .

We can first vary the Lagrangian (2.1) with respect to the two metrics and , which yields the following two covariant Einstein field equations [90]

 M2gY(μρGν)ρg−m2M2eff3∑n=0α4−nβnY(μρUν)ρ(n) =Y(μρ(Tν)ρbar+Tν)ρg)+α√−geff√−g(αY(μρTν)ρgeff+βTμνgeff), (2.9a) M2fX(μρGν)ρf−m2M2eff3∑n=0αnβ4−nX(μρVν)ρ(n) =X(μρTν)ρf+β√−geff√−f(βX(μρTν)ρgeff+αTμνgeff), (2.9b)

where and are the Einstein tensors for the and metrics. The tensors and are defined from the following sums of matrices [8]555Here we adopt a slight change of notation with respect to [8]. Our notation is related to the one of [8] by and .

 U(n)=n∑p=0(−)pen−p(X)Xp, (2.10a) V(n)=n∑p=0(−)pen−p(Y)Yp, (2.10b)

by raising or lowering indices with their respective metrics, thus and . By the same property which makes (2.5) to be symmetric, one can show that and are indeed automatically symmetric.

In the right sides of (2.9) the stress-energy tensors and are defined with respect to the metric , is defined with respect to and with respect to . Thus, for pressureless baryonic and dark matter fluids we have , and in terms of the corresponding scalar densities and normalized four velocities, while represents in fact the stress-energy tensor of the vector field , i.e.

 Tμνgeff=M2eff[WXFμρFννρ+λWgμνeff], (2.11)

where is the derivative of with respect to its argument defined by (2.3).

The equation of motion for the baryons follows a geodesic law , whereas the equations of motion for the dark matter particles are non-geodesic,

 agμ =rguνgFμν, (2.12a) afμ =−rfuνfFμν, (2.12b)

where the accelerations and are defined with their respective metrics. On the other hand, the equation for the vector field yields

 (2.13)

with denoting the covariant derivative associated with . This equation is obviously compatible with the conservation of the currents, since by (2.2) we have also and . The equation (2.13) can be also written as

 ∇νgeffTgeffμν=−(jνg−jνf)Fμν. (2.14)

Combining the equations of motion (2.12) with (2.14) we obtain the conservation law

 √−g∇νgTgμν+√−f∇νfTfμν+√−geff∇νgeffTgeffμν=0. (2.15)

Alternatively, such global conservation law can be obtained from the scalarity of the total matter action under general diffeomorphisms.

## 3 Linear perturbations

In this section we will be interested in computing the linear perturbations of our model, which will be extensively used in the following section for the study of the polarisation mechanism. Following the ideas of [87] we will assume that the two fluids of dark matter particles are slightly displaced from the equilibrium configuration by displacement vectors and . This makes the dark matter medium to act as an analogue of a relativistic plasma in electromagnetism. The dark matter currents will be slightly displaced from the equilibrium current as well, with satisfying . Furthermore, we will perturb the two metrics around a general background in the following way666By which we mean, for instance, .

 gμν =(¯gμν+hμν)2, (3.1a) fμν =(¯fμν+ℓμν)2, (3.1b)

with main interest in background solutions where . In the latter case the potential interactions take the form

 √−geff=√−¯g4∑n=0(α+β)4−nen(k), (3.2)

where , and traces are defined with the common background metric .

Working at first order in the displacement of the dark matter particles and assuming that their gradients are of the same order as the metric perturbations, the dark matter currents can then be expanded as [87]

 jμg =jμ0+2∇geffν(j[ν0⊥μ]ρyρg)+O(2), (3.3a) jμf =jμ0+2∇geffν(j[ν0⊥μ]ρyρf)+O(2), (3.3b)

with the projector operator perpendicular to the four velocity of the equilibrium given as , and the remainder indicating that the expansion is at first order. The desired plasma-like solution is achieved after plugging our ansatz (3.3) into the equations of the vector field (2.13), which results in

 WXFμν=−2M2effj[μ0ξν]+O(2), (3.4)

where we have defined the projected relative displacement as . Note that is necessarily a space-like vector. In Sec. 4 we shall see that the spatial components of this vector, which can be called a dipole moment, define in the non-relativistic limit the polarisation field of the dark matter medium as , where is the coordinate density associated with . But for the moment we only need to notice that is a first order quantity, therefore the field strength is itself a first order quantity, and hence the stress-energy tensor (2.11) with respect to is already of second order,777We are anticipating the form of given by (4.23) and which implies and .

 Tμνgeff=O(2). (3.5)

For this reason, in the case of our desired plasma-like solution the contributions of the self-interactions of the internal vector field will be at least third order in perturbation at the level of the action. For more detail see the comprehensive derivations in [87].

Let us first expand the Lagrangian (2.1) to first order in perturbation around a common background . We assume that there is no matter in the background, so the matter interactions do not contribute at that order. Using (3.2) we obtain (modulo a total divergence)

 L=L¯g+√−¯g[−(M2ghμν+M2fℓμν)G¯gμν+m2M2eff(α+β)3[k]]+O(2), (3.6)

where is the Einstein tensor for the background metric , and is the background value of the Lagrangian. In order for the background in our interest to be a solution of the theory, we have to impose that linear perturbations vanish. We find

 G¯gμν=m2M2effM2gα(α+β)3¯gμν=m2M2effM2fβ(α+β)3¯gμν. (3.7)

This is only possible if

 αM2g=βM2f. (3.8)

With this choice we guarantee that we are expanding around the background that is a solution to the background equations of motion.

As next we shall compute the action to second order in perturbations. Let us emphasize again that thanks to the plasma-like ansatz (3.3) the self-interaction of the internal vector field does not contribute to this order, see (3.5). Our quadratic Lagrangian above the background simply becomes

 L=L¯g +√−¯g[−M2ghρσ¯Eμνρσhμν−M2fℓρσ¯Eμνρσℓμν+hμν(Tμνbar+Tμνg)+ℓμνTμνf +12m2M2eff(α+β)2([k]2−[k2])]+O(3), (3.9)

where is the Lichnerowicz operator on the background given in the general case by

 −2¯Eρσμνhρσ =□¯g(hμν−¯gμνh)+∇¯gμ∇¯gνh−2∇¯g(μHν)+¯gμν∇¯gρHρ −2C¯gμρσνhρσ−23R¯g(hμν−14¯gμνh), (3.10)

where , , with denoting the Weyl curvature of the background metric. Of course this expression is to be simplified using the background equations (3.7). Also bear in mind that for having a common background . The field equations obtained by varying with respect to and yield

 −2M2g¯Eμνρσhρσ+Tμνbar+Tμνg+α(α+β)2m2M2% eff(¯gμνk−kμν)=O(2), (3.11a) −2M2f¯Eμνρσℓρσ+Tμνf+β(α+β)2m2M2eff(¯gμνk−kμν)=O(2), (3.11b)

where we recall that and . Of course those equations can also be recovered directly from the general Einstein field equations (2.9).

In the present paper we shall mostly make use of that combination of the two Einstein field equations (3.11) which corresponds to the propagation of a massless spin-2 field. We subtract the two equations from each other so as to cancel the mass term, resulting in

 −2¯Eμνρσ(M2gαhρσ−M2fβℓρσ)+1α(Tμνbar+Tμνg)−1βTμνf=O(2). (3.12)

However we shall also use the contribution of the mass term in the linearised Bianchi identities associated with (3.11). Thus, last but not least we can act on (3.11) with , yielding

 ∇¯gν(kμν−¯gμνk)=O(2). (3.13)

In deriving this relation we use the result that, for instance, which comes from the equations of motion (2.12) and the fact that there is no matter in the background. We shall see in the next section how important is this relation for the present model to recover the looked-for phenomenology of dark matter at galactic scales.

## 4 Polarisation mechanism and MOND

We next consider the Newtonian or non-relativistic (NR) limit of our model. In the previous section we derived the field equations at linear order around a vacuum background metric , which necessarily obeys the equations (3.7). Thus, the background is a de Sitter solution, with cosmological constant given by the graviton’s mass, . When computing the Newtonian limit applied to describe the physics of the local universe, e.g. the solar system or a galaxy at low redshift, we shall neglect the cosmological constant and shall approximate the de Sitter metric by a flat Minkowski background .

The best way to implement (and justify) this approximation is to perform an expansion on small scales, say . When we solve for the gravitational field in the solar system or a galaxy embedded in the de Sitter background, we get terms like (de Sitter-Schwarzschild solution), where is the mass of the Sun or the galaxy. On small scales the term will dominate over , so we can ignore the influence of the cosmological constant if the size of the system is . In the present case we shall require that the size of our system includes the very weak field region far from the system where the MOND formula applies. In that case the relevant scale is , where is the MOND acceleration. So our approximation will make sense provided that . For a galaxy with mass the MOND transition radius is . With the value of the cosmological constant we find , so we are legitimate to neglect the influence of the cosmological constant. For the solar system, while and the approximation is even better.

Actually we shall make the approximation only on the particular massless spin-2 combination of the two metrics given by (3.12), namely

 −2¯Eμνρσ(M2gαhρσ−M2fβℓρσ)=−1α(Tμνbar+Tμνg)+1βTμνf+O(2), (4.1)

where the Lichnerowicz operator (3) is now the flat one, and also in the constraint equation (3.13). We shall not need to consider the other, independent (massive) combination of the two metrics. In the NR limit when we parametrize the two metric perturbations and by single Newtonian potentials and such that

 h00=Ugc2+O(1c4), ℓ00=Ufc2+O(1c4), (4.2a) h0i=O(1c3), ℓ0i=O(1c3), (4.2b) ℓij=δijUfc2+O(1c4). (4.2c)

We shall especially pay attention to the potential since it represents the ordinary Newtonian potential felt by ordinary baryonic matter. Similarly we parametrize the internal vector field by a single Coulomb scalar potential such that

 A0 =ϕc2+O(1c4), (4.3a) Ai =O(1c3). (4.3b)

It is now straightforward to show that (4.1) reduces in the NR limit to a single scalar equation for a combination of the Newtonian potentials and ,888We shall adopt the usual boldface notation for ordinary three-dimensional Euclidean vectors. Also, from now on we no longer write the neglected remainders .

 (4.4)

where is the ordinary Laplacian, and , and are the ordinary Newtonian (coordinate) densities, satisfying usual continuity equations such as . Notably, all non-linear corrections in (4.1) are negligible. In principle, we should a priori allow some parametrized post-Newtonian (PPN) coefficients in the spatial components of the metrics (4.2), say and . However, when plugging (4.2) into (4.1) we determine from the components of (4.1) that in fact and . Similarly we find that the equations (2.13) governing the internal vector field reduce to a single Coulomb type equation,

 ∇⋅[WX∇ϕ]=1M2eff(rgρ∗g−rfρ∗f). (4.5)

Notice that the constants and in the effective metric to which is coupled the internal field intervene in the NR limit only in the expression

 X=|∇ϕ|22λ(α+β)4. (4.6)

Next, the equations of the baryons which are geodesic in the metric reduce to

 dvbardt=∇Ug, (4.7)

while the equations of motion (2.12) of the dark matter particles become

 dvgdt =∇Ug+rg∇ϕ, (4.8a) dvfdt =∇Uf−rf∇ϕ. (4.8b)

At this stage the important point is to recall that the two sectors associated with the metrics and do not evolve independently but are linked together by the constraint equation (3.13). We now show that this constraint provides a mechanism of gravitational polarisation in the NR limit and yields the MOND equation for the potential felt by baryons. In flat space-time this constraint is

 ∂ν(kμν−ημνk)=O(2), (4.9)

where and . Plugging (4.2) into (4.9) we readily obtain a constraint on the gravitational forces felt by the DM particles, namely

 ∇(αUg+βUf)=0. (4.10)

Note that the constraint (4.10) comes from a combination between the and components of the metrics (4.2). Using (4.10) we express the equations of motion (4.8) solely in terms of the potential ruling the baryons,

 dvgdt =∇Ug+rg∇ϕ, (4.11a) dvfdt =−αβ∇Ug−rf∇ϕ. (4.11b)

As we see, the gravitational to inertial mass ratio of the particles when measured with respect to the metric is .

We look for explicit solutions of (4.11) in the form of plasma-like oscillations around some equilibrium solution. We immediately see from (4.11) the possibility of an equilibrium (i.e. for which ) when we have the following relation between constants,

 αβ=rfrg. (4.12)

The equilibrium holds when , i.e. when the Coulomb force annihilates the gravitational force. To describe in the proper way the two DM fluids near or at equilibrium we use the NR limits of the relations (3.3), with an appropriate choice of the equilibrium configuration. From (4.11) we note that is constant, hence we define the equilibrium in such a way that the two displacement vectors and with respect to that equilibrium obey , and in particular we choose the equilibrium velocity to be . We shall now define the relative displacement or dipole moment vector by .999Here and denote the spatial components of the displacement four vectors and considered in Sec. 3, see (3.3)–(3.4). Hence (3.3) imply

 rgρ∗g =ρ∗0−βα+β∇⋅P, (4.13a) rfρ∗f =ρ∗0+αα+β∇⋅P, (4.13b)

where is the common density of the two DM fluids at equilibrium in the absence of external perturbations (far from any external mass), and is the polarisation. The velocity fields of the two fluids are

 (4.14)

where denotes the convective derivative, which reduces here to the ordinary derivative since .

By inserting (4.13) into (4.5) we can solve for the polarisation resulting in

 P=−M2effWX∇ϕ. (4.15)

Furthermore, by combining (4.14) and (4.11) and making use of the previous solution for the polarisation (4.15), we readily arrive at a simple harmonic oscillator describing plasma like oscillations around equilibrium, namely

 d2ξdt2+ω2gξ=α+ββ∇Ug. (4.16)

In this case the plasma frequency is given by

 ωg= ⎷α+ββrgρ∗0M2effWX. (4.17)

Next, consider the equation for the Newtonian potential in the ordinary sector. Combining (4.10) with (4.4) we readily obtain

 ΔUg=−12(M2g+α2β2M2f)[ρ∗bar+ρ∗g−αβρ∗f]. (4.18)

To recover the usual Newtonian limit we must impose, in geometrical units,

 M2g+α2β2M2f=18π. (4.19)

This condition being satisfied, the right-hand side of (4.18) can be rewritten with the help of (4.13) and the relation (4.12). We obtain an ordinary Poisson equation but modified by polarisation effects,

 ∇⋅[∇Ug−4πPrg]=−4πρ∗bar. (4.20)

We still have to show that the polarisation will be aligned with the gravitational field , i.e. we grasp a mechanism of gravitational polarisation. This is a consequence of the constitutive relation (4.15) taken at the equilibrium point, neglecting plasma like oscillations. At this point we have . Thus,

 P=M2effrgWX∇Ug. (4.21)

Finally (4.20) with (4.21) takes exactly the form of the Bekenstein-Milgrom equation [64]. To recover the correct deep MOND regime we must impose that when

 1−4πM2effr2gWX=|∇Ug|a0, (4.22)

where is the MOND acceleration scale. This is easily achieved with the choice

 W(X)=X−23(α+β)2X3/2+O(X2), (4.23)

together with fixing the constants and to the values

 M2eff=r2g4π,λ=a202. (4.24)

We thus recovered gravitational polarisation and the MOND equation (in agreement with the dielectric analogy of MOND [82]) as a natural consequence of bimetric gravity since it is made possible by the constraint equation (4.9) linking together the two metrics of bigravity.101010In the previous model [87] we had to assume that some coupling constant in the action tends to zero. With bimetric gravity and its nice potential interactions we see that no particular assumption is required. The polarisation mechanism is possible only if we can annihilate the gravitational force by some internal force, here chosen to be a vector field, and therefore assume a coupling between the two species of DM particles living in the and sectors.

Let us recapitulate the various constraints we have found on the parameters in the original action (2.1). These are given by (4.12) for having a polarisation process, (4.19) for recovering the Poisson equation and (4.24) for having the correct deep MOND regime. In addition, we recall (3.8) which was imposed in order to be able to expand the two metrics around the same background. Strictly speaking the relation (3.8) is not necessary for the present calculation because we neglected the influence of the background, and (3.8) may be relaxed for some applications. Finally, we note that the result can be simplified by absorbing the remaining constant together with the sum into the following redefinitions of the vector field and the effective metric: and . When this is done, and after redefining also the graviton mass , we see that we can always choose and without loss of generality.

Finally the fully reduced form of the action (2.1) reads

 Lfinal =√−g