Cosmology of bigravity with doubly coupled matter

We study cosmology in the bigravity formulation of the dRGT model where matter couples to both metrics. At linear order in perturbation theory two mass scales emerge: an hard one from the dRGT potential, and an environmental dependent one from the coupling of bigravity with matter. At early time, the dynamics is dictated by the second mass scale which is of order of the Hubble scale. The set of gauge invariant perturbations that couples to matter follow closely the same behaviour as in GR. The remaining perturbations show no issue in the scalar sector, while problems arise in the tensor and vector sectors. During radiation domination, a tensor mode grows power-like at super-horizon scales. More dangerously, the only propagating vector mode features an exponential instability on sub-horizon scales. We discuss the consequences of such instabilities and speculate on possible ways to deal with them.


Introduction and Summary
The physical mechanism responsible for the present day acceleration of our universe is unknown. The simplest explanation is a positive cosmological constant; however the large amount of tuning required to fit the data seems excessive. Alternatives based on modifications of Einstein General Relativity (GR), which become observationally relevant at large scales, are being actively explored nowadays [1]. Among them, Massive Gravity [2] has received special attention: from an effective field theory perspective, it is one of the most natural options to investigate when renouncing to the diffeomorphism invariance of GR. In such a scenario, the graviton mass introduces a new energy scale that can be related with the scale of dark energy. In its simplest incarnation, to build a massive deformation of GR a reference non-dynamical metric is needed. Besides the aether-like nature of the reference metric, an unattractive feature from a theoretical perspective, there are various motivations to go beyond massive gravity and enter in the realm of bigravity theories [3,4,5,6]. For example, spatially flat homogenous Friedmann-Robertson-Walker (FRW) solutions do not exist [7] in Lorentz invariant ghost free massive gravity, and even allowing for open FRW solutions [8] strong coupling [9] and ghostlike instabilities [10] develop. Flat FRW solutions exist [11,12] in the case of Lorentz breaking models [13,14,15,16,17].
In bigravity one can find various branches of regular cosmological solutions describing flat FRW cosmologies [18,19]. A branch where the gravity modification is equivalent to an effective cosmological constant suffers of strong coupling [20]. A more promising branch is unstable at early time, being characterized by an exponential growth of fluctuations [20,21]. Also, singular FRW-type solutions exist [18] which exhibit only mild instabilities (power-law growth of vector and tensor modes [22,23]). On the other hand, they correspond to bouncing universes characterized by a naked curvature singularity, which makes their physical relevance questionable. Recently, it has been proposed to extend the theory of massive (bi)gravity by considering a more general coupling to matter -called doubly matter coupling -in which the physical metric coupled to the matter energy momentum tensor is an appropriate linear combination of the two metrics [24,25] (see also [26,27] for different approaches). This scenario, although problematic for massive gravity [28], is potentially interesting in the bigravity setup since a qualitatively new branch of FRW cosmological solutions exists [29], hence its cosmological perturbations deserve to be investigated. Such a theory is not ghost-free [30,31], but there exist physically interesting situations where the Boulware-Deser (BD) ghost does not represent an immediate phenomenological problem. This is the case if its mass is above the cut-off scale Λ c of the theory under consideration. In addition, Λ c might be parametrically larger than the strong coupling scale where the effects of the graviton mass term become important and interesting. Moreover, the ghost does not manifest itself at linear order in an expansion in fluctuations around particularly symmetric configurations (as for example FRW cosmologies). In this work, after studying the two branches of cosmological solutions at the homogeneous level, we focus on the dynamics of linearized cosmological fluctuations around the new homogenous backgrounds allowed by the doubly matter coupling. No hints of BD ghost mode are found at linear level in fluctuations, and the theory propagates the seven degrees of freedom as expected for a healthy bigravity theory. The dynamics of scalar fluctuations is healthy, and no instabilities are found in this sector. The dynamics of tensor and vector fluctuations is richer, but it shows problematic behaviours. The tensor sector exhibits a power-like instability at superhorizon scales during the radiation domination era. We argue that such instability is not extremely serious, and can be tamed by an appropriate choice of initial conditions, possibly motivated by inflation. Much worse is the behaviour of vector fluctuations. In this case, during radiation domination, we find a gradient instability at subhorizon scales, which leads to an exponential growth of small scale fluctuations, rapidly driving the theory outside the regime of validity of perturbation theory. Hence, this serious instability rules out the cosmological configurations that we consider. Nevertheless, we speculate on possible extension of the bigravity theory under consideration, that might be able to cure such instability problems.

The theory under consideration 2.1 Scalar Field
Before plunging in the study of bigravity, it is interesting to understand in a simplified setting the peculiar feature of the non-minimal coupling of matter to gravity proposed by Ref. [24]. We will show how the consistency of the effective description of matter as a (perfect) fluid necessarily requires the dynamical character of the second metric, selecting bigravity as the only consistent formulation. Take as matter a scalar field φ that couples with gravity not simply by the metric g µν but trough a combination of g µν and a non-dynamical flat metric f µν We denoted by the time derivative with respect to t. In order to restore diffeomorphism (diff) invariance, the non-dynamical metric f can be written using four Stuckelberg fields; in the unitary gauge we have Φ 0 = z(t) and Φ i = x j δ i j . The scalar field φ, with a potential F (φ), couples to gravity according to Thus, φ is minimally coupled to G µν defined by where Y µ ν = ( √ X) µ ν and X µ ν = g µσ f σν and α , β are two arbitrary constants. Notice that φ is not minimally coupled to g µν . Setting β = 0 and α = 1, we recover the standard minimal coupling. It is clear from the ADM canonical analysis that the total Hamiltonian is not anymore linear in the lapse N and shift N i of the dynamical metric g µν , then the number of propagating degrees of freedom (DoF) will be more than 3. The seemingly innocent action (2) actually represents a modification of gravity.
Let us consider FRW homogeneous cosmological solutions where and The Energy Momentum Tensor (EMT) for the scalar is diagonal, fluid-like and can be written as which has a peculiar dependence on N and z . The expression for the EMT reduces to the standard one when β → 0, α → 1. Contrary to the case of a scalar field minimally coupled to gravity, the time-time and the spatial components of the Einstein equations and the equation of motion for φ are all independent. Indeed, taking the time derivative of the time-time component of the Einstein equations and using the equation of motion of φ, one can solve for a ; then inserting this expression in the spatial components of the Einstein equations one gets the following constraint 1 Thus, unless β = 0, Einstein equations require that p φ = 0. The same constraint follows from the requirement that the scalar EMT is conserved. Of course such constraint has no counterpart in GR, where the EMT for φ is automatically conserved when φ satisfies its equation of motion. Hence: -in a FRW background the dynamics of a scalar field in (2) is not equivalent to a perfect fluid; -the pressure has to vanish.
Such difficulties were not taken into account in [32]. In the following we will show that both issues are absent when the non-dynamical Stuckelberg metric (1) is promoted to a full-fledged dynamical one, see also [28]. The reason for such a behaviour can be traced back to the non-dynamical nature of the metric f µν used in the new matter coupling. The conservation of the energy momentum tensor defined as the response of the matter action S matt to a diffeomorphism variation of the dynamical metric gives such a strong constraint that in general cannot be satisfied unless very special conditions like (7) are met. This is not very surprising and it is typical of theories with non-dynamical object [33,34].
When instead the metric f µν gets dynamical, inserting its own Ricci scalar in the action, for the following FRW parametrization and a new possibility of realizing (11) opens up. Within this new way, the scalar field dynamics can be still captured by the perfect fluid description and no spurious constraint is required. Therefore, for the rest of this paper, we will parametrize the matter content of the Universe through a perfect fluid and the metric f µν entering in (3) will become dynamical in the bigravity formulation.

Bigravity and Matter Coupling
Consider the action of massive bigravity with the dRGT potential as interaction between the two dynamical metrics g µν and f µν In the presence of two metrics, it is not a priori clear to what metric matter couples to. In general, the BD ghost revives in the presence of doubly coupled matter [30,31]. In [24] a new matter coupling was proposed where matter is minimally coupled to the effective metric G µν given by (3). Although the BD ghost persists even with this special doubly coupled bigravity model, it was shown that the BD ghost does not appear in the decoupling limit [24,30]. As for the scalar field non-minimally coupled to gravity of section 2.1, the theory described by the action (12) propagate more than the seven DoF expected in the bigravity formulation of dRGT. The extra scalar mode is most probably a ghost, however, the theory is still acceptable if the mass of such a mode is above the ultraviolet cutoff Λ c . Clearly this point deserves further investigation. The matter EMT is defined as the response of the matter action to a variation of g (f ) The modified Einstein equations can be written as where Q ( Q) are the effective energy-momentum tensors induced by the interaction term for the two metrics.
Let us introduce the BD ghost free potential [35,36,37] where the V n are the symmetric polynomials of Y with τ n = tr(Y n ). We have that where (V ) µ ν = ∂V /∂Y ν µ and q = det X = det(f )/ det(g).

FRW Ansatz and Conservation Laws
Let us consider homogeneous FRW bi-diagonal metrics with flat spatial slices in conformal time, so the form of g and f is as follows: The effective metric gets the following form Consistency of the equations of motion requires the following Bianchi-type constraints When T µν is the EMT of a perfect fluid, i.e. T µν = (p + ρ) u µ u ν + p G µν with u α G αβ u β = −1, the previous equations can be combined to give where H = a /a and H ω = ω /ω and the equation of state w, p = w ρ. Notice that condition (23) corresponds to the conservation of the matter EMT T µν with respect to the metric G µν . The constraint (24) can be realised in two inequivalent ways:

Branch 1
In this branch eq (24) is realised through the vanishing of the square bracket, i.e. the pressure of the fluid is determined by the massive potential. Notice that at early time, when ρ m 2 M 2 pl , consistence requires that w ρ ≈ 0, even though the Universe should be radiation dominated at that epoch. In this branch we have to deals with the very same issue found in section 2.1 and the description of matter as a perfect fluid is not consistent with the one of a scalar field. This branch was studied in [32], when f µν is a flat non-dynamical metric, in presence of a scalar field. Though, contrary to the case of massive gravity with minimally coupled matter, flat FRW solutions exist, the non-physical requirement of w = 0 makes the present branch not very appealing as discussed in Ref. [28].

Branch 2
In this case Notice that the limit β → 0 exists and we recover the very same branch of FRW cosmology in bigravity with standard matter minimally coupled to the metric g µν . Contrary to branch 1, since condition (25) is matter independent, the scalar field dynamics is equivalent to a perfect fluid; in this sense matter has the standard effective description. It is also interesting to note that (25) is equivalent to requiring that the matter's action and the interaction part of the action are separately diff invariant, namely and equivalently for the metric f . Clearly, the branch two is the most interesting one and from now on we will focus on it.

Background solutions for Branch 2
Introducing the ratio of the two scale factors ξ = ω/a, the tt-component of the modified Einstein equations for g reads The relation between ξ and a is determined by using (25) and (27) in the tt-component of the modified Einstein equations for f ; the result is an algebraic equation The case α → 0 and β → 0 has been already studied in [20] and instabilities were found in the scalar sector. Throughout this paper we assume that the scale of the graviton mass m is of the order of the present Hubble scale, i.e m 2 ∼ M −2 pl ρ Λ ∼ H 2 0 . Such a choice is the most natural one if massive gravity has anything to do with the present acceleration of the Universe. Then, according to the evolution equation for the matter energy density (23) that gives at early times, provided that ξ 1/a, we can always consider m 2 M 2 pl /ρ 1 as a dimensionless expansion parameter. Therefore, at early time, solutions of equation (28) can be classified according to the different regimes of ξ, i.e ξ 1, ξ 1, and ξ ∼ 1.
ξ 1 In this case the solution for the matter energy density is given by At the leading order we have ρ ∼ ρ Λ . This rules out the small ξ regime at early times.
ξ 1 In this case we have that at early times the matter energy density is given by The behaviour is equivalent to a cosmological constant; thus, also the large ξ regime is not suitable for early time cosmology.
ξ 1 For this last case we have two solutions for ξ, and at the leading order they read When ξ − α β the spatial components of the effective metric (21) are singular and moreover the early time cosmology is dominated by a cosmological constant, i.e. H 2 ∝ m 2 a 2 . When instead ξ β κ α , at the leading order we have and the corrections are of the order ∼ O(m 2 M 2 pl /ρ). Therefore, up to a renormalisation of the Newtonian constant 2 , the early time cosmology is very similar to GR at the background level.
Notice that the ξ ∼ 1 regime is incompatible with the β → 0 limit. 3 Indeed, when β is very small, ξ ∼ O m 2 M 2 pl /ρ and we turn back to the small ξ regime studied in [20]. Thus the early time cosmology for ξ β κ α exists only when β m 2 M 2 pl /ρ. In order to have this early time cosmology until the solution reaches the de Sitter point, β needs to be O(1), assuming that all other parameters are also O(1).
In the next section we will analyse cosmological perturbations around this background.

Cosmological Perturbations
Perturbations around a FRW background can be classified according to representations of the SO(3) group, namely scalar, vector and tensor modes. It is convenient to use the gauge invariant formulation following [20,21] and for the benefit of the reader the relevant expression are collected in Appendix A.
In order to avoid cluttering in the main text, the general expressions of the perturbed Einstein equations for the scalar and vector sector are given in Appendices B.1 and B.2, see also [39] for a different derivation.
The general feature that emerges is that comparing with the ghost free bigravity theory with minimally coupled matter [20], the new matter coupling gives rise to an additional environmental dependent mass; namely, the new equations can be obtained from the old ones by the replacement where J is a function of time that depends on the field under consideration. The bottom line is that two mass scales are present: a hard mass m 2 coming from the interaction potential for the metrics, see (12); and, as a consequence of the non-minimal coupling, a second "soft" running mass scale. In the physical branch, where ξ ∼ 1, at early times the soft mass always dominates on the hard one, since ρ/M 2 pl m 2 . Thus at early time the potential V plays no role and can be simply neglected at the leading order in ∼ O m 2 M 2 pl /ρ . Of course this is the case only when β = 0, otherwise the regime ξ ∼ 1 simply does not exist and the m 2 corrections become the leading part. The perturbed equations of motion can be fully solved in this approximation.
Notice that the new matter coupling is driven by the matter EMT T µν , so the contribution to the ij component of the perturbed Einstein equations is always proportional to w; as a result such a contribution vanishes when the Universe is dominated by non-relativistic matter where w = 0.
In the following we will often use the Fourier transform with respect to x i of the various perturbations, the corresponding 3-momentum will be denoted k i and k 2 = k i k i δ ij . To keep notations as simple as possible we will use the same name for the field and its Fourier transform.

Scalar Sector
In the scalar sector, see Appendix B.1, the fields E, B 1 and Ψ 1/2 are non dynamical and can be expressed in terms of Φ 1/2 , basically one of the Bardeen potentials for g and f , that satisfy two second order equations; thus 2 scalar DoF propagate. As already pointed out in section 2.1, by using canonical analysis an additional scalar is expected, however as matter of fact such a mode does not propagate around a homogeneous background. Nevertheless it is expected to appear in a less symmetric background and/or at higher order in perturbation theory [30]. The consequences of the missing BD mode, at leading order in cosmological perturbations around FRW cosmologies, deserve further study.
In order to solve for the two propagating fields, it is convenient to use the following combination of themselves: It is Φ + that couples to matter thus Φ + is the relevant metric perturbation for observations. Φ + turns out to decouple completely (at the leading order) and satisfies the same equation as in GR. Instead the equation of motion for Φ − is more complicated; at early times it can be expanded in powers of τ /τ U 1, where τ U is the age of the universe in conformal time. In what follows, we shall consider the leading order in such an expansion.
Consider first the case of a radiation dominated Universe. For Φ + we have the very same equation of GR, i.e.
For Φ − , outside the horizon where k τ 1 we get and for the modes inside the horizon, namely k τ 1, we have Both scalars show no instability and their non-decaying modes are constant outside the horizon. Inside the horizon instead they both oscillate with a different decaying amplitude and frequency. The corresponding energy density perturbation is at the leading order and therefore behaves like in GR. Outside the horizon it turns out to be constant, instead inside the horizon it has a fixed amplitude oscillating behaviour. In the case of matter dominated Universe both the fields Φ + and Φ − obey the same equation of GR at the leading order, namely It is worth to stress that the gauge-invariant density perturbation δρ gi /ρ is the real observable quantity and is given by δρ gi ρ = κ α 2 6 (κ α 2 + β 2 ) 12 + k 2 τ 2 Φ + .
Again the behaviour of δρ gi /ρ is similar to GR: outside the horizon it is frozen and inside the horizon it grows like τ 2 leading to the formation of structures.
In order to study the other Bardeen potentials, it is useful to introduce also the combinations In general it turns out that It is interesting to note that in bigravity the double diff invariance is broken down to the diagonal diff invariance when an interaction between the two metric is introduced. The + sector is protected by diagonal diffs and indeed in the tensor sector corresponds to massless spin 2 modes, while the − gives rise to the massive ones, see also [39].
Remarkably, the instability that was present in the field Φ 2 for modes inside the horizon with matter minimally coupled to g µν , found in [20], is not present. For such a behaviour the new matter coupling is instrumental, indeed the ξ ∼ 1 regime for the background emerges and the very fast gradient instabilities disappear. However, as we will see soon, instabilities in the tensor and vector sectors emerge.

Tensor Sector
The first hint of the problematic behaviour of fluctuations manifests itself when studying the tensor fluctuations. For this sector, the final expression for the perturbed Einstein equations is rather simple where f 1 ,ρ and d are defined in Appendix B.1. Also in this sector, it is useful to introduce the following combination of fields with the indices and the T T symbol understood. Again it is h + that is relevant for observed gravitational waves, since it is the combination appearing inside G µν . In particular, for the physical background solution, we get for general w The combination h + behaves as in GR and represents spin 2 masseless modes protected by diagonal diff invariance. Instead the combination h − features a contribution to its mass proportional to w (see the parenthesis in the first line of eq (49)). This contribution is entirely due to the new matter coupling. If w is positive, the mode h − is tachyonic. If w is negative, this mode acquires a positive mass squared. Let us discuss some physical implications of the evolution equation (49) for the mode h − . During matter domination, w = 0 and the new contributions to the effective mass to the mode h − vanish: there are no instabilities associated with the new coupling of gravity with matter.
In the radiation dominated era, w = 1/3: the mode h − acquires a tachyonic mass, and an instability is expected. Indeed, radiation dominated super-horizon solutions for h + and h − , in terms of the scale factor a ∝ τ , are where C 1 , C 2 , C 3 , C 4 are integration constants. The quantities C 1 , C 2 control the healthy evolution of the mode h + , and contribute also to the mode h − due to the source term in the second line of eq (49).
The integration constants C 2 , C 3 correspond to decaying modes, that can be neglected, while C 4 controls a growing mode. As a consequence, this system features a mild power-like instability at superhorizon scales, but only during radiation dominated era. Hence, in order to ensure that this set-up is under perturbative control, we need to impose that the amplitude of h − does not exceed unity during radiation domination. Let us choose units in which the scale factor after the reheating period, the beginning of radiation domination, is a rh = 1. Focus on the contribution of the growing mode only, and take into account the evolution of the scale factor during radiation domination. We find that at the end of the radiation period, at matter-radiation equality, we have to satisfy the inequality where the temperature of transition from radiation to matter domination is T eq 3 eV, while the reheating temperature when the radiation era starts is model dependent, and depends on the reheating temperature after inflation. Taking a representative value T rh = 10 8 GeV, we find an upper bound for the integration constant C 4 : Interestingly, such small values for the quantity C 4 , can be motivated by the inflationary phase that precedes radiation domination. Indeed, it is inflation that sets the initial condition for the amplitude of tensor fluctuations, that evolve in the radiation dominated era. We sketch here an argument to explain this fact. Inflation is a phase of quasi-de Sitter expansion, where the parameter w −1. During this phase, the mode h − acquires a positive mass squared in its evolution equation (49) with D 1 . . . D 4 integration constants. The quantities D 1 , D 2 control the evolution of the mode h + : neglecting the decaying mode D 2 , the constant mode h + matches continuously between the inflationary and radiation dominated era, hence we set D 1 = C 1 . The quantities D 3 , D 4 control the specific properties of the mode h − , so they determine the initial conditions for the mode C 4 at the beginning of radiation era, after the end of inflation. Notice that D 3 , D 4 are decaying modes. We can write the matching relation between solutions (51) and (55) where we make the crude assumption of a sudden transition from inflation to radiation domination, so that the end of inflation occurs at a rh = 1. Assuming that during inflation the size of h − remains bounded, and that there are no cancelations among the different terms in eq (55), at the beginning of inflation we have the condition Notice that in our units in which a rh = 1, the value of a in , the scale factor at the beginning of inflation, is very small. So in these units it is natural to choose a very small value for D 4 to satisfy condition (57). Indeed, saturating the previous equality, and using the relation a rh /a in = e N ef , with N ef the e-fold number (that we take for definiteness N ef = 60), recalling that in our units a rh = 1, we find that comfortably satisfies inequality (53).
Of course the previous arguments are based on linear perturbation theory. Non-linear effects couple different modes: higher order fluctuations in the scalar sector are expected to feed the initial value of the tensor amplitude at the beginning of radiation era. On the other hand, in estimating their effects, one would have to take in careful consideration the coupling factors between scalars and h − , that can be suppressed by factors of the graviton mass. So it is possible that even taking into account non-linearities, the growth of tensor fluctuations can still be maintained under control: it would be interesting to investigate more in detail this topic.

Vector Sector
The dynamics of fluctuations in the vector sector is unfortunately much more problematic. Using the conservation of the matter EMT T µν with respect to the effective metric G µν (equation (78) in Appendix A), the velocity perturbations δv i can be completely solved. From equations (88-91) in Appendix B.2, one can show that all vectors can be expressed in terms of a single combination V 12 that satisfies a second order equation; thus, only a single transverse vector propagates. As for the tensors, the qualitative features of vector dynamics depends on the equation of state of the fluid constituting the matter EMT. When the fluid equation of state p = w ρ is such that w ≤ 0, no instabilities are found and the system is healthy. Instead, we find serious gradient instabilities when w > 0. In particular, in the case of a radiation dominated universe we have that δv = δv 0 , where δv 0 is an arbitrary constant, and at the leading order in τ /τ U So we find the 'wrong sign' in front of the gradient term for the propagating mode V 12 , leading to a gradient instability in the sub-horizon limit kτ 1. In this regime, we have an exponential growth for V 12 , and the leading contribution to the solution of (61) is Also when selecting tuned initial conditions, such exponential instability is so severe that even a small initial amplitude of vector generated by non-linear effects (as the ones mentioned at the end of the previous section) will be rapidly amplified to a level that drives perturbation theory out of control. Notice that once more in the + sector no instability is present. Indeed, the combination represents a decreasing mode. Thus instabilities are present only in the − sector. To conclude, the dynamics of vector fluctuations ruin the cosmology of bigravity doubly coupled to matter. Possible ways-out, to be investigated in the future, could be to modify our ansatz for the EMT. Indeed in this work we mainly focussed on an EMT of perfect fluid form: it might be that other choices of matter content can lead to better behaved system. As an example, one could further include vector degrees of freedom, also non-minimally coupled to gravity (as recently explored in massive gravity or related systems in [40,41,42]), and study cosmological configurations in the scenario of doubly coupled matter. Acting as additional source, they might be able to fix the problems we found with vector fluctuations. We leave this interesting issues to a future investigation.

Conclusions
Massive gravity has been the subject of an extensive investigation. The main phenomenological motivation is to explain the present acceleration of the universe. However devising a satisfactory model is not an easy task. The simplest ghost free Lorentz invariant version of massive gravity [35] with an auxiliary non-dynamical metric has no flat homogenous solution [7] and only Lorentz breaking 5 models [13,14,15,16,17] can support such configurations [11]. Sticking with Lorentz invariant models, to overcome the limitations of a non-dynamical auxiliary metric, one can promote it to a dynamical one in the contest of bigravity. As far as homogeneous configurations are concerned, in bigravity theories the situation drastically improves with respect to massive gravity, and one finds branches of flat FRW solutions. Leaving aside the ones with a curvature singularity at late time with c < 0 (see eq. (20)) -as discussed in [18] and more recently in [43,22,23] -regular background solutions unfortunately features an exponential instability in the scalar sector at early time, quickly invalidating perturbation theory already during radiation domination [20]. Things do not change when both metrics are coupled with two different matter sectors [21]. In the presence of two metrics there is a certain degree of ambiguity on how to couple matter with gravity [24] and actually one can consider a sort of democratic coupling of the two metrics with matter, see eq. (3). Though the new coupling reintroduce the Boulware-Deser ghost [30], one can argue that its mass is above the cutoff and does not affect the low energy physics in the spirit of effective field theory. One of the interesting features of the new coupling is that it gives rise to an effective background dependent soft mass, see for instance (35) for the case of an FRW background. As a result even taking the "hard" mass m in the deforming potential to zero, we still have a massive gravity theory thanks to the environmental soft mass, that for an FRW background is proportional to the Hubble parameter. This is the feature that opens up a new dynamical regime compared with minimally coupled case such that ξ = ω/a stays constant at early time and then flows toward a de Sitter attractor responsible for the present acceleration phase, see figure 1. The next step is to study the behaviour of cosmological perturbations. Things get better in the scalar sector where no instability is found, however troubles develop in the tensor and specially in the vector sector. Among the tensor modes the + combination (see eq.(47)) which is protected by diagonal diffs has the same dynamics as in GR, while the second independent tensor mode develops a power-law growth until the matter domination. Such a growth is naturally counterbalanced by a sufficient low primordial production during inflation. In the vector sector the situation is more worrisome: while again the + combination of vector fields, see eq.(63), has only a decreasing mode, the − combination shows an exponential instability when w > 0 at subhorizon scales. The main effect is to loose theoretical control of the theory at the perturbative level already during the radiation domination. On the other hand, we stress that in all cases studied the anomalous growth is only present in sectors that do not couple directly with observed matter; moreover, we speculated on possible extension of the bigravity theory under consideration, that might be able to cure such instability problems.
The fields in (70) can be of course expressed in terms of the ones in (69). In the matter sector, since matter is minimally coupled to the effective metric G µν , we define the gauge invariant perturbed pressure and density in the following way where and For matter, together with pressure and density perturbation, there is also the perturbed 4-velocity u µ that consists of a scalar part v and a vector part δz i u µ =ū µ + δu µ , u µ u ν G µν = −1 .
The conservation of the matter EMT T µν with respect to the effective metric G µν leads to the following differential relation for vector matter perturbations In the vector sector we have 4 fields and 1 gauge transformation; thus, we can form 3 independent gauge invariant vector perturbations